Hey Guys,

I hope you are great. I am completely new to Fenics and this is my first time posting to the group. I appreciate that if you can help me.

I am working on some phase field problem. For start, I need to solve a time-dependent Ginzburg-Landau equation with a highly nonlinear source term:

(1/lmda) du/dt = lmbda*betta*grad^2(u) - 2*A*lmbda*(u+2u**3-3u**2)+12*lmbda*delG*(u**2 - u**3)

where all the coefficients are given below in the code. In fact, the RHS is the driving force for propagation of interface. u is in fact a scalar parameter. The initial condition for u is pow(1+exp(x[0] - 0.5), -1).

I have started from Cahn-Hilliard equation, then try the nonlinear time-dependent burgers equation, However I was not able to solve the problem. The latest code that I have written is given below for your consideration:

from fenics import *

# Class for interfacing with the Newton solver

class GinzburgLandau(NonlinearProblem):

def **init**(self, a, L):

NonlinearProblem.**init**(self)

self.L = L

self.a = a

def F(self, b, x):

assemble(self.L, tensor=b)

def J(self, A, x):

assemble(self.a, tensor=A)

# Create mesh and build function space

mesh = UnitSquareMesh(80, 80)

V = FunctionSpace(mesh, ‘P’, 2)

# Define trial and test functions

du = TrialFunction(V)

v = TestFunction(V)

# Define functions and initial condition

u = Function(V) # current solution

u_n = Function(V) # solution from previous converged step

u_init = Expression(‘pow((1+exp((x[0] - 0.5) / 2.95e-9)), -1)’, degree=2)

u.interpolate(u_init)

u_n.interpolate(u_init)

# Weak statement of the equations

lmbda = 2600.0

betta = 2.59e-10

delG = -1073.3e6

A = 7924.0e6

dt = 0.1e-14

k = Constant(dt)

L0 = u*v*dx - u_n*v*dx + lmbda*betta*k*dot(grad(u), grad(v)) dx*A

L1 = (2k

*lmbda*(u + 2

*u**3 - 3*u

**2) + 12**2 - u**3))

*k*lmbda*delG*(u*v*dx

L = L0 + L1

# Compute directional derivative about u in the direction of du (Jacobian)

a = derivative(L, u, du)

# Create nonlinear problem and Newton solver

problem = GinzburgLandau(a, L)

solver = NewtonSolver()

solver.parameters[“linear_solver”] = “lu”

# solver.parameters[“preconditioner”] = “ilu”

solver.parameters[“convergence_criterion”] = “incremental”

solver.parameters[“relative_tolerance”] = 1e-12

# solver.parameters[“absolute_tolerance”] = 1e-13

# solver.parameters[“maximum_iterations”] = 25

# solver.parameters[“relaxation_parameter”] = 1.2

file_u = File(‘ginzburglandau/u.pvd’)

num_steps = 500

t = 0.0

for n in range(num_steps):

```
# Update current time
t += dt
print('t = ', t)
solver.solve(problem, u.vector())
file_u << u
# Update the previous solution
u_n.assign(u)
```

Unfortunately, the interface does not propagate with time and it stops in the middle of the square. I appreciate that if you can show me a start point that I can get to the solution of my problem.

I know you are really busy, but it would be very kind of you if you can respond it.

Looking forward to hearing form you.

Bests,

Benhour