### Assign mixed element function values on subdomain

157

views

0

Short question: Given a

Details:

The following solves the lid-driven cavity with mixed finite elements (tested with FEniCS 2017.2.0):

Here's the resulting velocity vector field for reference:

Assume that I need to assign the velocity local to some subdomain by the lid, while preserving the other guess values from our previous solution. I could for example define the subdomain with

I remember seeing some related posts that play with degree of freedom maps for simpler cases; but I haven't seen a mixed finite element example, and asking the user to figure out the mapping for this case seems a bit complicated. I sure hope I'm missing something more straightforward.

`fenics.Function`

from a `fenics.FunctionSpace`

made from a `fenics.MixedElement`

, how does one assign a subset of the function values restricted to a subdomain, i.e. without changing values outside of the subdomain?Details:

The following solves the lid-driven cavity with mixed finite elements (tested with FEniCS 2017.2.0):

```
import fenics
mesh = fenics.UnitSquareMesh(20, 20)
P2P1 = fenics.MixedElement([
fenics.VectorElement('P', mesh.ufl_cell(), 2),
fenics.FiniteElement('P', mesh.ufl_cell(), 1)])
W = fenics.FunctionSpace(mesh, P2P1)
solution = fenics.Function(W)
u, p = fenics.split(solution)
v, q = fenics.TestFunctions(W)
mu = fenics.Constant(0.01)
inner, dot, grad, div, sym = \
fenics.inner, fenics.dot, fenics.grad, fenics.div, fenics.sym
F = (dot(v, dot(grad(u), u)) - div(v)*p \
+ 2.*mu*inner(sym(grad(v)), sym(grad(u)))
- q*div(u))*fenics.dx
near = fenics.near
class Lid(fenics.SubDomain):
def inside(self, x, on_boundary):
return on_boundary and near(x[1], 1.)
class FixedWalls(fenics.SubDomain):
def inside(self, x, on_boundary):
return on_boundary and (near(x[0],0.) | near(x[0],1.) | near(x[1],0.))
bcs = [
fenics.DirichletBC(W.sub(0), (1., 0.), Lid()),
fenics.DirichletBC(W.sub(0), (0., 0.), FixedWalls())]
fenics.solve(F == 0., solution, bcs)
```

Here's the resulting velocity vector field for reference:

```
import matplotlib
%matplotlib inline
u, p = fenics.split(solution)
fenics.plot(u)
```

Now I want to use most of this solution as the initial guess for a different boundary value problem. For example, I could change the lid velocity.

Assume that I need to assign the velocity local to some subdomain by the lid, while preserving the other guess values from our previous solution. I could for example define the subdomain with

```
class TopLayer(fenics.SubDomain):
def inside(self, x, on_boundary):
return x[1] > 0.9
```

**How do I change the values of the solution only in this subdomain?**For now I only need to assign a constant value, and it's okay to create a new function (from the same function space) rather than assigning the values in place.I remember seeing some related posts that play with degree of freedom maps for simpler cases; but I haven't seen a mixed finite element example, and asking the user to figure out the mapping for this case seems a bit complicated. I sure hope I'm missing something more straightforward.

Community: FEniCS Project

### 1 Answer

3

@klunkean came up with a hack that solves this problem in the discussion here: https://www.allanswered.com/post/mqpbr/#bamzr

So continuing from the script in my question, we can set the velocity and pressure in the subdomain with

which appears to give the right result:

So continuing from the script in my question, we can set the velocity and pressure in the subdomain with

```
hack = fenics.DirichletBC(W.sub(0), fenics.Constant((0., 0.)), TopLayer())
hack.apply(solution.vector())
hack = fenics.DirichletBC(W.sub(1), fenics.Constant(0.), TopLayer())
hack.apply(solution.vector())
```

which appears to give the right result:

`fenics.plot(u)`

```
fenics.plot(p)
fenics.plot(mesh)
```

To avoid hacking; perhaps the better answer lies within the `fenics.DirichletBC.apply`

code.

Edit: Removed comment about this not working with adapted meshes. One just has to call `leaf_node()`

in all of the right places.

The BCs are homogenised if you solve a nonlinear problem.

written
10 weeks ago by
Nate

Perhaps I shouldn't have added my comment about

`fenics.NonlinearVariationalProblem`

without showing exactly what I meant by "corrupt". The same issue should arise with linear problems. It's not the values of the BC's that get "corrupted", but rather the data structure gets ruined somehow and an error is thrown about the BC's being in the wrong space.
written
10 weeks ago by
Alexander G. Zimmerman

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

`fenics.NonlinearVariationalProblem`

, rather than calling`fenics.solve`

as done above, this seems to somehow corrupt the boundary conditions of the problem. This can be avoided by creating new temporary copies of the function space and solution on which to call`DirichletBC.apply`

, e.g.