# How to calculate Euler-Lagrange equations in the linear elasticity equation?

**URL:** https://community.freefem.org/t/how-to-calculate-euler-lagrange-equations-in-the-linear-elasticity-equation/1678
**Category:** General Discussion
**Created:** [April 17, 2022, 9:26am UTC](https://community.freefem.org/t/how-to-calculate-euler-lagrange-equations-in-the-linear-elasticity-equation/1678 "2022-04-17T09:26:53Z")
**Posts on this page:** 2
**Page:** 1

<div class="post-metadata">

### Author: ![sosome](https://avatars.discourse-cdn.com/v4/letter/s/2bfe46/32.png) [@sosome](https://community.freefem.org/u/sosome)
#### Post date: [April 17, 2022, 9:26am UTC](https://community.freefem.org/t/how-to-calculate-euler-lagrange-equations-in-the-linear-elasticity-equation/1678/1 "2022-04-17T09:26:53Z")

</div>

Dear all,

I am trying to euler-lagrange equation in the linear elasticity equation, but I don’t know how to calculate energy form.  
 ![energy_form](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/4/48a01580e1c7d9349d3e9b6ca398a7942f5bb1f4.png)  
σ = stress, ε = strain, Ω = boundary, Γ = domain geometry

```auto
//Parameters
real Rho = 8000.; //Density
real E = 210.e9; //Young modulus
real Nu = 0.27; //Poisson ratio

real Gravity = -9.81;	//Gravity

//Mesh
real nn = 10; //Mesh quality
real L = 20.; //Beam length
real H = 1.; //Beam height
int Fixed = 1; //Beam fixed label
int Free = 2; //Beam free label
border b1(t=0., L){x=t; y=0.; label=Free;};
border b2(t=0., H){x=L; y=t; label=Fixed;};
border b3(t=0., L){x=L-t; y=H; label=Free;};
border b4(t=0., H){x=0.; y=H-t; label=Fixed;};

int nnL = max(2., nn*L);
int nnH = max(2., nn*H);
mesh Th = buildmesh(b1(nnL) + b2(nnH) + b3(nnL) + b4(nnH));

//Fespace
func Pk = P2;
fespace Uh(Th, [Pk, Pk]);
Uh [ux, uy];

//Macro
real sqrt2 = sqrt(2.);
macro Epsilon(ux, uy) [dx(ux), dy(uy), (dy(ux)+dx(uy))/sqrt2] //
macro Divergence(ux, uy) (dx(ux) + dy(uy)) //

//Problem
real Mu = E/(2.*(1.+Nu));
real Lambda = E*Nu/((1.+ Nu)*(1.-2.*Nu));

varf vElasticity ([ux,uy], [vx, vy])
	= int2d(Th)(
		  Lambda * Divergence(vx, vy) * Divergence(ux, uy)
		+ 2. * Mu * (
			  Epsilon(vx, vy)' * Epsilon(ux, uy)
		)
	)
	+ int2d(Th)(
		  Rho * Gravity * vy
	)
	+ on(Fixed, ux=0, uy=0)
	;

matrix<real> Elasticity = vElasticity(Uh, Uh, solver=sparsesolver);
real[int] ElasticityBoundary = vElasticity(0, Uh);
ux[] = Elasticity^-1 * ElasticityBoundary;

//Movemesh
Th = movemesh(Th, [x+ux, y+uy]);
[ux, uy] = [ux, uy];

//Plot
plot([ux, uy], value=true, cmm="u");

```

[https://modules.freefem.org/modules/elasticity/](https://modules.freefem.org/modules/elasticity/)

---

<div class="post-metadata">

### Author: ![frederichecht](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/frederichecht/32/15_2.png) [@frederichecht](https://community.freefem.org/u/frederichecht)
#### Post date: [April 18, 2022, 8:49am UTC](https://community.freefem.org/t/how-to-calculate-euler-lagrange-equations-in-the-linear-elasticity-equation/1678/2 "2022-04-18T08:49:07Z")

</div>

This is trivial

```auto

real energie =
	= int2d(Th)(
	2* Lambda * Divergence(ux, uy) * Divergence(ux, uy)
		+ 4 * Mu * (
			  Epsilon(ux, uy)' * Epsilon(ux, uy)
		)
	)
	+ int2d(Th)(
		  Rho * Gravity * uy
	)

```