Mesh adaptivity in fenics

Hello guys,
I would like to control the meshsize [h] based of value of the solution [u] i.e f(u) be the function/ criterion evaluated at each node. Then in 2d case with linear elements I will have 3 values for each element f_i(u)

One example for f(u) would be f(u):= u ||grad(u)|| / 25. Solution u is scalar quantity.

for each cell :
if h > max(f_i(u)) :
Perform refinement

This is used to solve nonliner solid mechanics problem with load step.

This is what I would like achieve this. But i don’t think demo example Auto adaptive Poisson equation won’t be helpful since it requires global criterion. However if I were to refine the mesh manually using refine function then I have to manually updated the spaces and also transfer the solution to the new mesh.

Is there a better approach to achieve this local adaptivity using some inbuilt functions of fenics?

1 Like

Hello guys,
As with the second approach to manually adapting the mesh. The code would look something like:

from dolfin import *
import numpy as np

# Sub domain for Dirichlet boundary condition
class DirichletBoundary(SubDomain):
    def inside(self, x, on_boundary):
        return abs(x[0] - 1.0) < DOLFIN_EPS and on_boundary

# Create mesh and define function space
mesh = UnitSquareMesh(32, 32)
meshfile = File("mesh.pvd")
meshfile << mesh

V = FunctionSpace(mesh, "CG", 1)

# Define boundary condition
g = Expression("t",t=1.0,degree=1)
bc = DirichletBC(V, g, DirichletBoundary())

# Define variational problem
u = Function(V)
v = TestFunction(V)
f = Expression("x[0]*sin(x[1])", degree=2)
F = inner((1 + u**2)*grad(u), grad(v))*dx - f*v*dx


#cell local refinement criterion
def meshAdaptivityMarker(mesh, uTemp):
	cellMarkers = MeshFunction("bool",mesh,2)
	cellMarkers.set_all(False)
	for cell in cells(mesh):
		coordinates = np.reshape(cell.get_coordinate_dofs(),(-1,mesh.geometry().dim()))
		cellRegPar = np.array([uTemp(x) for x in coordinates])
		minRegPara= np.nanmin(cellRegPar)
		if(cell.h() > minRegPara/147.0):
			cellMarkers[cell]=True
	return cellMarkers


# Save solution in VTK format
file = File("nonlinear_poisson.pvd")

# This for loops through various load step
for i in [1.0,2.0,3.0]:
	g.t = i
# Compute solution
	solve(F == 0, u, bc,
      	solver_parameters={"newton_solver":{"relative_tolerance":1e-6}})
	file << u
	cellMarkers = meshAdaptivityMarker(mesh, u)
	mesh = refine(mesh,cellMarkers)
	V = FunctionSpace(mesh, "CG", 1)
	u = Function(V)
	v = TestFunction(V)
	bc = DirichletBC(V, g, DirichletBoundary())
	F = inner((1 + u**2)*grad(u), grad(v))*dx - f*v*dx
	meshfile << mesh

However this approach is not optimal. So I would be very grateful for any advise to improve this or suggest any other approach.

1 Like

Hi Djay92,

Have you found a solution for your problem?
I am wondering if I could use such a method to deal with large mesh deformation (?)

Thank you in advance,