Mixed function space with dolfinx - how to implement an integral constraint

[lclaire]

Hello,

I would like to define a variational problem that uses a mixed function space on a 2D domain in dolfinx. Here is how I define the 2D domain:

L = 1
W = 0.1
domain = mesh.create_rectangle(MPI.COMM_WORLD, [np.array([0, -W/2]), np.array([L, W/2])], [2, 2], cell_type=mesh.CellType.triangle)

Next, I would like the mixed function space to handle a 3D displacement field, with

Vy = ufl.VectorElement('Lagrange', domain.ufl_cell(), degree=1, dim = 3)

and 3 Lagrange multipliers, uniform over the domain, like so

Vp = ufl.VectorElement('Real', domain.ufl_cell(), dim=3)

However, it seems that I cannot assemble these two components:

V = dolfinx.fem.FunctionSpace(domain, ufl.MixedElement([Vy, Vp]))

yields an error:

File ~/anaconda3/envs/fenicsx-0.5.1/lib/python3.9/site-packages/basix/ufl_wrapper.py:1269, in convert_ufl_element(element)
   1266 if family_name == "DPC":
   1267     discontinuous = True
-> 1269 family_type = _basix.finite_element.string_to_family(family_name, element.cell().cellname())
   1270 cell_type = _basix.cell.string_to_type(element.cell().cellname())
   1272 variant_info = {
   1273     "lagrange_variant": _basix.LagrangeVariant.unset,
   1274     "dpc_variant": _basix.DPCVariant.unset
   1275 }

File ~/anaconda3/envs/fenicsx-0.5.1/lib/python3.9/site-packages/basix/finite_element.py:87, in string_to_family(family, cell)
     84 if family in families:
     85     return families[family]
---> 87 raise ValueError(f"Unknown element family: {family} with cell type {cell}")

ValueError: Unknown element family: Real with cell type triangle

Many thanks in advance for your help!

[nate]

I believe this is because elements of type Real are not supported by dolfinx. There’s a great deal of discussion of the best way to do this. The main problem is that populating an entire row of a matrix kills performance.

Perhaps @mscroggs or @dokken know more.

Implement global support dofs

Missing Real finite element

[lclaire]

thank you for your answer. I will try to find an alternative formulation for my problem that does not involve elements of type Real.

My first intention was to use the element of type Real to enforce an integral constraint over my domain via a Lagrange multiplier, for example, suppose I want to impose the average displacement u over the 2D domain \Omega to be zero \int_\Omega u(x) \mathrm{d}x=0. As this is a scalar constraint, it would be sufficient for the corresponding Lagrange multiplier to be defined over 1 DOF only, I am not sure how to implement that within dolfinx though.

[bagla0]

I am also struggling to implement Real element in dolfinx function space for implementing lagrange multiplier similar to Periodic homogenization of linear elastic materials — Numerical tours of continuum mechanics using FEniCS master documentation.

Any example or mwe where real element implemented in fenicsx would be highly helpful.

I have seen posts of work still going to update dolfinx with real element. Let me know if any alternative way is available based on above example.

[francesco-ballarin]

See recent posts on this topic:

• How can I use Block Forms, Multiphenicsx, and Lagrange Multipliers to Solve Poisson's Equation with Pure Neumann Boundary Conditions? - #4 by francesco-ballarin (go to the tutorial listed in the accepted solution, and notice that the case without domain decomposition does not actually require multiphenicsx)
• Singular poisson equation in DOLFINx (click on the link in the second post)

[bagla0]

In domain without decomposition, I am getting difficulty to relate the lagrange multiplier.
In reference tutorial multiphenicsx/tutorials/08_singular_poisson/tutorial_create_nullspace.ipynb at main · multiphenics/multiphenicsx (github.com), only one test function \delta v is used.

However, in implementing Lagrange multiplier \lambda

@bagla0

Periodic homogenization of linear elastic materials — Numerical tours of continuum mechanics using FEniCS master documentation.

we have two test functions \delta \lambda for real space and \delta v for vector space. Can you tell me with analytical expressions how I can relate the tutorial with implementing constraints using lagrange multiplier deined in peroidic homo tutorial?

[francesco-ballarin]

By following the tutorials above, you wouldn’t need the lagrange multiplier at all, so no \lambda or \delta\lambda at all. Rather, you would must inform the solver that the displacement has the following null space (sentence copied from comet-fenics):

The above problem is not well-posed because of the existence of rigid body translations

[bagla0]

