A bug appears when using locate_dofs_topological(...)

General
1

Hi all! I have a question when I try to find the boundary dofs:
I consider a rectangular domain with length L and height d. I give the minimal code below

import dolfinx
from dolfinx import nls
import ufl 
import numpy as np
from mpi4py import MPI
from petsc4py.PETSc import ScalarType
import dolfinx.nls.petsc
from dolfinx import geometry
import matplotlib.pyplot as plt
from dolfinx import fem, io, mesh
import math
from ufl import (Measure, SpatialCoordinate, TestFunctions, TrialFunctions,
                 div, exp, inner)

# Define domain parameters
L = 0.1  # total length
delta = 2.0  # aspect ratio
d=L *delta  # thickness

# Create domain
my_domain = dolfinx.mesh.create_rectangle(MPI.COMM_WORLD, [np.array([0.0, 0.0]), np.array([L, d])], [5,9], cell_type=dolfinx.mesh.CellType.quadrilateral)

# Define mixed mesh
u_el = ufl.VectorElement('Lagrange', my_domain.ufl_cell(), degree=1, dim=2)
phi_el = ufl.FiniteElement('Lagrange', my_domain.ufl_cell(), degree=1)
mixed_el = ufl.MixedElement([u_el, phi_el]) # mixed element
V = dolfinx.fem.FunctionSpace(my_domain, mixed_el) 

#find the facets for the left and right boundaries
fdim = my_domain.topology.dim - 1
facets_left = mesh.locate_entities_boundary(my_domain, fdim, lambda x: np.isclose(x[0], 0.0))
facets_right = mesh.locate_entities_boundary(my_domain, fdim, lambda x: np.isclose(x[0],L))
Q, _ = V.sub(0).collapse()

#find the dofs
dofs_left_u = fem.locate_dofs_topological((V.sub(0), Q), fdim, facets_left)
dofs_right_u = fem.locate_dofs_topological((V.sub(0), Q), fdim, facets_right)

print(dofs_left_u)
print(dofs_right_u)

with the results

[array([  0,   1,   2,   3,  12,  13,  24,  25,  39,  40,  57,  58,  78,
        79,  96,  97, 114, 115, 132, 133], dtype=int32), array([  0,   1,  12,  13,  24,  25,  36,  37,  48,  49,  60,  61,  72,
        73,  84,  85,  96,  97, 108, 109], dtype=int32)]
[array([ 72,  73,  74,  75,  93,  94, 111, 112, 129, 130, 147, 148, 159,
       160, 168, 169, 174, 175, 177, 178], dtype=int32), array([ 10,  11,  22,  23,  34,  35,  46,  47,  58,  59,  70,  71,  82,
        83,  94,  95, 106, 107, 118, 119], dtype=int32)]

My understanding is that the above is the index of dofs at the left and right boundary. But I realize that the number “73” and “74” appear both at the right and left boundary. Could you please explain why and how I fix this error?
Thanks in advance for your help!

2

Note that each call returns two arrays. The first array are the indices in V, while the second array are the corresponding indices in Q.
This is done to map the bc from a function in Q to the Mixed space V

1 Like