Re: [Getfem-users] Stiffness matrix for linear elasticity

[Renard Yves]

Roman Putanowicz <address@hidden> a écrit :

> Dear Yves,
> > This formula implements in fact a weak term. This represents
> >
> > lambda div(u).div(v) + mu grad(u):grad(v) + mu grad(u)^T:grad(v)
> > y
> > where A:B is the Frobenius (scalar) product of matrices.
> Thank you for your comment, it was much helpful. Somehow I couldn't
> see the forest for the trees.
>
> Now I can link the implementation and the formal derivations.
> Below I include the explanation of this link I have built for myself.
> It is rather informal so please correct me if I to brutally violate
> some math rules :-)
>
>
> For linearised elasticity (assuming epsilon = 1/2 (grad(u) + grad^T(u))
> we have:
>
> Cauchy stress:
>   sigma =   lambda div(u) I + 2 mu (1/2 (grad(u) + grad^T(u)))
>
>   sigma =   lambda div(u) I + mu grad(u) + mu grad^T(u)
>
> The force balance equation (neglecting inertia term and body forces)
>   div(sigma) = 0
>
> Weak form of weighted residual equation
>   int{div(sigma):v} ->  int{sigma:grad(v)} + boundary term
>
> Considering the term under the integral sign
>
>   sigma:grad(v) ->
>   lambda div(u) I:grad(v) + mu grad(u):grad(v) + mu grad^T(u):grad(v)
>
> but as I:grad(v) = div(v)
>
> Then we get:
>
>   lambda div(u).div(v) + mu grad(u):grad(v) + mu grad^T(u): grad(v)

This is correct, yes.


> Well now I see the correspondence with
>
>   lambda.t(:,i,i,:,j,j) + mu.t(:,i,j,:,i,j) + mu.t(:,j,i,:,i,j)
>
> where: t = comp(vGrad(#1).vGrad(#1))
>
> Somehow I couldn't get this using the index notation (maybe was not patient
> enough :).
>
> I decided to write this (now obvious :) elaborate derivation because
> the above way the to calculate the stiffness matrix for linearised elasticity
> is a bit different from what one can find in the elementary books on FEM.
>
> Is it true that this particular way of calculation is driven
> a) by the availability of tensor t = comp(vGrad(#1).vGrad(#1))
> b) avoidance of problems with incompressible materials ?

In fact, it just seems to us that this was the simpler way to  
implement linear elasticity.

For incompressible materials, lambda goes to infinity and this does  
not avoid the
standard problems (locking phenomena). The best way in this case is to  
introduce a
a mixed problem with the pressure as a separate unknown.


Yves.