How to use correctly split with MixedElements?

variational formulation
1

Hello everyone,
Recently I have encountered a problem when using the split function with MixedElements. When I use w0.split() the result is very similar to the initial condition but the solver gives an error “Unable to successfully call PETSc function ‘MatSetValuesLocal’.” and when I use split(w0) the solver works but there are perturbations with respect to the initial condition when I plot it.

Attached is a snippet of the code:

mesh = IntervalMesh(900, 0, 1.5e-4) 
V1 = FiniteElement("CG", mesh.ufl_cell(), 2)
V2 = FiniteElement("CG", mesh.ufl_cell(), 2)
ME = FunctionSpace(mesh,MixedElement([V1, V2]))
dw = TrialFunction(ME)
q, p = TestFunctions(ME)
w = Function(ME)
w0 = Function(ME)
class InitialConditions(UserExpression):
    def __init__(self, **kwargs):
        super().__init__(**kwargs)
    def eval(self, values, x): #Returns values for a function of dimension two
        if x[0]<xi=4*1.5e-6/6:
            values[0]=1.0
            values[1]=1.0
        else:
            values[0]=0.0
            values[1]=0.0
    def value_shape(self): 
         return (2,)
u_init = InitialConditions()
w0.interpolate(u_init)
c, eta = split(w)
c0, eta0 = split(w0)
fig, ax = pl.subplots()
figure()
plot(eta0)
ax.set_xlim(0.0000975, 0.0001025)
c0, eta0 = w0.split()
figure()
plot(eta0)
ax.set_xlim(0.0000975, 0.0001025)

I would be very grateful if anyone knows why this is happening.

Thanks,
Jesús

2

So, in general, you want to use split(w) to get components for defining a variational form, then w.split() after solving, to write output for visualization. For example, you might follow a workflow like

# (Define mixed function space W)

w = Function(W)
q,p = split(w)

# (Define nonlinear residual F in terms of q and p)

solve(F==0,w,...)

q_viz, p_viz = w.split()
File("q.pvd") << q_viz
File("p.pvd") << p_viz

The variables q and p from split(w) are essentially just convenient labels for w[0] and w[1], while q_viz and p_viz from w.split() are separate Function objects. (The fact that this is a potential source of confusion has been a known issue for some time.)

5 Likes
3

I tried it. It seems that w.split() and split(w) retrun the same results. It is w.split(True) who returns the different result.