Dolfinx: dirichlet BCs for mixed function spaces

dolfinx
1

I have two functions, u1 and u2, and I am trying to apply a Dirichlet boundary condition to just u1. I define the function space and functions as follows:

from mpi4py import MPI
from dolfinx import fem, io, nls, log, mesh, plot
import numpy as np
import ufl
from petsc4py.PETSc import ScalarType

domain = mesh.create_interval(MPI.COMM_WORLD, 100, [0, L])
P1 = ufl.FiniteElement('CG', domain.ufl_cell(), 1)
element = ufl.MixedElement([P1, P1])
V = fem.FunctionSpace(domain, element)

# test functions
v1, v2 = ufl.TestFunctions(V)
# fields to solve for
u = fem.Function(V)
u1, u2 = u.split()

And then I apply the boundary condition:

# function determining if a node is on the tray top
def on_top_boundary(x):
    return(np.isclose(x[0], L))

# determine boundary DOFs
boundary_dofs = fem.locate_dofs_geometrical(V.sub(0), on_top_boundary)
# apply dirichlet BC to boundary DOFs
bc = fem.dirichletbc(ScalarType(C_0), boundary_dofs, V.sub(0))

But this returns the error:
RuntimeError: Cannot tabulate coordinates for a FunctionSpace that is a subspace.

So instead I use the mixed space instead of a subspace:

# function determining if a node is on the tray top
def on_top_boundary(x):
    return(np.isclose(x[0], L))

# determine boundary DOFs
boundary_dofs = fem.locate_dofs_geometrical(V, on_top_boundary)
# apply dirichlet BC to boundary DOFs
bc = fem.dirichletbc(ScalarType(C_0), boundary_dofs, V)

But this returns: RuntimeError: Cannot locate dofs geometrically for mixed space. Use subspaces.

It seems like I’m being sent in circles. Any thoughts on what to try next?

2

How about:

from mpi4py import MPI
from dolfinx import fem, io, nls, log, mesh, plot
import numpy as np
import ufl
from petsc4py.PETSc import ScalarType

L = 1
domain = mesh.create_interval(MPI.COMM_WORLD, 100, [0, L])
P1 = ufl.FiniteElement('CG', domain.ufl_cell(), 1)
element = ufl.MixedElement([P1, P1])
V = fem.FunctionSpace(domain, element)

# test functions
v1, v2 = ufl.TestFunctions(V)
# fields to solve for
u = fem.Function(V)
u1, u2 = u.split()
u.name = "u"
# function determining if a node is on the tray top
def on_top_boundary(x):
    return(np.isclose(x[0], L))

V0, submap = V.sub(0).collapse()

# determine boundary DOFs
boundary_dofs = fem.locate_dofs_geometrical((V.sub(0),V0), on_top_boundary)


# apply dirichlet BC to boundary DOFs
C_0 = 5
bc = fem.dirichletbc(ScalarType(C_0), boundary_dofs[0], V.sub(0))

u.x.array[bc.dof_indices()[0]] = bc.value.value

with io.XDMFFile(domain.comm, "u.xdmf", "w") as xdmf:
    xdmf.write_mesh(domain)
    xdmf.write_function(u.sub(0))
    xdmf.write_function(u.sub(1))

As a mixed space can consist of spaces that cannot tabulate degrees of freedom (such as N1Curl), one has to use the collapsed subspace and its map to the parent space to get the degrees of freedom.

2 Likes
3

Sorry to bother

For the same case, if I want to apply Dirichletbc for u2, just replace the V.sub(0) to V.sub(1)? Should I change any index in boundary_dofs[0]?

Thanks

4

Dear @dokken ,

This script is working for me, but if instead of using a scalar function space I use a vector function space, I get the following error:

RuntimeError: Cannot tabulate coordinates for a mixed FunctionSpace.

Below there is a minimal version of the code.

from mpi4py import MPI
from dolfinx import fem, io, nls, log, mesh, plot
import numpy as np
import ufl
from petsc4py.PETSc import ScalarType

L = 1
domain = mesh.create_interval(MPI.COMM_WORLD, 100, [0, L])
P1 = ufl.VectorElement('CG', domain.ufl_cell(), 1)
element = ufl.MixedElement([P1, P1])
V = fem.FunctionSpace(domain, element)

# test functions
v1, v2 = ufl.TestFunctions(V)
# fields to solve for
u = fem.Function(V)
u1, u2 = u.split()
u.name = "u"
# function determining if a node is on the tray top
def on_top_boundary(x):
    return(np.isclose(x[0], L))

V0, submap = V.sub(0).collapse()

# determine boundary DOFs
boundary_dofs = fem.locate_dofs_geometrical((V.sub(0),V0), on_top_boundary)


# apply dirichlet BC to boundary DOFs
C_0 = 5
bc = fem.dirichletbc(ScalarType(C_0), boundary_dofs[0], V.sub(0))

I’ve just changed the line

P1 = ufl.FiniteElement('CG', domain.ufl_cell(), 1)

by

P1 = ufl.VectorElement('CG', domain.ufl_cell(), 1)

Is there a solution for this issue?

Thanks in advance.