ArityMismatch error when applying nonlinear operator "Power"

Hi all!

I am encountering an ArityMismatch error: “ArityMismatch: Applying nonlinear operator Power to expression depending on form argument v_1.”

Im trying to solve the ginzburg landau equations, written as two coupled expressions for the real and imaginary part of the order parameter. Here is my code:

from fenics import *

# Create mesh
mesh = UnitSquareMesh(32, 32)

# Define function spaces
P1 = FiniteElement("P", triangle, 1)
element = MixedElement([P1, P1])
V = FunctionSpace(mesh, element)

# Define trial and test functions
rho, theta = TrialFunctions(V)
v, w = TestFunctions(V)

# Define constants
alpha = Constant(1.0)
beta = Constant(1.0)
gamma = Constant(0.1)

# Define boundary conditions
def left_boundary(x, on_boundary):
    return on_boundary and near(x[0], 0.0)

def right_boundary(x, on_boundary):
    return on_boundary and near(x[0], 1.0)

bc_rho_left = DirichletBC(V.sub(0), Constant(1.0), left_boundary)
bc_theta_left = DirichletBC(V.sub(1), Constant(0.0), left_boundary)
bc_rho_right = DirichletBC(V.sub(0), Constant(1.0), right_boundary)
bc_theta_right = DirichletBC(V.sub(1), Constant(pi), right_boundary)

bc = [bc_rho_left, bc_theta_left, bc_rho_right, bc_theta_right]

# Define variational forms
F_rho = alpha * rho * v * dx - beta * rho**3 * v * dx + gamma * dot(grad(rho), grad(v)) * dx - gamma * rho * dot(grad(theta), grad(theta)) * v * dx
F_theta = dot(grad(theta), grad(w)) * dx + 2 * gamma * rho * dot(grad(rho), grad(theta)) * w * dx

F = F_rho + F_theta

# Compute solution
u = Function(V)
solve(F == 0, u, bc)

# Split solution into components
rho, theta = u.split()

# Plot solution
plot(rho)
plot(theta)

# Save solution to file
File("solution.pvd") << u

# Show plots
interactive()

I’m relatively new to Fenics and I really don’t know what’s going wrong here. I’d appreciate some help!
Thanks,
Jaime

Your problem is non-linear, and thus you should not use a TrialFunction, but a function to define rho and theta, i.e.

should be replaced by

u = Function(V)
rho, theta = split(u)
1 Like