### Interaction of MixedFunction's sub(), [], derivative and diff

331
views
0
10 months ago by
I would expect all the vectors assembled in the attached code snippet to be different from zero (ideally the same result, possibly with different permutations and an additional zero for the case with the MixedElement). However, the result varies with different permutations of .sub(), as_vector, diff and derivative. What I'm missing, and why some of the results are equal to zero? I'm particularly bothered by the apparent need to use as_vector when dealing with a MixedElement.

File attached: Test.py (2 KB)

Community: FEniCS Project

0
10 months ago by
Not all the vectors are supposed to be zero. This can be seen from the following 1D example:
from fenics import *

mesh = UnitIntervalMesh(1)
V = FunctionSpace(mesh, 'P', 1)

u = Function(V)
v = TestFunction(V)

f = u.dx(0)

F = derivative(f, u, v)
print(as_backend_type(assemble(F * dx)).vec().view())

F = diff(f, u) * v
print(as_backend_type(assemble(F * dx)).vec().view())​
For
$f(u) = \frac{\partial u}{\partial x},$
you have in the first case
$F = \delta [ f (u) v ] [v] = \frac{\partial}{\partial k} \left. \int \frac{\partial}{\partial x} (u + k v) v \,\mathrm{d} x \right|_{k=0} = \int \frac{\partial v}{\partial x} v \,\mathrm{d} x,$
which is non-zero, as can be seen from the shape functions. In the second case, you have
$F = \frac{\partial^2 u}{\partial x^2} v,$
which leads to zero values, since $u \equiv 0$.
Thank you, I was convinced that all the vectors of my example should have been different from 0. Anyway: why the following code gives a null vector?
from fenics import *
mesh = UnitSquareMesh(1, 1)
UF3_ELEMENT = VectorElement("CG", mesh.ufl_cell(), 1, 2)
R3_ELEMENT = VectorElement("R", mesh.ufl_cell(), 0, 1)
UF3_6_ELEMENT = MixedElement(UF3_ELEMENT, R3_ELEMENT)
UF3_6_SPACE = FunctionSpace(mesh, UF3_6_ELEMENT)
UF3_6 = Function(UF3_6_SPACE)
VU3_6 = TestFunction(UF3_6_SPACE)
F = grad(UF3_6.sub(0))
delta_F = derivative(F, UF3_6, VU3_6)
LinearForm = inner(delta_F, Identity(2)) * dx
print(as_backend_type(assemble(LinearForm)).vec().view())
written 10 months ago by Marco Morandini
I'm not sure but it looks like F is somehow not a function of UF3_6 anymore. The following seems to work
from fenics import *

mesh = UnitSquareMesh(1, 1)

P1 = VectorElement('P', mesh.ufl_cell(), 1)
R0 = FiniteElement('R', mesh.ufl_cell(), 0)
mixed_element = MixedElement(P1, R0)
V = FunctionSpace(mesh, mixed_element)

u = Function(V)
(v0, _) = TestFunctions(V)

u0 = u.sub(0)

print(as_backend_type(assemble(l)).vec().view())​