# Dolfinx, ‘Failure to conjugate test function in complex Form’ in Hyperelasticity code

Hello,
I am new to FEniCSx and now studying with demo code.
Among those, I met some error by just following the code on the website :
https://jsdokken.com/dolfinx-tutorial/chapter2/hyperelasticity.html

``````import numpy as np
import ufl

from petsc4py import PETSc
from mpi4py import MPI
from dolfinx import fem, mesh, plot
L = 20.0
domain = mesh.create_box(MPI.COMM_WORLD,[[0.0,0.0,0.0], [L, 1, 1]], [20, 5, 5], mesh.CellType.hexahedron)
V = fem.VectorFunctionSpace(domain, ("CG", 2))
#We create two python functions for determining the facets to apply boundary conditions to

def left(x):
return np.isclose(x, 0)

def right(x):
return np.isclose(x, L)

fdim = domain.topology.dim -1
left_facets = mesh.locate_entities_boundary(domain, fdim, left)
right_facets = mesh.locate_entities_boundary(domain, fdim, right)
#Next, we create a marker based on these two functions

# Concatenate and sort the arrays based on facet indices. Left facets marked with 1, right facets with two
marked_facets = np.hstack([left_facets, right_facets])
marked_values = np.hstack([np.full_like(left_facets, 1), np.full_like(right_facets, 2)])
sorted_facets = np.argsort(marked_facets)
facet_tag = mesh.meshtags(domain, fdim, marked_facets[sorted_facets], marked_values[sorted_facets])
#We then create a function for supplying the boundary condition on the left side, which is fixed.

u_bc = np.array((0,) * domain.geometry.dim, dtype=PETSc.ScalarType)
#To apply the boundary condition, we identity the dofs located on the facets marked by the MeshTag.

left_dofs = fem.locate_dofs_topological(V, facet_tag.dim, facet_tag.find(1))
bcs = [fem.dirichletbc(u_bc, left_dofs, V)]
#Next, we define the body force on the reference configuration (B), and nominal (first Piola-Kirchhoff) traction (T).

B = fem.Constant(domain, PETSc.ScalarType((0, 0, 0)))
T = fem.Constant(domain, PETSc.ScalarType((0, 0, 0)))
#Define the test and solution functions on the space

v = ufl.TestFunction(V)
u = fem.Function(V)
#Define kinematic quantities used in the problem

# Spatial dimension
d = len(u)

# Identity tensor
I = ufl.variable(ufl.Identity(d))

# Right Cauchy-Green tensor
C = ufl.variable(F.T * F)

# Invariants of deformation tensors
Ic = ufl.variable(ufl.tr(C))
J  = ufl.variable(ufl.det(F))
#Define the elasticity model via a stored strain energy density function
#, and create the expression for the first Piola-Kirchhoff stress:

# Elasticity parameters
E = PETSc.ScalarType(1.0e4)
nu = PETSc.ScalarType(0.3)
mu = fem.Constant(domain, E/(2*(1 + nu)))
lmbda = fem.Constant(domain, E*nu/((1 + nu)*(1 - 2*nu)))
# Stored strain energy density (compressible neo-Hookean model)
psi = (mu / 2) * (Ic - 3) - mu * ufl.ln(J) + (lmbda / 2) * (ufl.ln(J))**2
# Stress
# Hyper-elasticity
P = ufl.diff(psi, F)
#Comparison to linear elasticity

#To illustrate the difference between linear and hyperelasticity, the following lines can be uncommented to solve the linear elasticity problem.

# P = 2.0 * mu * ufl.sym(ufl.grad(u)) + lmbda * ufl.tr(ufl.sym(ufl.grad(u))) * I
#Define the variational form with traction integral over all facets with value 2. We set the quadrature degree for the integrals to 4.

# Define form F (we want to find u such that F(u) = 0)
F = ufl.inner(ufl.grad(v), P)*dx - ufl.inner(v, B)*dx - ufl.inner(v, T)*ds(2)
#As the varitional form is non-linear and written on residual form, we use the non-linear problem class from DOLFINx to set up required structures to use a Newton solver.

problem = fem.petsc.NonlinearProblem(F, u, bcs)
``````

The error is as follows:

``````ArityMismatch                             Traceback (most recent call last)
Cell In, line 87
84 F = ufl.inner(ufl.grad(v), P)*dx - ufl.inner(v, B)*dx - ufl.inner(v, T)*ds(2)
85 #As the varitional form is non-linear and written on residual form, we use the non-linear problem class from DOLFINx to set up required structures to use a Newton solver.
---> 87 problem = fem.petsc.NonlinearProblem(F, u, bcs)

File /usr/local/lib/python3.10/dist-packages/ufl/algorithms/check_arities.py:177, in check_form_arity(form, arguments, complex_mode)
175 def check_form_arity(form, arguments, complex_mode=False):
176     for itg in form.integrals():
--> 177         check_integrand_arity(itg.integrand(), arguments, complex_mode)
...
...
File /usr/local/lib/python3.10/dist-packages/ufl/algorithms/check_arities.py:170, in check_integrand_arity(expr, arguments, complex_mode)
168 for arg, conj in arg_tuples:
169     if arg.number() == 0 and not conj:
--> 170         raise ArityMismatch("Failure to conjugate test function in complex Form")
171     elif arg.number() > 0 and conj:
172         raise ArityMismatch("Argument {0} is spuriously conjugated in complex Form".format(arg))

