Dear firedrakers, assume I want to do the following solve of a mixed system: V1 = FunctionSpace(mesh, 'RT', 1) V2 = FunctionSpace(mesh, 'DG', 0) W = V1 * V2 lmbda = 1 u, p = TrialFunctions(W) v, q = TestFunctions(W) f = Function(V1) g = Function(V2) a = (p*q - q*div(u) + lmbda*inner(v, u) + div(v)*p)*dx L = (dot(u,f)+p*g)*dx u = Function(W) solve(a == L, u, solver_parameters=...) But now I take out the line L = (dot(u,f)+p*g)*dx and assume that somewhere I calculate instead: L1 = dot(TestFunction(V1),f)*dx L2 = TestFunction(V2)*g*dx L1 and L2 taken together contain the same information as L, but I can’t write L=L1+L2, because L1 and L2 were not defined using the mixed function spaces (assume that I calculated L1 and L2 somewhere else in the code, and I only use the mixed function spaces in this subroutine). So given L1 and L2, how can I build L, which can be used in the solve? What I currently do is this: m_u = dot(TestFunction(V1),TrialFunction(V1))*dx m_p = TestFunction(V2)*TrialFunction(V2)*dx f_u = Function(V1) f_p = Function(V2) solve(m_u=L1,f_u) solve(m_p=L2,f_p) L = (dot(u,f_u)+p*f_g)*dx which I believe does the right thing, but it involves mass solves. Is there a better way? Can I copy the data over directly? Thanks a lot, Eike -- Dr Eike Hermann Mueller Research Associate (PostDoc) Department of Mathematical Sciences University of Bath Bath BA2 7AY, United Kingdom +44 1225 38 5803 e.mueller@bath.ac.uk http://people.bath.ac.uk/em459/