Define function from an expression with a variable parameter

Question

Dear Fenics users,

I need to define a 2-d vector function from an analytical expression with a variable parameter. I am currently using the class Expression with

class MyFunctionExpression(Expression):
    def eval(self, values, x):
        values[0] = 1+((x[0]**2 +x[1]**2)**0.5 - 0.5)**2
        values[1] = 1+((x[0]**2 +x[1]**2)**0.5 - 0.5)**2
    def value_shape(self):
        return (2,)

V = VectorFunctionSpace(MyMesh, "Lagrange", 1)
f = MyFunctionExpression(element=V.ufl_element())
f_u = interpolate(f,V)

Now I need to make the expression in MyFunctionExpression dependent on a parameter p that I will change during computation... like in

ParamExpression = Expression('param',param=Constant(0.),degree = 0)

that I update with

ParamExpression.user_parameters["param"]  = Constant(1.)

I don't know how to adapt this last syntax with my previous definition of MyFunctionExpression...

Many thank in advance for your help !

Claire

Answer 1

Hi, consider

from dolfin import *

class MyFunctionExpression(Expression):
    def __init__(self, param, **kwargs):
        self.param = param

    def eval(self, values, x):
        values[0] = self.param+((x[0]**2 +x[1]**2)**0.5 - 0.5)**2
        values[1] = self.param+((x[0]**2 +x[1]**2)**0.5 - 0.5)**2

    def value_shape(self):
        return (2,)

mesh = UnitSquareMesh(10, 10)
V = VectorFunctionSpace(mesh, "Lagrange", 1)
f = MyFunctionExpression(param=1, element=V.ufl_element())

for p in (1, 2, 3):
    f.param = p
    fu = interpolate(f, V)
    print p, '->', sqrt(assemble(inner(fu, fu)*dx)) 

Answer 2

For this simple case I would write expression as cpp string, i.e.

import dolfin as d

expr = d.Expression(("p + pow(pow(pow(x[0], 2) + pow(x[1], 2), 0.5) - 0.5, 2)",
                     "p + pow(pow(pow(x[0], 2) + pow(x[1], 2), 0.5) - 0.5, 2)"),
                     degree = 4, p = 1.0)

mesh = d.UnitSquareMesh(32, 32)
V = d.VectorFunctionSpace(mesh, "CG", 1)

expr_interp_1 = d.interpolate(expr, V)

expr.p = 0.0
expr_interp_0 = d.interpolate(expr, V)