Adaptive mesh in time-dependent tasks

Hi everyone!

I’m new in FEniCS and have some troubles in using of the AdaptiveNonlinearVariationalSolver. From some point my solution stop to update with new values. Here is my code

######################################################################################
# Script for calculation of particle shrinkage because of boundary energy minimization
# Differential equation has look 
# \dot{\phi} = - M \left( \alpha\Delta\phi - \omega \phi (1 - \phi)  (1 - 2\phi) \right)
######################################################################################

from fenics import *
import numpy as np
import time

## Input parameters for calculation

start = time.time()

L = 100 # Domain size in x and y, microns
alpha = 2 # Gradient energy pre-factor
omega = 1 # Double-well potential heigth
M_phi = 6 # Phase-field mobility
r_0 = 50 # initial particle radius
n_delta = 1 # numbers of elements per width of diffusive zone
d_t = 1 # time increment
t_end = 100 # time of ending of solution

# For mesh building and initiation of phi_0 field it is crucial to calculate boundary width, i.e
delta = 4.0*np.sqrt(alpha/omega)
print('Diffusive zone width:',delta)

# uniform mesh for this task
N_el = np.int64(L * n_delta / delta) # number of elements in calculation
print('Number of elements in 1D:',N_el)
mesh = RectangleMesh(Point(0,0),Point(L,L),N_el,N_el) # Generation of mesh
V = FunctionSpace(mesh,'P',1) #Function space on created mesh

# Boundary conditions
tol = 1E-10 #tolerance value for boundary searching
def boundary_Dir(x, on_boundary):
    return (near(x[0],L,tol) and on_boundary) or (near(x[1],L,tol) and on_boundary) # Boundary of the Dirichlet condition
# Boundary on Neumann conditions should be not added because they have no influence of energy
bc = DirichletBC(V,Constant(1), boundary_Dir) #initialization of BC
# Initial conditions
phi_0 = Expression('0.5*(1 + tanh(0.0625*delta*(sqrt(x[0]*x[0]+x[1]*x[1])-r_0)))',
                    degree=1,r_0=r_0,delta=delta) # particle on the edge of corner in the mesh
phi_0_int = interpolate(phi_0,V) #interpolation of existing function on the mesh

# Defining of functions
phi = Function(V)
var_phi = TestFunction(V)

# Defining of variational task
phi_i = phi_0_int # phi_i - order parameter at previous test 
F_int = M_phi*alpha*d_t*dot(grad(phi),grad(var_phi)) # Fraction of boundary energy
F_dw = M_phi*omega*d_t*phi*(1-phi)*(1-2*phi) # Fraction of double-well potential
F_td = phi - phi_i  #Fraction of time discretization
F = ( F_int + (F_dw + F_td)*var_phi) * dx # Full variational formulation

# Defining a Jacobian of non-linear form
dphi = TrialFunction(V)
J = derivative(F,phi,dphi)

# Defining of goal function and tolerance value for mesh adaption
M = phi * dx 
tol_opt = 1e-4 * L**2

t = 0

vtkfile = File('particle-shrinkage_adaptive/solution.pvd') # initialization of file variable
vtkfile << (phi_i,t)

t = d_t

# Calculation itself
while t <= t_end :
    
    #Solution of equation
    problem = NonlinearVariationalProblem(F, phi, bc,J)
    solver = AdaptiveNonlinearVariationalSolver(problem, M)
    solver.solve(tol_opt)

    print(phi_i.vector() == phi.leaf_node().vector())
    if (phi_i.vector() == phi.leaf_node().vector()).all():
        phi_i.assign(phi.root_node())
    else:
        phi_i.assign(phi.leaf_node()) #Update data from current step
    print('time=',t)
    vtkfile << (phi_i,t) # writing solution in file
    t += d_t

end = time.time()
print('Time for calculation = ',(end-start)/60,' min')

Have you found the solution to your problem?
Because with a Nonlinear solver for time dependent, i find an accurate solution but with the adaptive solver, the solution is not updated after some timesteps.