Writing symmetric tensor for paraview visualization

Hello,

I would like to write a 2D symmetric tensor as a vector (XX, YY, 0, XY, 0, 0) with pvd files to use the tensor glyph representation of paraview.

Something like:

from fenics import *


mesh = UnitSquareMesh(2,2)
V = TensorFunctionSpace(mesh, 'P' ,1)

tensor = as_tensor([[1,0], [0,1]])
tensor = project(tensor, V)

vector = as_vector([tensor[0,0], tensor[1,1], 0, tensor[0,1], 0, 0])

file = File('test.pvd') 
file << vector #error

I guess the problem here is that vector is not projected in any FunctionSpace
Moreover the symmetry parameter in TensorFunctionSpace seems broken (see Writing symmetric tensor function fails).

Anyone has a workaround ?

Hi,
if you just export the tensor directly it should be working no ?
Best

I tried to export it directly but I can’t apply the Tensor Glyph filter on the tensor data.
Tensor Glyph apparently expect a format (XX, YY, ZZ, XY, YZ, XZ) (https://kitware.github.io/paraview-docs/latest/python/paraview.simple.TensorGlyph.html)

What happens if you define V as a Function Space of 3x3 tensor using shape=(3,3) ?

Adding shape=(3,3) leads to a shape mismatch :

Shapes do not match: <Argument id=140284654524520> and <ListTensor id=140283568815184>.
---------------------------------------------------------------------------
UFLException                              Traceback (most recent call last)
<ipython-input-1-9e086fc1b63d> in <module>
      6
      7 tensor = as_tensor([[1,0], [0,1]])
----> 8 tensor = project(tensor, V)
      9
     10 vector = as_vector([tensor[0,0], tensor[1,1], 0, tensor[0,1], 0, 0])

~/Dev/miniconda3/envs/cellfem/lib/python3.7/site-packages/dolfin/fem/projection.py in project(v, V, bcs, mesh, function, solver_type, preconditioner_type, form_compiler_parameters)
    127     Pv = TrialFunction(V)
    128     a = ufl.inner(w, Pv) * dx
--> 129     L = ufl.inner(w, v) * dx
    130
    131     # Assemble linear system

~/.local/lib/python3.7/site-packages/fenics_ufl-2019.1.0-py3.7.egg/ufl/operators.py in inner(a, b)
    167     if a.ufl_shape == () and b.ufl_shape == ():
    168         return a * Conj(b)
--> 169     return Inner(a, b)
    170
    171

~/.local/lib/python3.7/site-packages/fenics_ufl-2019.1.0-py3.7.egg/ufl/tensoralgebra.py in __new__(cls, a, b)
    159         ash, bsh = a.ufl_shape, b.ufl_shape
    160         if ash != bsh:
--> 161             error("Shapes do not match: %s and %s." % (ufl_err_str(a), ufl_err_str(b)))
    162
    163         # Simplification

~/.local/lib/python3.7/site-packages/fenics_ufl-2019.1.0-py3.7.egg/ufl/log.py in error(self, *message)
    170         "Write error message and raise an exception."
    171         self._log.error(*message)
--> 172         raise self._exception_type(self._format_raw(*message))
    173
    174     def begin(self, *message):

UFLException: Shapes do not match: <Argument id=140284654524520> and <ListTensor id=140283568815184>.

But i found a workaround with MixedFunctionSpace:

from fenics import *


mesh = UnitSquareMesh(2,2)

U = FiniteElement("CG", mesh.ufl_cell(), 1)
ME  = FunctionSpace(mesh, U*U*U*U*U*U)

tensor = as_tensor([[1,0], [0,2]])

vector = as_vector([tensor[0,0], tensor[1,1], 0, tensor[0,1], 0, 0])
vector = project(vector, ME)

# file = File('test.xdmf') 
file = XDMFFile('test.xdmf')
file.write(vector)
file.close()

And finally the visualization using the filter Tensor Glyph:

OK. in my above suggestion I was thinking about changing to:

tensor = as_tensor([[1, 0, 0], [0, 1, 0], [0, 0, 0]])