Set domain definition from gmsh created mesh

Question

Hi, I am solving a bimaterial domain problem. For that i am creating mesh from gmsh and want use domain information in the FEniCS. Is there any elegant way to use the domain information. Here i am attaching .geo file.

Thanks for your help

Point(1) = {0, 0, 0, 1.0};
Point(2) = {1, 0, 0, 1.0};
Point(3) = {1, 1, 0, 1.0};
Point(4) = {0, 1, 0, 1.0};
Point(5) = {0, 0.5, 0, 1.0};
Point(6) = {1, 0.5, 0, 1.0};
Line(1) = {1, 2};
Line(2) = {2, 6};
Line(3) = {6, 3};
Line(4) = {3, 4};
Line(5) = {4, 5};
Line(6) = {5, 1};
Line(7) = {5, 6};
Line Loop(8) = {4, 5, 7, 3};
Plane Surface(9) = {8};
Line Loop(10) = {7, -2, -1, -6};
Plane Surface(11) = {10};

Answer 1

This is a gmsh question, rather than fenics.
An example, if you add:

Physical Surface('top') = 9;
Physical Surface('bottom') = 11;

to your .geo dolfin-convert will create an extra region.xml containing the subdomain information which you can load like this in python:

domains = MeshFunction("size_t", mesh, 'file_region.xml')

Comments

Thanks for the help and time.

I have to set different value of material properties in both the domain, how can we do that I am little bit confuse here.
for example in top domain E = 10 and bottom domain E = 1.

---

Have a look at the tutorial

---

I have created function as given in tutorial but I am unable to print value at any point it is giving error: Unable to evaluate expression. can you correct me where I am going wrong ?

k_0 = 3.0
k_1 = 10.0
class K(Expression):
    def __init__(self, materials, k_0, k_1, **kwargs):
        self.materials = materials
        self.k_0 = k_0
        self.k_1 = k_1

    def eval_cell(self, values, x, cell):
            if self.materials[cell.index] == top:
                values[0] = self.k_0
            else:
                values[0] = self.k_1

kappa = K(materials, k_0, k_1, degree=1)
print kappa(0.5,0.5) 

---

I am not entirely sure why this does not behave like a normal expression.
As a quick fix, you could do something like this:

kappa = K(materials, k_0, k_1, degree=0)
V = FunctionSpace(mesh,'DG',0)
kappa = project(kappa,V)
plot(kappa,interactive=True)

Note that 'top' is not the marker in your materials meshfunction, it is 1 or 2 in this case

---

Thank for the help, Now it is working. :)