dolfin team mailing list archive

[Johan Hake <johan.hake@xxxxxxxxx>]

Couldn't you just assemble the boundary condition directly from the solution?

  neuman_vector = assemble(u*vv*ds(bottom),\
                  exterior_facet_domains=mesh_function)

If this is what you want you can preassemble a boundary mass matrix:

  Mb = assemble(uu*vv*ds(bottom), exterior_facet_domains=mesh_function)

and then each times step:

  neuman_vector = Mb*u.vector()

Johan

On Saturday November 13 2010 12:58:06 David Maxwell wrote:
> New question #133939 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/133939
> 
> I'd like to extract Neumann data from a solution to an elliptic PDE.
> 
> For example, if u is a weak solution of -Laplacian(u) = f, then its Neumann
> data (\partial_n u) on the boundary is well defined as follows:
> 
> Given sufficiently regular v defined on the boundary, extend v to a
> function in H^1 in the domain.  Then
> 
> <\partial_n u, v> = \int_\Omega \nabla u \nabla v - f v
> 
> Since v belongs to L^2 of the boundary, I can interpret \partial_n as an
> element of L^2 of the boundary.
> 
> In truth, I really only want the Neumann data on part of the boundary.
> 
> I thought I might proceed as follows.  In the following, 'u' is the
> solution of the PDE in function space V, and mesh_function is a mesh
> function on the edges that equals 1 on the 'top' boundary, 2 on the
> remainder of the boundary, and 0 on all other edges.
> 
> My first attempt went something like:
> 
> uu = TestFunction(V)
> vv = TrialFunction(V)
> 
> inside=0; top=1; bottom=2
> 
> a = uu*vv*ds(bottom)
> rhs = inner(grad(u),grad(vv))*dx - vv*f*dx
> 
> top_bc = DirichletBC(V,Constant(0),mesh_function,top)
> interior_bc = DirichletBC(V,Constant(0),mesh_function,inside) # I thought I
> was being clever here!
> 
> u_n =
> VariationalProblem(a,rhs,bcs=[top_bc,interior_bc],exterior_facet_domains=m
> esh_function).solve()
> 
> This crashes on applying the boundary conditions using 'ident', so I tried
> setting the 'use_ident' parameter to False on the boundary conditions.
> 
> This still doesn't work because the interior_bc is (rightfully) too
> enthusiastic -- every vertex is joined to some interior edge and so I pick
> up the zero solution.
> 
> I tried some other tacts, but they are perhaps too foolish to describe
> here.
> 
> And I have the sense that I really ought to be solving a problem on a
> BoundaryMesh (rather than a highly constrained problem on the original
> mesh) but I don't know how to extend Functions on a boundary mesh to
> Functions on the original mesh.
> 
> Help!
> 
> David Maxwell
> 
> 
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
> 
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to     : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp