Dear all,
I would like to define a new differentiation operator as a customized UFL Operator. In particular the divergence in tangent spaces.
It can be defined by
div_tan(V,n) = div(V) - inner(dot(grad(V),n),n)
Implementing this in a naive way makes processing the form further very difficult, as I would have to look for matching pairs of terms and check if they make up the div_tan expression.
So, I very much rather have div_tan(V,n) as a new UFL-Operator with two operands.
Attached is a minimal example where this already works quite nicely. But when I try to identify the operation of the Form as in the UFL-Handbook, my new operator is always identified as "expression" and not by "div_tan".
A minimal example is attached, any help is much appreciated:
from dolfin import *
from ufl.core.ufl_type import ufl_type
from ufl.differentiation import CompoundDerivative
from ufl.assertions import ufl_assert
from ufl.checks import is_cellwise_constant
from ufl.corealg.multifunction import MultiFunction
@ufl_type(num_ops=1, inherit_indices_from_operand=0, is_terminal_modifier=True)
class Div_Tan(CompoundDerivative):
__slots__ = ()
def __new__(cls, f):
ufl_assert(not f.ufl_free_indices,"Free indices in the divergence argument is not allowed.")
# Return zero if expression is trivially constants
if is_cellwise_constant(f):
return Zero(f.ufl_shape[:-1]) # No free indices asserted above
return CompoundDerivative.__new__(cls)
def __init__(self, f):
CompoundDerivative.__init__(self, (f,))
@property
def ufl_shape(self):
return self.ufl_operands[0].ufl_shape[:-1]
def __str__(self):
return "div_tan(%s)" % self.ufl_operands[0]
def __repr__(self):
return "Div_Tan(%r)" % self.ufl_operands[0]
def div_tan(f):
return Div_Tan(as_ufl(f))
mesh = UnitSquareMesh(10,10)
FV = VectorFunctionSpace(mesh, "CG", 1)
V = Function(FV)
V.rename("V", "")
#This works as expected
Form1 = div(V)
print "Should say div:"
print "Form 1:", Form1
print "Should say div_tan:"
Form2 = div_tan(V)
print "Form 2:", Form2
#now see if we can identify the ufl operators of Form1 and Form2:
class TestOperator(MultiFunction):
def __init__(self):
MultiFunction.__init__(self)
def expr(self, o):
return "General Expression"#, o.ufl_operands
def div(self, o):
return "UFL-Build-in standard divergence"#, o.ufl_operands
def div_tan(self, o):
return "Custom tangent-divergence"#, o.ufl_operands
print "Now try to identify the forms:"
MyTest = TestOperator()
print MyTest(V)
print MyTest(Form1)
#Should report as "div_tan", but returns as "general expression" :(
print "The following should report as div_tan, but reports as general expression:"
print "I suspect the Multifunction environment needs to be told there is a new operator?"
print MyTest(Form2)