ArityMismatch: Adding expressions with non-matching form arguments ('v_0',) vs ()

General
1

Hi,

I am new to fenics and I was trying to put advection to Cahn-Hilliard equation as I wanted to couple it to hydrodynamics. For the time being, I am using constant velocity just to see if the code works. But it gives me error. Here is the code.

import random
from dolfin import *

Class representing the intial conditions

class InitialConditions(UserExpression):
def init(self, **kwargs):
random.seed(2 + MPI.rank(MPI.comm_world))
super().init(*kwargs)
def eval(self, values, x):
values[0] = 0.63 + 0.02(0.5 - random.random())
values[1] = 0.0
def value_shape(self):
return (2,)

Class for interfacing with the Newton solver

class CahnHilliardEquation(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)

Sub domain for Periodic boundary condition (periodic in x)

class PeriodicBoundary(SubDomain):

# Left boundary is "target domain" G
def inside(self, x, on_boundary):
    return bool(x[0] < DOLFIN_EPS and x[0] > -DOLFIN_EPS and on_boundary)

# Map right boundary (H) to left boundary (G)
def map(self, x, y):
    y[0] = x[0] - 1.0
    y[1] = x[1]

Create periodic boundary condition

pbc = PeriodicBoundary()

Model parameters

lmbda = 1.0e-02 # surface parameter
dt = 5.0e-06 # time step
theta = 0.5 # time stepping family, e.g. theta=1 → backward Euler, theta=0.5 → Crank-Nicolson

Form compiler options

parameters[“form_compiler”][“optimize”] = True
parameters[“form_compiler”][“cpp_optimize”] = True

Create mesh and build function space

mesh = UnitSquareMesh.create(100, 100, CellType.Type.quadrilateral)
P1 = FiniteElement(“Lagrange”, mesh.ufl_cell(), 1)
ME = FunctionSpace(mesh, P1*P1,constrained_domain=pbc)

Define trial and test functions

du = TrialFunction(ME)
q, v = TestFunctions(ME)

Define functions

u = Function(ME) # current solution
u0 = Function(ME) # solution from previous converged step

Split mixed functions

dc, dmu = split(du)
c, mu = split(u)
c0, mu0 = split(u0)

Create intial conditions and interpolate

u_init = InitialConditions(degree=1)
#u_init=Expression((‘tanh(1000*(x[1]-0.5)+1.0)0.01’,‘0’), degree=2)
#u_init=Expression(('tanh(1000(x[1]-0.5)+1.0)’,‘0’), degree=2)
u.interpolate(u_init)
u0.interpolate(u_init)

Compute the chemical potential df/dc

c = variable(c)
f = 100c**2(1-c)**2
dfdc = diff(f, c)

mu_(n+theta)

mu_mid = (1.0-theta)mu0 + thetamu

#add advection (advection is causing problems)
vel=Constant((0.0,0.1))
vel_=interpolate(vel,ME)

Weak statement of the equations

L0 = cqdx - c0qdx + dtdot(grad(mu_mid), grad(q))dx + dtdot(grad(c),vel_)dx
L1 = muvdx - dfdcvdx - lmbda*dot(grad(c), grad(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 = CahnHilliardEquation(a, L)
solver = NewtonSolver()
solver.parameters[“linear_solver”] = “lu”
solver.parameters[“convergence_criterion”] = “incremental”
solver.parameters[“relative_tolerance”] = 1e-6

Output file

file = File(“output.pvd”, “compressed”)

Step in time

t = 0.0
T = 50*dt
while (t < T):
t += dt
u0.vector()[:] = u.vector()
solver.solve(problem, u.vector())
file << (u.split()[0], t)

2

Please use 3x` formatting to make your code properly formatted, ie.

```python
#add code here
```
3

I fixed the issue. There was a typo. Thank anyway. I will take care of your suggestion next time.