1D mesh from 2D mesh subdomain (boundary)

On Fri, Nov 27, 2009 at 03:22:40PM -0000, Patrick Riesen wrote:
> New question #91944 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/91944
>
> Hello all,
>
> I'd like to extract the points of a 2d mesh boundary subdomain as a 1D mesh
> for separate computation of initial/boundary conditions.

You can do this simply by

  boundary_mesh = BoundaryMesh(mesh)

> Using the SubMesh class crashed in dolfin::TopologyComputation::compute_entities with a segfault.
>
> Is it supposed to work to get a 1D mesh from 2D with SubMesh?
> Or do i need to take another approach to extract a 1D mesh from 2D?
>
> thanks for your help

The problem with the above (BoundaryMesh) is that it does give you the
boundary mesh which is topologically 1D but it is embedded in R^2
which means you'll have trouble assembling over it since DOLFIN
assumes nD meshes to be embedded in R^n.

The SubMesh thing is intended for different use.

--
Anders

Ok, so BoundaryMesh(mesh) gives me the complete boundary. could i then extract from this a subboundary mesh?
or can i avoid the problem about assembling if i orientate my 2D mesh so that the subboundary i want to extract is a straight segment along an axis (i.e. x-axis)?

many thanks for the answer, i'll see what i can do :-)

best regards,
patrick

> You can, but you'll have to do it manually:
>
> 1. Extract the boundary mesh
>
> 2. Create a new mesh which is has both topological and geometric
> dimension 2.
>
> 3. Manually build the 2D mesh by copying data over from the boundary
> mesh.
>
> Use a MeshEditor for 2 and 3.

thank you for the reply, Anders.
I thought about something like this before i saw SubMesh/BoundaryMesh.
I'll use this approach.

many thanks,
patrick

> --
> Anders
>

Thanks Anders Logg, that solved my question.

I assume the above discussion regarding dolfin assembly over meshes embedded in higher degree space (ie 1d mesh in R^2),
is why the below code only prints the x[0] correctly, with x[1] incorrectly taking the value 0,0.5 or 1?

from dolfin import *

class MyExp(Expression):
    def eval(self, values, x):
        print x # <-----------------prints x[0], x[1]: but x[1] is not true value
        values[0] = 0

mesh = UnitSquare(10,10)
bmesh = Mesh(BoundaryMesh(mesh))
for i,j in [(i,j) for i in range(2) for j in range(2)]:
    bmesh.init(i,j)

V = VectorFunctionSpace(bmesh, "CG", 2)

exp = MyExp()
f = Function(V)
f.interpolate(exp)