Hello,
I am trying to solve a nonlinear complex-value problem with item of 1j*|u|^2*u. When I create the weak form in axi-symmetric for this item as:
(1j)*(ufl.inner(u,u))*u*ufl.conj(v)*xm[0]*ufl.dx
which bring error of “ArityMismatch: Adding expressions with non-matching form arguments (‘v_1’,) vs (‘conj(v_1)’,)”
Does any one know how to solve this issue.
Below is the full code:
# %%
import dolfinx
from mpi4py import MPI
import numpy as np
from dolfinx.fem import functionspace, Function,Constant
from dolfinx.mesh import create_unit_cube,create_unit_square,create_rectangle
import matplotlib.pyplot as plt
from ufl import TrialFunction, TestFunction, inner,grad
import ufl
from petsc4py import PETSc
import dolfinx.fem.petsc
from dolfinx.nls.petsc import NewtonSolver
comm = MPI.COMM_WORLD
comm_rank = comm.Get_rank()
numProc=comm.Get_size()
# %%
nEle=20
#msh=create_unit_cube(comm,nEle,nEle,nEle)
msh = create_rectangle(comm,points=((0.0, 0.0), (1.0, 2.0)), n=(16, 32))
xm = ufl.SpatialCoordinate(msh)
V = functionspace(msh, ("Lagrange", 1))
#A = Function(V, dtype=np.complex128)
u_c = Function(V, dtype=np.complex128)
u_c.interpolate(lambda x:(x[1]+(x[0]**2 + x[1])*1j))
#u = TrialFunction(V)
u=Function(V, dtype=np.complex128)
v = TestFunction(V)
f = Function(V, dtype=np.complex128)
f.interpolate(lambda x:(1j+5-(x[1]**2+(x[0]**2+x[1])**2)*(x[1]*1j-x[0]**2-x[1])))
k=Constant(msh,PETSc.ScalarType(((1,0),(0,0))))
# %%
a =(u.dx(1)*ufl.conj(v))*xm[0]*ufl.dx+(1j)*ufl.inner(u.dx(0),v.dx(0)) *xm[0]* ufl.dx+(1j)*(ufl.inner(u,u))*u*ufl.conj(v)*xm[0]*ufl.dx
L = ufl.inner(f, v)*xm[0] * ufl.dx
F=a-L
# %%
msh.topology.create_connectivity(msh.topology.dim-1, msh.topology.dim)
boundary_facets = dolfinx.mesh.exterior_facet_indices(msh.topology)
boundary_dofs = dolfinx.fem.locate_dofs_topological(V, msh.topology.dim-1, boundary_facets)
bc = dolfinx.fem.dirichletbc(u_c, boundary_dofs)
problem = dolfinx.fem.petsc.NonlinearProblem(F, u, bcs=[bc])
solver = NewtonSolver(comm, problem)
# Set Newton solver options
solver.atol = 1e-8
solver.rtol = 1e-8
solver.convergence_criterion = "incremental"
# We can customize the linear solver used inside the NewtonSolver by
# modifying the PETSc options
ksp = solver.krylov_solver
opts = PETSc.Options() # type: ignore
option_prefix = ksp.getOptionsPrefix()
opts[f"{option_prefix}ksp_type"] = "gmres"
opts[f"{option_prefix}pc_type"] = "lu"
ksp.setFromOptions()
r = solver.solve(u)
print("Number of iterations: %d, converged %d"%(r[0],r[1]))
# %%
u.name='numerical'
u_c.name='exact'
file_results= dolfinx.io.XDMFFile(MPI.COMM_WORLD, "./fenicsx_complex/example.xdmf","w" )
file_results.write_mesh(msh)
file_results.write_function(u,0)
file_results.write_function(u_c,0)
file_results.close()
# %%
max_error = (np.max(np.abs(u_c.x.array-u.x.array)))
idx= (np.argmax(np.abs(u_c.x.array-u.x.array)))
print(max_error)
print(u_c.x.array[idx],u.x.array[idx])
fig=plt.figure(1)
ax=plt.subplot(1,1,1)
ax.plot(u_c.x.array.imag[:], 'o', label='exact')
ax.plot(u.x.array.imag[:], '+', label='Fenicsx')
ax.legend()
fig=plt.figure(2)
ax=plt.subplot(1,1,1)
ax.plot(u_c.x.array.real[:], 'o', label='exact')
ax.plot(u.x.array.real[:], '+', label='Fenicsx')
ax.legend()