ArityMismatch: Failure to conjugate test function in complex Form

``````

to fix the same issue when I was studying with Heat equation demo code, by using ufl.inner. However, it seems that there is no such point in hyperelasticity code, so I am stuck.
Could you help me to solve this issue, and how to cope with such type of error in general?

Thanks a lot.

I would suggest reading: The Poisson problem with complex numbers — FEniCSx tutorial
As the ordering of test-functions in inner products, and derivatives (used to compute Jacobians) in non-linear systems, need special care when you want to run with complex numbers. Is there a specific reason for wanting to run the Hyper-elasticity code with complex numbers?

Yes, as you noticed, there would be no need for running with complex numbers.
I found that for my environment

``````print(PETSc.ScalarType)
``````

yielding the result,

``````<class 'numpy.complex128'>
``````

SO, I replaced the whole PETSc.ScalarType to ‘numpy.float32’,

``````L = 20.0
domain = mesh.create_box(MPI.COMM_WORLD,[[0.0,0.0,0.0], [L, 1, 1]], [20, 5, 5], mesh.CellType.hexahedron)
V = fem.VectorFunctionSpace(domain, ("CG", 2))

def left(x):
return np.isclose(x, 0)

def right(x):
return np.isclose(x, L)

fdim = domain.topology.dim -1
left_facets = mesh.locate_entities_boundary(domain, fdim, left)
right_facets = mesh.locate_entities_boundary(domain, fdim, right)

marked_facets = np.hstack([left_facets, right_facets])
marked_values = np.hstack([np.full_like(left_facets, 1), np.full_like(right_facets, 2)])
sorted_facets = np.argsort(marked_facets)
facet_tag = mesh.meshtags(domain, fdim, marked_facets[sorted_facets], marked_values[sorted_facets])

u_bc = np.array((0,) * domain.geometry.dim, dtype=np.float32)

left_dofs = fem.locate_dofs_topological(V, facet_tag.dim, facet_tag.find(1))
bcs = [fem.dirichletbc(u_bc, left_dofs, V)]

B = fem.Constant(domain, np.float32((0, 0, 0)))
T = fem.Constant(domain, np.float32((0, 0, 0)))

v = ufl.TestFunction(V)
u = fem.Function(V)

d = len(u)

I = ufl.variable(ufl.Identity(d))

C = ufl.variable(F.T * F)

Ic = ufl.variable(ufl.tr(C))
J  = ufl.variable(ufl.det(F))

E = np.float32(1.0e4)
nu = np.float32(0.3)
mu = fem.Constant(domain, E/(2*(1 + nu)))
lmbda = fem.Constant(domain, E*nu/((1 + nu)*(1 - 2*nu)))

psi = (mu / 2) * (Ic - 3) - mu * ufl.ln(J) + (lmbda / 2) * (ufl.ln(J))**2

P = ufl.diff(psi, F)

F = ufl.inner(ufl.grad(v), P)*dx - ufl.inner(v, B)*dx - ufl.inner(v, T)*ds(2)
problem = fem.petsc.NonlinearProblem(F, u, bcs)
``````

BUT still yielding the same error.
Do I miss something still?

Thanks for spending your valuable time for discussion.

You have installed DOLFINx with complex numbers (i.e. PETSc with complex numbers). PETSc is what DOLFINx interfaces to for linear algebra. You would need to have an installation of PETSc with real numbers (and DOLFINx built on top of that) to be able to solve any problem.

How did you install DOLFINx?

I used Spack to install dofinx, directly following the instruction in GitHub - FEniCS/dolfinx: Next generation FEniCS problem solving environment.
While dealing with this, I found a similar article,
Complex assembly
However, it seems that he used source installation, which did not work in my computer.
I tried to use

``````source /usr/local/bin/dolfinx-complex-mode
source /usr/local/bin/dolfinx-real-mode
``````

but couldn’t find the path, thus failed.
Could you give some instructions for dealing with this situations?

The petsc installation should not be complex by default, as one can see that it is set to false in the petsc spack package. spack/package.py at develop · spack/spack · GitHub

@IgorBaratta can you reproduce this?

I additionally found that my dolfinx is installed as ‘DOLFINX complex’.

Then, do I have to remove and re-install dolfinx?
As I remember I did not choose complex option, however, I might committed some mistakes in installation.

What happens if you use the `Python 3 (ipykernel)`? Does that have DOLFINx installed?

Dokken, You save my life.
Using ipykernel makes everything fine.

I don’t know how to thank you!!

However, how

(python3 ipykernel)

can use dolfinx, and why there are two kernels installed?
Also, would there be no difference between the two, except for complex and real?
As I am not a CS man, I didn’t understand the mechanism.
Anyway, my problem is solved, so I want to say thank you once again.

To me, this does not look like you have used `spack` to install DOLFINx, but that you are running a docker container (possibly `dolfinx/lab` or `ghcr.io/fenics/dolfinx/lab:nightly` as these would have these two python kernels installed.

The difference is how PETSc has been configured. PETSc is the linear algebra backend, which is used to represent matrices and solve the linear and non-linear problems with, once assembled in DOLFINx.
PETSc does not allow for both float types (real and complex) in the same installation, and thus DOLFINx/petsc4py has to be configured with the correct one for the application in hand

In that post, when computing the `L2_error`, shouldn’t it be `ufl.inner` instead of `ufl.dot` ?