The ufl multidomain infrastructure is now in next and will soon migrate to master. I would recommend from now on to avoid calls like assemble(form, mesh=mesh), and instead associate the mesh with the form properly.
The mesh to integrate over will be associated with the form if functions and geometric quantities in the integrand are all defined with the same mesh, so there is no ambiguity, or the mesh is attached to the measure with the "domain" keyword. Examples: from dolfin import * mesh = UnitSquareMesh(3, 3) V = FunctionSpace(mesh, "CG", 1) f = Function(V) a = f*dx assert a.domains()[0].data() == mesh n = FacetNormal(mesh) a = n[0]*ds assert a.domains()[0].data() == mesh a = 1*dx(mesh) assert a.domains()[0].data() == mesh a = 1*dx(domain=mesh) assert a.domains()[0].data() == mesh The domain_data keyword to the measure is recommended over the dx[markers] syntax: subdomain_markers = None # MeshFunction... dx1 = dx(domain=mesh, domain_data=subdomain_markers) a = 1*dx1 assert a.domains()[0].data() == mesh If your integrand contains some quantities defined on cells the old way, things should work out fine as long as there is only one mesh with the same properties in there: x = SpatialCoordinate(mesh) a = x[0]*f*dx assert a.domains()[0].data() == mesh x = SpatialCoordinate(mesh.ufl_cell()) a = x[0]*f*dx assert a.domains()[0].data() == mesh x = SpatialCoordinate(mesh.ufl_cell()) a = x[0]*dx(mesh) assert a.domains()[0].data() == mesh In UFL files there is no mesh, so we can use the ufl Domain class instead. The only use case for this so far is this little thing: D = Domain(triangle) M = 1*dx(D) but in a while we'll hopefully have more advanced multidomain features in place on top of this. Scream out if you have problems with the new code, I'm pretty sure it's robust now but there were a lot of subtleties involved here. Martin
_______________________________________________ fenics mailing list [email protected] http://fenicsproject.org/mailman/listinfo/fenics
