Any way to access numbering of the DoFs in the assembled mpc matrix

dolfinx
,
1

Hello everyone,

When assembling matrix using the dolfinx-mpc interface, the slave DoFs will be eliminated. So is there any routine to access the numbering of the remained master DoFs in the distributed matrix? Or is there any rule for the creation of the sparsity pattern, for example, master dofs will remain on their original ranks after elimination of the slave DoFs?

To be more specifical, what I want to do is to add a penalty term to certain diagonal elements to eliminate the rigid motion in an elastic problem. Because Dirichlet bcs cannot be applied to the same DoFs where mpc is applied in dolfinx-mpc, and Nitsche method is rather complicated to implement, I’m trying to directly manipulate the assembled matrix.

2

The matrix retains its original size, so all slave dofs are kept in their original position. There are only additions of ghost degrees of freedom to handle transferring data from slave to master dofs.
The MPC class can give you both slave and master indices, see

Please note that rigid Motions with MPC is handled in

and

3

Sorry, but I didn’t get the point. When applying mpc to the system, we will calculate K' = T^TKT, where K is the system stiffness matrix and T is the transformation matrix, am I right? So why does the matrix retain its original size?

4

I keeps the original size to avoid having to resize the problem matrix and do a whole lot of extra stuff with the sparsity pattern, similar to how Dirichlet conditions are handled in dolfinx. The slave rows are set do being the identity row.

5

Thank you very much, I fully understand now!

6

Just one more question. Can we use A.setNullspace to eliminate the rigid motion of a system and use A.setNearNullSpace to improve the gamg solver at the same time?

7

I would just use nearnullspace with gamg. I’ve not seen any positive effect using both.

8

Thank you, I’ll have a try.

9

Hi, dokken. Sorry to disturb you again. I have another question and I think it is relative to this post.

As discussed before, the rows corresponding to mpc slaves in the system stiffness matrix K will be set to zeros except for the diagonal elements. And the corresponding rhs elements will be set to zeros as well. Thus, mathematically speaking, elements of the nullspace vectors for K should also be zeros at the slave dofs. However, the rigid_motions_nullspace routine in dolfinx_mpc.utils will not set the slave elements to zeros. After a simple change like:

...

la.orthonormalize(nullspace_basis)
assert la.is_orthonormal(nullspace_basis, float(np.finfo(_x.x.array.dtype).eps))

# insert this code snippet
for b in nullspace_basis:
    mpc._cpp_object.homogenize(b.array)
    b.scatter_forward()

local_size = V.dofmap.index_map.size_local * V.dofmap.index_map_bs
basis_petsc = [
        PETSc.Vec().createWithArray(x[:local_size], bsize=gdim, comm=V.mesh.comm)  
        for x in basis
]

...

The nullspace vectors will be successfully set to zeros at the slave dofs. So why does the original code still works even though it does not homogenize the nullspace vectors? I’m a little confused.