Where are the strains evaluated when using the "grad" function?

Hello,

I have a quick question. Let’s say that I have obtained a displacement u. Now, I want to do postprocessing and plot the strains and stress. To do that, I use the ufl.grad function on the displacement u. I need to convert the (strain and stress) values obtained to a NumPy array. My question is, how can I get the strains as a NumPy array? Also, would the reported strains be evaluated at the Gauss points? My goal is to have NumPy arrays for the stress and strain with the corresponding Gauss points.

one row of data would ideally be
GaussP_x_coord, GaussP_x_coord, E11, E22, E12, S11, S22, S12

The total number of rows should be equal to the total number of Gauss points

Many thanks in advance

Project the stresses or strains into a quadrature space (i would do this per stress or strain component, as it is easy to then map them to the tabulated dof coordinates).

Thanks @dokken, I am relatively new to FEniCS, may you direct me to the function name or an example?

Consider:

import dolfin


mesh = dolfin.UnitSquareMesh(10, 10)
V = dolfin.VectorFunctionSpace(mesh, "Lagrange", 2)
u = dolfin.Function(V)
u.interpolate(dolfin.Expression(("2*x[0]*x[0]", "x[1]"), degree=2))

degree = 2
el = dolfin.FiniteElement("Quadrature", mesh.ufl_cell(),
                          degree=degree, quad_scheme="default")
Q = dolfin.FunctionSpace(mesh, el)
qh = dolfin.Function(Q)

p = dolfin.TrialFunction(Q)
q = dolfin.TestFunction(Q)
dx = dolfin.Measure("dx", domain=mesh, metadata={
    "representation": "quadrature",
    "quadrature_degree": degree})
a = dolfin.inner(p, q)*dx
stress = dolfin.sym(dolfin.grad(u))
L = dolfin.inner(stress[0, 0], q)*dx
dolfin.solve(a == L, qh)

coords = Q.tabulate_dof_coordinates()
vals = qh.vector().get_local()
for (coord, val) in zip(coords, vals):
    print(f"Coordinate {coord}, Value {val} Exact value {4 * coord[0]}")

1 Like

Many thanks @dokken, I am using FEniCS 2019, and I keep getting this warning
*** FFC: quadrature representation is deprecated! It will ***
*** likely be removed in 2018.2.0 release. Use uflacs ***
*** representation instead. ***
*** =========================

For the specific problem above uflacs does not work, so you can just ignore the warning.

1 Like

Thank you very much, @dokken !

Hi Dokken,

I have the same question on this. Does the dof coordinates mean the nodal coordinates? For linear element (constant strain element), is it possible to output the strain tensor element-wise? Thanks.

If you want the stresses per element I would project them into a DG-0 space in legacy dolfin, or use interpolate with DOLFINx expression in DOLFINx: Implementation — FEniCSx tutorial

How about the strain tensor? I think it is a bit different from the Von-mises case probably. I want to output the strain tensor for each element and the corresponding element number. Since for linear element, the stain is consider inside the element, but it is still a tenor or matrix. Can I still project then into a DG-0 space? I am quite new to Fenics.

Yes, you’d project or interpolate the stress tensor using an appropriate space. For example the stress would be in the nonconforming discontinuous space \sigma_h \in \left[\text{DG}^{h,p-1}\right]^{d\times d} if the displacement is in some standard conforming space \vec{u}_h \in \left[V^{h,p}\right]^d for spatial dimension d and polynomial approximation order p.

Thank you for the quick response. I am talking about the strain tensor. Is that just the same as stress tensor in terms of how to output?

Can fenics achieve this? Like el = dolfin.FiniteElement (“Quadrature”, mesh.ufl_cell(), degree=degree, quad_scheme=“default”). I don’t have dolfin installed.