# Hyperelastic problem with two different materials

hi all,

im trying to run an hyperplastic problem with 2 different materials and its always fails after few iterations ,
any idea what can i improve?

im not sure which parts of the code to share

this is the defenition of my mesh

``````E1, E2, nu = 1,10,0.2
mesh = UnitCubeMesh(24, 16, 16)

class Omega_0(SubDomain):
def inside(self, x, on_boundary):
return x[0] <= 1 / 2 + DOLFIN_EPS

class Omega_1(SubDomain):
def inside(self, x, on_boundary):
return x[0] > 1 / 2 - DOLFIN_EPS

subdomain_0 = Omega_0()
subdomain_1 = Omega_1()
materials = MeshFunction('size_t', mesh, mesh.topology().dim())
materials.set_all(0)
subdomain_0.mark(materials, 1)
subdomain_1.mark(materials, 2)

boundary_parts = MeshFunction('size_t', mesh, mesh.topology().dim() - 1)
x0 = AutoSubDomain(lambda x: near(x[0], 0))
x1 = AutoSubDomain(lambda x: near(x[0], 1))
y0 = AutoSubDomain(lambda x: near(x[1], 0))
z0 = AutoSubDomain(lambda x: near(x[2], 0))
x0.mark(boundary_parts, 1)
y0.mark(boundary_parts, 2)
z0.mark(boundary_parts, 3)
x1.mark(boundary_parts, 4)
class K(fe.UserExpression):
def __init__(self, materials, k_0, k_1, **kwargs):
super().__init__(**kwargs)
self.materials = materials
self.k_0 = k_0
self.k_1 = k_1

def eval_cell(self, values, x, cell):
if self.materials[cell.index] == 1:
values[0] = self.k_0
else:
values[0] = self.k_1

def value_shape(self):
return ()

E_global = K(materials, E1, E2, degree=0)

DG = FunctionSpace(mesh, "DG", 2)
materials_function = Function(DG)
E_values = [E1, E2]
materials_function = project(E_global, DG)

MU1, LAMBDA1 = Constant(E1/(2*(1 + nu))), Constant(E1*nu/((1 + nu)*(1 - 2*nu)))
MU2, LAMBDA2 = Constant(E2/(2*(1 + nu))), Constant(E2*nu/((1 + nu)*(1 - 2*nu)))

class MUMU(fe.UserExpression):
def __init__(self, materials, MU_0, MU_1, **kwargs):
super().__init__(**kwargs)
self.materials = materials
self.MU_0 = MU_0
self.MU_1 = MU_1

def eval_cell(self, values, x, cell):
if self.materials[cell.index] == 1:
values[0] = self.MU_0
else:
values[0] = self.MU_1

def value_shape(self):
return ()

MU = MUMU(materials, MU1, MU2, degree=0)

def __init__(self, materials, lam_0, lam_1, **kwargs):
super().__init__(**kwargs)
self.materials = materials
self.lam_0 = lam_0
self.lam_1 = lam_1

def eval_cell(self, values, x, cell):
if self.materials[cell.index] == 1:
values[0] = self.lam_0
else:
values[0] = self.lam_1

def value_shape(self):
return ()

LAMBDA = lamabada(materials, LAMBDA1, LAMBDA2, degree=0)
``````

of course im solving non linear problem

``````# Create nonlinear variational problem and solve
problem = NonlinearVariationalProblem(F, _u_p, bcs=bcs, J=J)
solver = NonlinearVariationalSolver(problem)
solver.parameters['newton_solver']['relative_tolerance'] = 1e-6
solver.solve()
``````

this is the error i get :
Newton iteration 5: r (abs) = 5.658e+22 (tol = 1.000e-10) r (rel) = 1.210e+20 (tol = 1.000e-06)
Traceback (most recent call last):
File “/home/mirialex/PycharmProjects/hyperelastic_model/new_from_internet.py”, line 228, in
solver.solve()
RuntimeError:

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at

*** Remember to include the error message listed below and, if possible,
*** include a minimal running example to reproduce the error.

*** -------------------------------------------------------------------------
*** Error: Unable to solve linear system using PETSc Krylov solver.
*** Reason: Solution failed to converge in 0 iterations (PETSc reason DIVERGED_PC_FAILED, residual norm ||r|| = 0.000000e+00).
*** Where: This error was encountered inside PETScKrylovSolver.cpp.
*** Process: 0

*** DOLFIN version: 2019.2.0.dev0
*** Git changeset: ubuntu
*** -------------------------------------------------------------------------

As your problem is not reproducible, i.e you haven’t supplied the necessary imports or most of the variable definitions, it is hard to Give you any guidance.

What I would suggest as a first approach is to let the two materials have almost the same properties, and see if it converges.
If it converges gradually increase the difference and find the point of failure.