Hi everyone,
I'm using Fenics to solve a system of coupled Cahn-Hilliard equations and it works fine. However, for a certain (meaningful and important) range of input parameters, the solution becomes unstable and diverges. I have tried load ramping and bubble stabilization, to no avail. Therefore I want to use the adaptive solver. However, I get the following error:
/usr/local/lib/python2.7/dist-packages/dolfin/fem/adaptivesolving.py in __init__(self, problem, goal)
101 # Generate error control object
--> 102 ec = generate_error_control(self.problem, goal)
/usr/local/lib/python2.7/dist-packages/dolfin/fem/adaptivesolving.py in generate_error_control(problem, goal)
145 # Generate UFL forms to be used for error control
--> 146 (ufl_forms, is_linear) = generate_error_control_forms(problem, goal)
/usr/local/lib/python2.7/dist-packages/dolfin/fem/adaptivesolving.py in generate_error_control_forms(problem, goal)
214 # Get DOLFIN's error control generator to generate all forms
--> 215 generator = DOLFINErrorControlGenerator(primal, goal, u)
/usr/local/lib/python2.7/dist-packages/dolfin/fem/errorcontrolgenerator.py in __init__(self, F, M, u)
---> 50 ErrorControlGenerator.__init__(self, __import__("dolfin"), F, M, u)
/usr/local/lib/python2.7/dist-packages/ffc/errorcontrol/errorcontrolgenerators.pyc in __init__(self, module, F, M, u)
71 assert(len(self.lhs.arguments()) == 2)
72 assert(len(self.rhs.arguments()) == 1)
---> 73 assert(len(self.goal.arguments()) == 0)
74 assert(len(self.weak_residual.arguments()) == 1)
75
AssertionError:
There appears to be a problem with the rank of my goal functional. Below you can find my (shortened) code. It is based on the Fenics documentation and on this (old) example.
# Build the function space
BasicElement = FiniteElement('Lagrange', mesh.ufl_cell(), 1)
BubbleElement = FiniteElement('Bubble', mesh.ufl_cell(), 2)
Element = BasicElement + BubbleElement
CompleteElement = MixedElement([Element, Element, Element, Element, Element])
Space = FunctionSpace(mesh, CompleteElement)
# Create test & trial functions
w = TestFunction(Space)
u = TrialFunction(Space)
# Add Dirichlet boundary conditions
bc = []
# Build the equation and the goal functional
w0, w1, w2, w3, w4 = split(w)
xv, muv, U, xa, mua = split(u)
equation =
goal = xv*dx # xv is the variable that has the most complex behaviour and strongest variation over the domain
# Set the initial solution
sol = Function(Space) # function that contains the solution
initial_conditions = []
initial_conditions.append()
init_sol = project(Expression(initial_conditions, degree=1), Space)
sol.assign(init_sol)
# Solve
F = equation
R = action(F, sol)
DR = derivative(R, sol) # Gateaux derivative
problem = NonlinearVariationalProblem(R, sol, bc, DR)
solver = AdaptiveNonlinearVariationalSolver(problem, goal)
solver_tol = 1.0e-6
solver.solve(solver_tol)
Does anyone know what I am doing wrong? Keep in mind that it works fine with a non-adaptive solver and that the problem also persists when disabling the bubble elements.
Thank you for your answers.
Best regards,
dvlaethe
