Solver Nolinear

Linear Algebra
1

** Solver problem using the Logaritm of Left Cauchy-Green tensor (spatial configuration)**

hello guys, I’m trying to implement code for large deformations using the Hecky’s tensor ,
I’m trying to use a non-linear solver model, given the conditions of the problem, however, I don’t know how to configure it, the model perform two iterations and finish. Would anyone have any ideas to help me?

Version: Dolfinx 0.6.0 Python
Windows Subsystem for Linux (WSL) in Visual Studio Code

from dolfinx import fem
import ufl
from petsc4py.PETSc import ScalarType
from dolfinx import fem,log, nls ,mesh
from dolfinx.io import gmshio
from mpi4py import MPI
import numpy as np
from petsc4py import PETSc


*Geometry , boundary conditions and Scaled variable:*

carregamento= -50000
E = 78e6
poisson = 0.3
lambda_ = E*poisson / ((1+poisson)*(1-2*poisson))
G= E / 2*(1+poisson)


#Geometry
domain= mesh.create_rectangle(MPI.COMM_WORLD, [np.array([0,0]), np.array([2,1])],
                  [100,100], cell_type=mesh.CellType.triangle)

V = fem.VectorFunctionSpace(domain, ("CG", 1))

#boundary conditions 

def clamped_boundary(x):
    return np.isclose(x[0], 0)

fdim = domain.topology.dim - 1
boundary_facets = mesh.locate_entities_boundary(domain, fdim, clamped_boundary)
u_D = np.array([0,0], dtype=ScalarType)
bc = fem.dirichletbc(u_D, fem.locate_dofs_topological(V, fdim, boundary_facets), V)



#Aproximação das funções teste e funcao incognita
v = ufl.TestFunction(V)
u= fem.Function(V)

#tensor de cauchy-Green left 
d = len(u)
I = (ufl.Identity(d)) 
F =(I + ufl.grad(u)) 
C= (F * F.T )

#jacobian 
Jhc_e= ufl.det(F)
J= Jhc_e

#Hencky's Strain
N, M = C.ufl_shape
T_hencky=ufl.as_tensor([[(1/2)*(ufl.ln(C[i,j])) for i in range(N)] for j in range(M)])


#Tensor de tension
T_tension= (1/J)*(2.0 * G * T_hencky+ lambda_ * ufl.tr(T_hencky) * I )

#Function load
f = fem.Constant(domain, ScalarType((0,carregamento )))

#Formulação variacional para Div T = 0, em Tn= f,u= 0 em x=0
a = ufl.inner(ufl.grad(v),T_tension) * ufl.dx
L = ufl.inner(f, v) * ufl.ds

#Weak form for  Div T = 0,  Tn= f, u= 0 in 
F_bilinear =  + ufl.inner(ufl.grad(v),T_tension) * ufl.dx - ufl.inner(f, v) *ufl.ds

#Formulação variacional para Div T = 0, em Tn= f,u= 0 in facet_tags.find(16)

problem = fem.petsc.NonlinearProblem(F=F_bilinear ,u=u,bcs=[bc])


#solver no linear 
solver = nls.petsc.NewtonSolver(MPI.COMM_WORLD, problem)
solver.convergence_criterion = "incremental"
solver.rtol = 1e-15
solver.report = True


#solver linear
ksp = solver.krylov_solver
opts = PETSc.Options()
option_prefix = ksp.getOptionsPrefix()
opts[f"{option_prefix}ksp_type"] = "cg"
opts[f"{option_prefix}pc_type"] = "gamg"
opts[f"{option_prefix}pc_factor_mat_solver_type"] = "mumps"

ksp.setFromOptions()


#vizualization 
log.set_log_level(log.LogLevel.INFO)
n, converged = solver.solve(u)
assert(converged)
print(f"Number of interations: {n:d}")

Failed: 2023-08-22 20:53:30.274 ( 2.108s) [main ] petsc.cpp:677 INFO| PETSc Krylov solver starting to solve system. 2023-08-22 20:53:30.334 ( 2.168s) [main ] petsc.cpp:677 INFO| PETSc Krylov solver starting to solve system. 2023-08-22 20:53:30.339 ( 2.173s) [main ] NewtonSolver.cpp:36 INFO| Newton iteration 2: r (abs) = 0 (tol = 1e-10) r (rel) = -nan(tol = 1e-15) 2023-08-22 20:53:30.340 ( 2.173s) [main ] NewtonSolver.cpp:255 INFO| Newton solver finished in 2 iterations and 0 linear solver iterations. Number of interations: 2 2023-08-22 20:53:30.500 ( 2.334s) [main ] loguru.cpp:526 INFO| atexit

2

I would start with:

Change the convergence criterion to “residual”,
This gives you slightly more information about the initial steps of the Newton solver

Use a direct solver (ksptype=preonly, pc_type=lu) to make sure that a direct solver can solve your problem.

When using gamg you should attach the near nullspace of the problem, see https://github.com/FEniCS/dolfinx/blob/18e0b139d37ddb764f755fd8267db93d044bc038/python/demo/demo_elasticity.py#L169

3

well, i try change this solver.

solver.convergence_criterion = "residual"
solver.rtol = 1e-6
solver.max_it = 100
solver.report = True

solver.atol = 1e-6
solver.solver_type = "gmres"
solver.preconditioner_type = "lu"
solver.initial_guess = None  
solver.line_search = True
solver.jacobian_update = "approximate"
solver.error_on_nonconvergence = True


but i keep getting this error: 
2023-08-23 13:53:49.124 (   6.927s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2023-08-23 13:53:49.143 (   6.946s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 50: r (abs) = *-nan* (tol = 1e-06) r (rel) = -nan(tol = 1e-06)

this *-nan*  term does not represent a number, so I should consider that there is some mismatch error between the input values ​​perhaps?
Could you explain to me the main causes of this term ?

’ I’m a beginner in the fenics project, I apologize if I ask too obvious questions