Hi @dokken and @francesco-ballarin
Thank you for your previous support.
I was understanding the domain without decomposition tutorial (neumann BC using nullspace) and relating with mine in Constrains using null space in dolfinx - dolfinx - FEniCS Project. I have general queries related to its application:

1. Does bringing nullspace vector means we are identifying the rigid body motion related to coefficient matrix A (obtained from weak form)?

2. For the application of constraining the linear system of equation (Ax=b) other than rigid body motion constraints, does nullspace still works in alternative to lagrange multiplier method for dolfinx?

3. Can you provide/comment on the nullspace vector for Periodic homogenization of linear elastic materials — Numerical tours of continuum mechanics using FEniCS master documentation
   where,

the fluctuating field was averaged to zero over entire domain \int v dx =0 (2 constraints for 2 dimensions) and third constraint ( \int x_1 v_2- x_2 v_1 dx =0) where v_1 and v_2 are two components of fluctuating field with x_1, x_2 as spatial coordinates.

For this case, I tried the nullspace vectors previously as below,

@bagla0

# Refering neumann BC example without domain decomposition
# Assemble the discrete system

# Nullspace Implement
c1 = dolfinx.fem.Function(V)
c1.interpolate(lambda x: np.array([[1],[0]])) # Nullspace related to Constraint $v1$ averaged to zero.
C1 = c1.vector
#
c2 = dolfinx.fem.Function(V)
c2.interpolate(lambda x: np.array([[0],[1]])) # Nullspace related to Constraint $v[1]$ averaged to zero.
C2 = c2.vector
#
c3 = dolfinx.fem.Function(V)
c3.interpolate(lambda x: np.array((-x[1],x[0])) # Nullspace related to Constraint $-x2*v[0]+v[1]*x1$ averaged to zero.
C3 = c3.vector

# Orthonormalize vectors
dolfinx.la.orthonormalize([C1, C2, C3])
nullspace = petsc4py.PETSc.NullSpace().create(vectors=[C1, C2, C3], comm=mesh.comm)

# Set the nullspace
A.setNullSpace(nullspace)

I made the above attempt previously for nullspace vectors, which you said wrong. Can you here provide the correct nullspace vectors through telling how we look for it (MWE)?

4. What is the process of identifying the nullspace vectors (line-2)?
   
   
   
   Suppose, if we would be having pure Dirichlet BC in poisson equation, then, what would be the nullspace vector in 17. Singular Poisson — FEniCS Project?

I am sorry if I am unable to put right questions. Kindly guide me with whether nullspace works here. I understand some questions might be repeated, sorry for that.

Thanking You

[francesco-ballarin]

Does bringing nullspace vector means we are identifying the rigid body motion related to coefficient matrix A (obtained from weak form)?

Yes, the null space you set up should be associated to the rigid body motions.

For the application of constraining the linear system of equation (Ax=b) other than rigid body motion constraints, does nullspace still works in alternative to lagrange multiplier method for dolfinx?

We cannot claim that it will work in any case until we know what your actual problem is.

Can you provide/comment on the nullspace vector for Periodic homogenization of linear elastic materials — Numerical tours of continuum mechanics using FEniCS master documentation

I fear that we are running in circles. I already quoted

The above problem is not well-posed because of the existence of rigid body translations

from that same tutorial, so let’s see if the author @bleyerj has any further suggestion.

Suppose, if we would be having pure Dirichlet BC in poisson equation, then, what would be the nullspace vector in 17. Singular Poisson — FEniCS Project?

That the case in my tutorial: it would be the constant vector, as given a solution u also the sum between u and any constant would still be a solution.

[bagla0]

Hello @francesco-ballarin,

I am stuck in previous problem of non-availability of real function space in dolfinx. My problem constraints are over constraining the non-singular coefficient matrix. For which, I know that nullspace works only for singular matrices.

I tried with multiphenicsx here but it clones Lagrange space for applying Lagrange multiplier method.

1. This gives twice the dof (Lagrange original space+cloned space).
2. It requires boundary condition to be applied only on boundary, however, my constraints require \int _V u dx over entire volume. I don’t find a way to apply constraints over volume.
   Can you suggest any way to manually bring real function space in dolfinx?

Ques-Is there any existing method like Penalty method for applying general constraints over volume in dolfinx? I might have less knowledge to provide the above observation, if it’s wrong then, kindly correct me. Any help is greatly appreciated.

[astqx]

Hey, I have a similar problem. Is there an update on this?

[bagla0]

@astqx Try using nullspace vector to constrain the rigid body motion. Use above quoted discussion links and try this.