Use "replace()" on a form to solve the problem on a mesh.

Asked by Simon Funke


I would like to implement some kind of multigrid scheme, which solves problems on different meshes to accelerate the solution procedure.

My strategy is following.
The user specifies the weak formulation on the coarsest mesh.
The multigrid scheme then takes this setup and iteratively:
1) applies the "refine" function to refine the mesh;
2) creates new functionspaces on that new mesh
3) uses the "replace" function to update and solve the weak formulation on the fine mesh.

To demonstrate what I mean have a look at this code.
It solves the Poisson equation first on the coarse mesh and then on the fine mesh.
Unfortunately it doesnt work (see error message below):

==================== Python code, python-dolfin version 1.0.0-7 on Ubuntu 12.10 ===============================
from dolfin import *

# Solve the coarse resolution problem
mesh = UnitSquare(2, 2)

V = FunctionSpace(mesh, "CG", 1)
u = TrialFunction(V)
v = TestFunction(V)
sol = Function(V)

s = Constant(1.0)
F = inner(grad(u), grad(v))*dx - s*v*dx
bc = DirichletBC(V, Constant(0.0), "on_boundary")

solve(lhs(F) == rhs(F), sol, bc)

# Solve the fine resolution problem
mesh_fine = refine(mesh)

V_fine = FunctionSpace(mesh_fine, "CG", 1)
u_fine = TrialFunction(V_fine)
v_fine = TestFunction(V_fine)
sol_fine = Function(V_fine)

# Note 1: I get "Expecting the solution variable u to be a member of the trial space.":
F_fine = replace(F, {u:u_fine, v:v_fine})

# Note 2: It works fine if I redefine the form (but this is what I try to avoid):
# F_fine = inner(grad(u_fine), grad(v_fine))*dx - s*v_fine*dx

bc_fine = DirichletBC(V_fine, Constant(0.0), "on_boundary")
solve(lhs(F_fine) == rhs(F_fine), sol_fine, bc_fine)

This is the error message I get:

*** -------------------------------------------------------------------------
*** Error: Unable to define linear variational problem a(u, v) = L(v) for all v.
*** Reason: Expecting the solution variable u to be a member of the trial space.
*** Where: This error was encountered inside LinearVariationalProblem.cpp.
*** ------------------------------------------------------------------------

Any idea what the issue is and how to fix it?
I suspect that the form F_fine still carries some information about the coarse mesh dimensions that did not get updated with the replace function.



P.S. Two other questions on the side:
Is there a way to "replace" the mesh of a FunctionSpace (e.g. V_fine = replace(V, mesh_fine))?
Is there a way to "replace" the function space for the boundary condition (e.g. bc_fine = replace(bc, V_fine))?

Question information

English Edit question
DOLFIN Edit question
No assignee Edit question
Solved by:
Simon Funke
Last query:
Last reply:
Revision history for this message
Marie Rognes (meg-simula) said :

Hi Simon,

This is clearly a bug: it should work. Could you please file a bug report? The really strange aspect is that it runs if you change both the mesh and the polynomial degree, or the polynomial degree, but not if you just change the mesh.

dolfin-devs: Could this be a caching or instant issue?

Revision history for this message
Marie Rognes (meg-simula) said :

Also, to answer your other questions:

1. No, there is no way to replace the mesh of a FunctionSpace. But you might want to take a look at the adapt(...) functionality, see

For instance, you can do something like the following pseudo-code:

  new_mesh = mesh.adapt()
  adapt(V, new_mesh)

2. No, but again see adapt (as for 1.)

Revision history for this message
Simon Funke (simon-funke) said :

This bug is now fixed in the trunk.

For the discussion see: