Normal vector of just one of the boundaries

Dear FEniCSx community,

I am trying to compute a stress vector in one of the boundaries of my geometry. In order to achieve this, I need to define the stress tensor and the normal vector associated to that boundary, and take the dot product of both variables.

The problem is that if I define the normal vector as usual, i.e., n=FacetNormal(mesh), I get the normal vector of all the boundaries. Therefore, if I take the dot product of the stress tensor and the normal vector, I get the stress vector over all the boundaries.

Is it possible to define the normal vector of, for example, the left boundary?

Here is a minimal example:

# Import
import dolfinx
from ufl import Identity, div, dot, ds, dx, inner, lhs, grad, nabla_grad, rhs, sym, FacetNormal
from ufl import  FiniteElement, VectorElement, TrialFunctions, TestFunctions, MixedElement
from mpi4py import MPI

# Mesh
mesh = dolfinx.UnitSquareMesh(MPI.COMM_WORLD, 10, 10)
n = FacetNormal(mesh)

# Function spaces
v_cg2 = VectorElement('CG', mesh.ufl_cell(), 2)
p_cg1 = FiniteElement('CG', mesh.ufl_cell(), 1)
mel = MixedElement([v_cg2, p_cg1])
V = dolfinx.FunctionSpace(mesh, v_cg2)
Q = dolfinx.FunctionSpace(mesh, p_cg1)

u_n = dolfinx.fem.Function(V)
p_n = dolfinx.fem.Function(Q)

# Stress tensor
mu = 0.01
def sigma(u, p):
    return 2*mu*sym(nabla_grad(u)) - p*Identity(u.geometric_dimension())

# Stress vector
stress = dot(sigma(u_n,p_n),n)

Many thanks in advance.

What do you want to use stress for later on in your code?
If you want to integrate the stress over a subset of a boundary, you can just restrict the integration measure ds

The problem is that stress over one of the boundaries is what I need (it is the gradient of my cost function, so I just need to transform it into numpy array and use scipy.minimize). Therefore, I do not need to integrate it.

How are you going to transform this into a numpy array?

I would suggest using a projection into a suitable function space, which you can do as:

which shows you how to get a projection of the facet normal over a certain boundary marked with a mesh tag.
You can of course adapt this to take in your stress and project in to a suitable space.

Thank you very much for your reply.
I have adapted the function facet_normal_approximation as you suggested.

Now, when I compute the stress vector as stress = dot(sigma(u_n,p_n),nn) I get the object <class 'ufl.tensoralgebra.Dot'>.

Do you know how I can transform it into a numpy array?

Many thanks again.

As I asked previously, what do you want that numpy array to contain? Is it values at degrees of freedom?

Yes, that is exactly what I need.

I would:

  1. Project the expression stress into a suitable function space
  2. Get the dof coordinates of the suitable space by using “FunctionSpace.tabulate_dof_coordinates()”
  3. Find what dofs are on your chosen boundary by using locate_dofs_topological or locate_dofs_geometrical.
  4. Extract the dof coordinates at these dofs and the corresponding values from the projected function
1 Like

@JLorenteMacias Could you please share your full code for calculating the stresses in this unit square…?