Get stresses from displacement for the thermo-elastic evolution problem (full coupling)

Dear all,

How can I output the stresses in the VTK format for the thermo-elastic evolution problem?

I have augmented the original code for the problem with a definition

def sigma(v, Delta_T):
    return (lmbda * ufl.tr(eps(v)) - kappa * Delta_T) * ufl.Identity(gdim) + 2.0 * mu * eps(v)

and added the following to the time-stepping loop

    stress = sigma(v.sub(0), 0) / 1e6    # convert to MPa
    V_stress = fem.functionspace(domain, ("DG", 0))
    stress_expr = fem.Expression(stress, V_stress.element.interpolation_points())
    stresses = fem.Function(V_stress, name="Stress")

The stresses are then written to VTK files analogously to the temperature variation and displacements.

Which Delta_T should I use to get the thermal stresses?

The stresses as written to the VTK files are identically zero. How to get the correct ones?

Thanks,
Nikolai

You need to call stresses.interpolate(stress_expr) to move the data into the output function before saving.

Hello Jorgen,

Adding stresses.interpolate(stress_expr) causes the following error:

Traceback (most recent call last):
  File "/mnt/c/Trainings/FEniCSx/thermoelasticity_full.py", line 242, in <module>
    stresses.interpolate(stress_expr)
  File "/usr/lib/petsc/lib/python3/dist-packages/dolfinx/fem/function.py", line 459, in interpolate
    _interpolate(u0)
  File "/usr/lib/python3.10/functools.py", line 889, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
  File "/usr/lib/petsc/lib/python3/dist-packages/dolfinx/fem/function.py", line 455, in _
    self._cpp_object.interpolate(e0._cpp_object, cells0, cells1)  # type: ignore
RuntimeError: Function value size not equal to Expression value size.

Per the error message, you’re attempting to interpolate a tensor-valued expression onto a scalar valied function space. Naturally, that doesn’t work.

In that same demo, you’ll find an approach to creating arbitrarily shaped function spaces. E.g.:

Sel = basix.ufl.element( "DG", domain.basix_cell(), 0, shape=(gdim,gdim) )  # Stress element
V_stress = fem.functionspace(domain,Sel)
1 Like