Hi,

I am trying to solve a 2D contact problem based on the 3D code available here:

But I can’t figure out why my code is not running. My code is the following

```
Ny=30 # Refinamento da malha na direção y
Nx=30 # Refinamento da malha na direção x
mesh = fn.RectangleMesh(fn.Point(0, -1), fn.Point(1, 0),Nx,Ny) # gera a malha
mesh.coordinates()[:,0]=mesh.coordinates()[:,0]**2
mesh.coordinates()[:,1]=-mesh.coordinates()[:,1]**2
fn.plot(mesh, title="Malha")
class Top(fn.SubDomain):
def inside(self, x, on_boundary):
return fn.near(x[1], 0.) and on_boundary
def symmetry_x(x, on_boundary):
return fn.near(x[0], 0.) and on_boundary
def bottom(x, on_boundary):
return fn.near(x[1], -1.) and on_boundary
# exterior facets MeshFunction
facets = fn.MeshFunction("size_t", mesh, 1)
facets.set_all(0)
Top().mark(facets, 1)
ds = fn.Measure('ds', subdomain_data=facets)
R = 0.5
d = 0.02
obstacle = fn.Expression("-d+pow(x[0],2)/2/R", d=d, R=R, degree=2)
V = fn.VectorFunctionSpace(mesh, "CG", 1)
V2 = fn.FunctionSpace(mesh, "CG", 1)
V0 = fn.FunctionSpace(mesh, "DG", 0)
u = fn.Function(V, name="Displacement")
du = fn.TrialFunction(V)
u_ = fn.TestFunction(V)
bc =[fn.DirichletBC(V, fn.Constant((0., 0.)), bottom),
fn.DirichletBC(V.sub(0), fn.Constant(0.), symmetry_x)]
E = fn.Constant(10.)
nu = fn.Constant(0.3)
mu = E/2/(1+nu)
lmbda = E*nu/(1+nu)/(1-2*nu)
def eps(v):
return fn.sym(fn.grad(v))
def sigma(v):
return lmbda*fn.tr(eps(v))*fn.Identity(2) + 2.0*mu*eps(v)
def ppos(x):
return (x+np.abs(x))/2.
pen = fn.Constant(1e4)
form = fn.inner(sigma(u), eps(u_))*fn.dx + pen*fn.dot(u_[1], ppos(u[1]-obstacle))*fn.ds(1)
J = fn.derivative(form, u, du)
problem = fn.NonlinearVariationalProblem(form, u, bc, J=J)
solver = fn.NonlinearVariationalSolver(problem)
solver.parameters["newton_solver"]["linear_solver"] = "cg"
solver.parameters["newton_solver"]["preconditioner"] = "ilu"
solver.solve()
import matplotlib.pyplot as plt
p = fn.plot(sigma(u)[1,1], mode='color')
plt.colorbar(p)
plt.title(r"$\sigma_{yy}$",fontsize=26)
plt.show()
```

I tried different nonlinear solvers but the result is the same. What is wrong?