Hi,

I’m solving an advection-diffusion-equation in a mixed system from which I obtain the concentrations for three substances, c1, c2 and c3. From these concentrations I want to calculate saturation functions with the following formula:

Where Ct is a threshold, so a constant, and phi is the saturation function, and C is the concentration c1,c2 or c3. I tried this (for c3 in this example):

```
element = FiniteElement('P', triangle, 1)
Vphi = FunctionSpace(mesh, element)
phi = Function(Vphi)
ptest = TestFunction(Vphi)
thr = Constant(0.002)
form = ((c3**2)/((c3**2)+(thr**2)))*ptest*dx
phi = assemble(form)
```

This works, but in this case phi is of the type ‘dolfin.cpp.la.Vector’ which I want to transfer back to a Function. The other option I tried is:

```
thr = 0.002
c3_vec = c3.vector().get_local()
phi_vec = (c3**2)/((c3**2)+(thr**2))
phi.vector().set_local(phi_vec)
```

This does work normally, but since c3 comes from a mixed functionspace I can’t convert it to c3_vec. The error I get is: ‘Unable to access vector of degrees of freedom. cannot access a non-const vector from a subfunction.’

So my question is: 1. Which of the two methods is correct, and 2. How do I solve the described problems?

Thanks in advance