Set Expression values according to gmsh physical regions

Hi, welcome to the FEniCS community!

It is pretty hard to guess what you would like to achieve without a MWE. I recommend you to read https://fenicsproject.discourse.group/t/read-before-posting-how-do-i-get-my-question-answered/21.

Assuming that you import the tags:

mesh = Mesh("MSL_XY.xml")
subdomains = MeshFunction("size_t", mesh, "MSL_XY_physical_region.xml")
boundaries = MeshFunction("size_t", mesh, "MSL_XY_facet_region.xml")

and redefine the measures as

dx = Measure('dx', domain=mesh, subdomain_data= subdomains)
ds = Measure('ds', domain=mesh, subdomain_data=boundaries)
dS = Measure('dS', domain=mesh, subdomain_data=boundaries)

Then you can use this to define your expression for each subdomain, i.e.,

def eps1(u):
    return 0.5*(grad(u) + grad(u).T)

def eps2(u):
    return 0.6*(grad(u) + grad(u).T)

I1 = inner(eps1(u),eps1(v))*dx(1)
I2 = inner(eps2(u),eps2(v))*dx(2)

where dx(1) and dx(2) are the volumetric measures for subdomains tagged as 1 and 2, respectively. Similarly for ds and dS.

You could also define a custom expression:

class Mu(UserExpression):
    """
    This class defines a variable mu that takes different values
    on the unitsquare domain (0,1)^2 at x<0.5 and x>0.5. The mesh is 
    assumed to be conforming in x=0.5
    """

    def __init__(self, p1, p2, **kwargs):
        super().__init__(**kwargs)
        self.p1 = p1
        self.p2 = p2

    def eval(self, value, x):
        if x[0] < 0.5. + DOLFIN_EPS:
            value[0] = self.p1(x)
        else:
            value[0] = self.p2(x)

    def value_shape(self):
            return ()


p1 = Expression('x[1]', degree=1)
p2 = Expression('x[0]+x[1]', degree=1)

mu = Mu(p1, p2, degree=2)
2 Likes