Associating mixed element function space dofs with entries in Jacobian matrix

dolfinx
1

Hello,

I have the following mixed element system including its Jacobian J.

from dolfinx import mesh
from ufl import SpatialCoordinate, FiniteElement, MixedElement, TestFunctions, dx, split, inner, grad, derivative
from dolfinx.fem import FunctionSpace, Function, petsc, form
from mpi4py import MPI

msh = mesh.create_interval(MPI.COMM_WORLD, 2, points=(0, 1))
x = SpatialCoordinate(msh)
CG1_elem = FiniteElement("CG", msh.ufl_cell(), 1)
ME_elem = MixedElement([CG1_elem, CG1_elem])
ME = FunctionSpace(msh, ME_elem)
u = Function(ME)
u0, u1 = split(u)
u0_test, u1_test = TestFunctions(ME)
F = (-inner(grad(u0), grad(u0_test))) * dx
F += (-inner(grad(u1), grad(u1_test))) * dx

problem = petsc.NonlinearProblem(F, u)
J = derivative(F, u)
A = petsc.assemble_matrix(form(J))
A.assemble()
A_dense = A.convert("dense")
print(A_dense.getDenseArray())

I now want to see for each of the matrix elements of A_dense (or A), to which subspace it belongs. E.g. by creating another matrix of the same dimension as the Jacobian, but with entries either 0 or 1, where 0 would indicate that this entry belongs to u0 and 1 that it belongs to u1. I suspect that this can be done using the function dolfinx.fem.dofmap. However, I was not able to get it to work. I need this because I want to be able to see the (potential) block structure of my system. Thank you in advance for your help.

Cheers

domi

2

How do you consider whether an entry in A “belongs” to a specific space? Such nomenclature could be defined for the diagonal, but what about the off-diagonal?

3

Hello nate,

my MWE results in the following matrix:

[-2.  2.  0.  0.  0.  0.]
[ 2. -4.  0.  0.  2.  0.]
[ 0.  0. -2.  2.  0.  0.]
[ 0.  0.  2. -4.  0.  2.]
[ 0.  2.  0.  0. -2.  0.]
[ 0.  0.  0.  2.  0. -2.]

However, when I use

F += 5*(-inner(grad(u1), grad(u1_test))) * dx

I end up with

[ -2.   2.   0.   0.   0.   0.]
[  2.  -4.   0.   0.   2.   0.]
[  0.   0. -10.  10.   0.   0.]
[  0.   0.  10. -20.   0.  10.]
[  0.   2.   0.   0.  -2.   0.]
[  0.   0.   0.  10.   0. -10.]

From this I conclude that the entries are related to the subspaces as follows:

[  u0   u0   0.   0.   0.   0.]
[  u0   u0   0.   0.   u0   0.]
[  0.   0.   u1  u1    0.   0.]
[  0.   0.   u1  u1    0.   u1]
[  0.   u0   0.   0.   u0   0.]
[  0.   0.   0.  u1    0.   u1]

Since my indirect approach becomes cumbersome when the mixed space holds more than two elements (and the mesh has more then two elements) I am looking for a direct way to relate the individual subspaces to the matrix entries.

4

With:

V1, ME_V1_map = ME.sub(0).collapse()
V2, ME_V2_map = ME.sub(1).collapse()

You can find the dofs in ME corresponding to the first and second subspace. They are encoded in ME_V1_map and ME_V2_map. Is that what you need?

Note that those are process local indices, so if you’re using MPI with more processes then you might have to change those to global indices, depending on your usecase.