### 1D boundary integrals and marking

156

views

3

Hello everyone,

I'm trying to implement a Robin boundary condition for a toy-problem in 1D.

I have a 1D domain and thus two boundaries:

The Left Boundary is marked as '1', and the Right boundary is marked as '2'

Now, the problem arises when I try to construct a weak form (this is just a general idea, not the full code):

Cheers, Petr.

PS. MWE: https://gist.github.com/Corwinpro/0c32046e56023e1e1dd26dc4b57360bb
Answer: ds should be returned in def Markboundaries(). Otherwise, it is just a local instance.

I'm trying to implement a Robin boundary condition for a toy-problem in 1D.

I have a 1D domain and thus two boundaries:

```
mesh = IntervalMesh(100, 0., 1.)
class boundaryLeft(SubDomain):
def inside (self, x, on_boundary):
return on_boundary and near(x[0], 0.)
class boundaryRight(SubDomain):
def inside (self, x, on_boundary):
return on_boundary and near(x[0],1.)
boundaryLeft = boundaryLeft();
boundaryRight = boundaryRight();
def MarkBoundaries(mesh):
sub_domains = MeshFunction("size_t", mesh, mesh.topology().dim()- 1)
boundaryLeft.mark(sub_domains,1)
boundaryRight.mark(sub_domains,2)
ds = Measure('ds', domain=mesh, subdomain_data=sub_domains)
```

The Left Boundary is marked as '1', and the Right boundary is marked as '2'

Now, the problem arises when I try to construct a weak form (this is just a general idea, not the full code):

```
p = TrialFunction(V)
q = TestFunction(V)
a = -inner(grad(p), grad(q))*dx + (q*p)*ds
```

As you can see, I have a boundary term "q*p*ds" here. If I use ds(1) or ds(2) or (ds(1) + ds(2)) instead of ds, it seems like the boundary term just doesn't appear in the weak form, i.e. the multiplication by ds(1),ds(2) just zeroes the boundary term.

UPD: I also tried VertexFunction, FacetFunction, etc. However, these don't work as well. Even if I just write:

`sub_domains.set_all(0)`

and then use ds(0) instead of ds, it zeroes the boundary term.

My question is:**Why do I get different answers (or different weak forms) if I use (ds(1) + ds(2)) or ds(0) instead of ds in my weak formulation ?**It seems like I do something wrong with the definition of MeshFunction, but I can't explain why and how to avoid this.

Cheers, Petr.

PS. MWE: https://gist.github.com/Corwinpro/0c32046e56023e1e1dd26dc4b57360bb

Community: FEniCS Project

written
7 months ago by
Corwinpro

Please login to add an answer/comment or follow this question.