On 27 May 2016, at 14:05, Lawrence Mitchell <lawrence.mitchell@imperial.ac.uk> wrote:
Hi Tomasz,
On 24 May 2016, at 08:47, David Ham <david.ham@imperial.ac.uk> wrote:
Hi Tomasz,
There are two things here. The first one is cell subdomains. There is a pull request in at PETSc right now which will make this work to some extent. We'll let you know once the PETSc people accept that pull.
This is now incorporated in Firedrake. If your Gmsh file has physical region ids then you can use those numbers in the dx measure to restrict integration to the given region.
e.g.
mesh = Mesh("...")
assemble(Constant(1)*dx(1, domain=mesh))
assembles a Constant over the region of the mesh tagged with the marker "1".
Using this information to construct the appropriate indicator functions is a little more work, but we can help with getting this right.
And now if you update again, you can initialise your indicator functions like so: Assume that the fluid domain is marked in the mesh with region 1, and the solid domain with region 2: from firedrake import * mesh = Mesh("my_mesh.msh") V = FunctionSpace(mesh, "DG", 0) fluid_indicator = Function(V) solid_indicator = Function(V) par_loop(""" for ( int i=0; i < f.dofs; i++ ) { f[i][0] = 1.0; }""", dx(1), {'f': (fluid_indicator, WRITE)}) par_loop(""" for ( int i=0; i < s.dofs; i++ ) { s[i][0] = 1.0; }""", dx(2), {'s': (solid_indicator, WRITE)}) Now fluid_indicator will be 1 in the fluid part of the domain and 0 elsewhere and vice-versa for the solid_indicator function. For the syntax of these par_loop calls, please see the relevant documentation: http://firedrakeproject.org/firedrake.html#module-firedrake.parloops Cheers, Lawrence