Nonlinear solves with mixed function spaces

Question

When solving a nonlinear equation using the solve(F == 0) interface, I'm not quite aware of how to properly arrange the terms.
In particular, two things are unclear to me:

• The expression up[0] for a Function of a mixed function space does not return the first component?
• split() only does a shallow copy, right?
• I eventually bump into an All terms in form must have same rank. error, with no idea what it's caused by.

The following simple example highlights the issues.

from dolfin import *

mesh = UnitSquareMesh(20, 20)

V = FunctionSpace(mesh, 'CG', 2)
W = MixedFunctionSpace([V, V])
P = FunctionSpace(mesh, 'CG', 1)

WP = W*P

up = Function(WP)

v = TestFunctions(WP)

## Doesn't work:
# Dot product requires non-scalar arguments, got arguments with ranks 0 and 1.
#F = dot(up[0], v[0]) * dx

# Alright then.
u, p = up.split()
F = dot(u, v[0]) * dx

# Fail:
# ufl.log.UFLException: All terms in form must have same rank.
solve(F == 0, up)

Answer 1

You've got Function.split(), static UFL split() and few __getitem__()s - I'm everytime a little bit confused. Nevertheless try this

up = Function(WP)
vq = TestFunction(WP) # Note TestFunction, not TestFunctions

u, p = as_vector((up[0], up[1])), up[2]
#u, p = split(up) # or use this

v, q = as_vector((vq[0], vq[1])), vq[2]
#v, q = split(vq) # or use this
#v, q = TestFunctions(WP) # or either this

F = dot(u, v) * dx

solve(F == 0, up)

You see, that splits and Test/TrialFunctions split the functions up to the first subspace level, while indexing using __getitem__() operator gets you down through all the subspace levels.

Comments

So split() on up doesn't do the same thing as the as_vector() construction?

---

u, p = split(up) works too. On the other hand u, p = up.split() does not work for some reason. It's UFL/DOLFIN magic.

---

Oh geez. Thanks for the hint!

---

So I guess none of this works anymore.
MWE:

from dolfin import *

mesh = UnitSquareMesh(20, 20)

W = VectorFunctionSpace(mesh, 'CG', 2)
P = FunctionSpace(mesh, 'CG', 1)

WP = MixedFunctionSpace([W, P])

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

up = Function(WP)

vq = TestFunction(WP)  # Note TestFunction, not TestFunctions
u, p = as_vector((up[0], up[1])), up[2]
#u, p = split(up)  # or use this
v, q = as_vector((vq[0], vq[1])), vq[2]
#v, q = split(vq)  # or use this
#v, q = TestFunctions(WP)  # or either this
F = dot(u, v) * dx - dot(f, v) * dx
solve(F == 0, up)

All of the suggested options result in error messages like

[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Arguments are incompatible!
[0]PETSC ERROR: Local column size 14009 cannot be larger than global column size 12894!