This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

Integrating a new UFL-Operator (Connectivity with MultiFunction?)

+1 vote

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)
asked Jun 30, 2016 by sschmidt FEniCS Novice (490 points)

1 Answer

+1 vote
 
Best answer

It seems that the expected handler name for Div_Tan generated by ufl doesn't match the one defined in TestOperator.

After renaming Div_Tan to match the ufl naming convention your code produces the right output for me. That is, the only change was writing

class DivTan(CompoundDerivative):

(alternatively, name the handler in TestOperator to div__tan)

answered Jun 30, 2016 by Magne Nordaas FEniCS Expert (13,820 points)
selected Jul 1, 2016 by johannr

Great!! Thank you so much!!

...