DOLFIN

computation of Neumann boundary data

Asked by David Maxwell on 2010-11-13

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=mesh_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

Question information

Language:: English Edit question

Status:: Answered

For:: DOLFIN Edit question

Assignee:: No assignee Edit question

Last query:: 2010-11-15

Last reply:: 2010-11-15

Link existing bug

Revision history for this message

Johan Hake (johan-hake) said on 2010-11-15:

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 : <email address hidden>
> Unsubscribe : https://launchpad.net/~dolfin
> More help : https://help.launchpad.net/ListHelp

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@lists.launchpad.net
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp

Revision history for this message

David Maxwell (damaxwell) said on 2010-11-15:

This would be too easy. I'm looking for a function defined on (a part of) the boundary, not just its average value.

-David

On Nov 15, 2010, at 8:33 AM, Johan Hake wrote:

> Your question #133939 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/133939
>
> Status: Open => Answered
>
> Johan Hake proposed the following answer:
> 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 : <email address hidden>
>> Unsubscribe : https://launchpad.net/~dolfin
>> More help : https://help.launchpad.net/ListHelp
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https://answers.launchpad.net/dolfin/+question/133939/+confirm?answer_id=0
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/dolfin/+question/133939
>
> You received this question notification because you are a direct
> subscriber of the question.

This would be too easy. I'm looking for  a function defined on (a part of) the boundary, not just its average value.

-David

On Nov 15, 2010, at 8:33 AM, Johan Hake wrote:

> Your question #133939 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/133939
> 
>    Status: Open => Answered
> 
> Johan Hake proposed the following answer:
> 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@lists.launchpad.net
>> Unsubscribe : https://launchpad.net/~dolfin
>> More help   : https://help.launchpad.net/ListHelp
> 
> -- 
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https://answers.launchpad.net/dolfin/+question/133939/+confirm?answer_id=0
> 
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/dolfin/+question/133939
> 
> You received this question notification because you are a direct
> subscriber of the question.

Revision history for this message

Johan Hake (johan-hake) said on 2010-11-15:

On Monday November 15 2010 09:57:36 David Maxwell wrote:
> Question #133939 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/133939
>
> Status: Answered => Open
>
> David Maxwell is still having a problem:
> This would be too easy. I'm looking for a function defined on (a part
> of) the boundary, not just its average value.

Ok.

Restricting your FunctionSpace to a geometrical subdomain is not supported
(yet) in DOLFIN.

You should, however, be able to hack the solution (not pretty, and maybe not
totaly correct wrt FEM analysis) by solving a reduces system if you are using
'CG' 1 elements.

lhs = assemble(uu*vv*ds(bottom), exterior_facet_domains=mesh_function)
rhs = assemble(inner(grad(u),grad(vv))*dx - vv*f*dx)

verts = VertexFunction('uint', mesh)

# Mark the vertices at the boundary with the bottom marker
inds = (verts.array()==bottom)

# Get the reduces system
rhs0 = rhs[inds] # Approximation of the rhs contribution
lhs0 = lhs[inds,inds] # Just getting rid of zeros
x = lhs0.copy()

solve(lhs0,x,rhs0)

# Put reduced solution on the full solution function
u.vector()[:] = 0.0
u.vector()[inds] = x

plot(u)
interactive()

Johan

> -David
>
> On Nov 15, 2010, at 8:33 AM, Johan Hake wrote:
> > Your question #133939 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/133939
> >
> > Status: Open => Answered
> >
> > Johan Hake proposed the following answer:
> > 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> = a_\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 : <email address hidden>
> >> Unsubscribe : https://launchpad.net/~dolfin
> >> More help : https://help.launchpad.net/ListHelp
>
> 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 : <email address hidden>
> Unsubscribe : https://launchpad.net/~dolfin
> More help : https://help.launchpad.net/ListHelp

On Monday November 15 2010 09:57:36 David Maxwell wrote:
> Question #133939 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/133939
> 
>     Status: Answered => Open
> 
> David Maxwell is still having a problem:
> This would be too easy. I'm looking for  a function defined on (a part
> of) the boundary, not just its average value.

Ok.

Restricting your FunctionSpace to a geometrical subdomain is not supported 
(yet) in DOLFIN.

You should, however, be able to hack the solution (not pretty, and maybe not 
totaly correct wrt FEM analysis) by solving a reduces system if you are using 
'CG' 1 elements.

lhs = assemble(uu*vv*ds(bottom), exterior_facet_domains=mesh_function)
 rhs = assemble(inner(grad(u),grad(vv))*dx - vv*f*dx)

verts = VertexFunction('uint', mesh)

# Mark the vertices at the boundary with the bottom marker
 inds = (verts.array()==bottom)

# Get the reduces system
 rhs0 = rhs[inds] # Approximation of the rhs contribution 
 lhs0 = lhs[inds,inds] # Just getting rid of zeros
 x = lhs0.copy()

solve(lhs0,x,rhs0)

# Put reduced solution on the full solution function
 u.vector()[:] = 0.0
 u.vector()[inds] = x

plot(u)
 interactive()

Johan

> -David
> 
> On Nov 15, 2010, at 8:33 AM, Johan Hake wrote:
> > Your question #133939 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/133939
> > 
> >    Status: Open => Answered
> > 
> > Johan Hake proposed the following answer:
> > 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> = a_\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@lists.launchpad.net
> >> Unsubscribe : https://launchpad.net/~dolfin
> >> More help   : https://help.launchpad.net/ListHelp
> 
> 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@lists.launchpad.net
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp

Can you help with this problem?

Provide an answer of your own, or ask David Maxwell for more information if necessary.

To post a message you must log in.

Ask a question

Edit question

DOLFIN

computation of Neumann boundary data

Question information

Related bugs

Related FAQ:

Can you help with this problem?

Subscribers