Stitch two CG2 solutions together (combined mesh is known)

Question

I have two meshes and their combination, say

from dolfin import *
mesh1 = RectangleMesh(0,0,1,1,3,3) #and
mesh2 = RectangleMesh(1,0,2,1,3,3)
mesh = RectangleMesh(0,0,1,2,6,3) #combination

I wish to combine the two solutions on a combined mesh

expr = Expression("sin(x[0])*cos(x[1])")
u1 = interpolate(expr,FunctionSpace(mesh1,'CG',2))
u2 = interpolate(expr,FunctionSpace(mesh2,'CG',2))
u = ?
meshplot = RectangleMesh(0,0,1,2,12,6)
uplot = interpolate(u,FunctionSpace(meshplot,'CG',1))
plot(uplot)

How should I construct u? It is easy for a CG1 field and I have played around enough to have an idea of how to reverse engineer the FEniCS CG2 ordering, but I was hoping for an easy solution...

Answer 1

Hi, in the following code functions living in spaces defined over parts of mesh are glued together by the $L^2$ projection. If you don't insist on using interpolation, perphaps it is a way to go

from dolfin import *

f = Expression('sin(x[0])*cos(x[1])')

# Combined mesh and space on it
mesh = RectangleMesh(-1, -1, 1, 1, 10, 10)
V = FunctionSpace(mesh, 'CG', 2)

# Get mesh parts
cell_f = CellFunction('size_t', mesh, 0)
for cell in cells(mesh):
    if cell.midpoint().x() > 0:
        cell_f[cell] = 1

mesh0 = SubMesh(mesh, cell_f, 0)
mesh1 = SubMesh(mesh, cell_f, 1)

# Create functions in parts
V0 = FunctionSpace(mesh, 'CG', 2)
u0 = interpolate(f, V0)

V1 = FunctionSpace(mesh, 'CG', 2)
u1 = interpolate(f, V1)

# Begin glueing
u = TrialFunction(V)
v = TestFunction(V)
a = inner(u, v)*dx

# Function v0 is u0 in mesh0 and 0 otherwise
L0 = inner(u0, v)*dx(0, subdomain_data=cell_f)
v0 = Function(V)
solve(a == L0, v0)

# Function v1 is u1 in mesh1 and 0 otherwise
L1 = inner(u1, v)*dx(1, subdomain_data=cell_f)
v1 = Function(V)
solve(a == L1, v1)

# Combine
v = Function(V, v0.vector())
v.vector()[:] += v1.vector()
plot(v)

# See about the error
w = interpolate(f, V)
plot(w)

v.vector().axpy(-1, w.vector())
print v.vector().norm('l2')

interactive()