JuAFEM shape gradient as a matrix

Specific Domains Modelling & Simulations
1

Hello everybody,

I hope I am asking the question in a correct way. I am working with JuAFEM and I have a example code given, which calculates the stiffness matrix elementwise.

@inbounds for q_point in 1:getnquadpoints(cellvalues_u) # loop over quadrature points
        dΩ = getdetJdV(cellvalues_u, q_point)
        for i in 1:n_basefuncs_u # loop over number of ansatz functions
            BuT = symmetric(shape_gradient(cellvalues_u, q_point, i))
            for j in 1:i 
                Bu = symmetric(shape_gradient(cellvalues_u, q_point, j))
                Ke[i, j] += BuT ⊡ C ⊡ Bu * dΩ 
            end
        end
    end

Is there a way to calculate the shape gradient such that I do not need to calculate the stiffness matrix elementwise but for each element in matrix form. In particular I do not want to access each element of the stiffness matrix (Ke[i,j]) but calulate the matrix itself i.e.

@inbounds for q_point in 1:getnquadpoints(cellvalues_u) # loop GPi
        dΩ = getdetJdV(cellvalues_u, q_point)
        Bu = ? # Calculate somehow the shape gradient Bu and BuT by a function
        BuT = transpose(Bu)
        Ke += BuT * C * Bu * dΩ 
 end
2

The entries of Bu are such that column i is the Voigt form of

grad_Ni = symmetric(shape_gradient(cellvalues_u, q_point, i))

So you could do something like

# This is for 2D where symmetric second order tensors have 3 entries
Bu = zeros(3, n_basefuncs_u)
for i in 1:n_basefuncs_u #
    Bu[:, i] = tovoigt(symmetric(shape_gradient(cellvalues_u, q_point, i)))
end

This notebook has an (old) example of three different ways of computing Ke in JuAFEM. I’m not sure if it still works.:

I don’t really understand why you would want to code like this though. What’s the advantage over the first way?

1 Like
3

Hi Kristoffer.carlsson,

thank you for your answer. I will try this out and hope it will work. The function “tovoigt” is it a function already implemented in JuAFEM?

Actually for normal FEM applications I guess there is no advantage. But I need the shape gradient matrices for later purposes.

4

I have to agree with Kristoffer. In Julia there is no advantage to using big matrices to formulate the elementwise finite element expressions. I did that in GitHub - PetrKryslUCSD/FinEtools.jl: Finite Element tools in Julia and I had to use buffers all over the place in order to avoid temporary allocations (which will eventually kill performance, if one is not careful).

JuAFEM does the right thing by calculating the element wise matrices entry by entry (most of the time, if I recall). StaticArrays really shine when used this way.

5

The function “tovoigt” is it a function already implemented in JuAFEM?

Sorry, it comes from Tensors.jl, you can add it and do using Tensors to get access to it.

1 Like