How to apply a constant inhomogeneous initial stress field (in linear elasticity)?

I/O
1

Hi,

I want to apply an initial (or pre / residual) stress field in a 3D continua FEM linear elastic problem.

So far, I have followed the official tutorial for linear elasticity, with the lhs a = ufl.inner(sigma(u), epsilon(v)) * ufl.dx

In comparison with https://fenicsproject.discourse.group/t/how-to-add-a-constant-stress-tensor-to-the-stress-expression-linear-elasticity/3050 , the initial stress field should be included in the right hand side like rhs = inner(my_constant_tensor,epsilon(v))*dx

However, there, the applied stress field is constant over the domain.

What is the best way to initialize “my_constant_tensor” in each element cell / material point? The stress data is a numpy array of dimensions (#material points, stress dimension).

Thank you in advance
Erik

2

See for instance:

3

If this data corresponds to data at a single gauss point per element, do

sig_0 = df.VectorFunctionSpace(mesh, 'DG', 0, dim = 3) # supposing 2d and voigt/mandel notation
sig_0.vector().set_local(data.flatten())
rhs = df.inner( voigt2tensor(sig_0), eps(v))*dx
1 Like
4

Thank you for your suggestion. However, I forgot to mention that I am using dolfinx, it appears to me that this solution holds only for dolfin.

5

Thank you.

The interpolate function seems to do the job. But I am wondering about the dimension x.shape[1] in your example. Printing it gives x.shape[1]=600, with 200 triangle cells. How do I obtain/set the number of integration points in this case?

6

Please provide a minimal reproducible example of what you have codes up so far.

7 
import dolfinx
from mpi4py import MPI
import ufl
import numpy as np
mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)

el = ufl.TensorElement("Lagrange", mesh.ufl_cell(), 1, shape=(2,2))

V = dolfinx.fem.FunctionSpace(mesh, el)
u = dolfinx.fem.Function(V)

def sig(x):
    print(x.shape[1])#=600
    sig_local = np.array([[1., 2.],[3., 4.]])
    values = np.repeat(sig_local, x.shape[1])
    return values.reshape(sig_local.shape[0]*sig_local.shape[1], x.shape[1])

u.interpolate(sig)
8

You are creating a CG-1 quantity, which has one degree of freedom per vertex.
If you want something that is constant per cell, use a DG 0 function space.

1 Like