Nonlinear finite elements/Rate form of hyperelastic laws

From testwiki
Revision as of 06:24, 11 April 2021 by 182.70.218.143 (talk) (For the spatial elasticity tensor a_{ijkl}, there was issues with indeices. Those have been corrected.)
(diff) ← Older revision | Latest revision (diff) | Newer revision β†’ (diff)
Jump to navigation Jump to search

Rate equations

Let us now derive rate equations for a hyperelastic material.

First elasticity tensor

We start off with the relation

𝑷=ρ0W𝑭

Then the material time derivative of 𝑷 is given by

First elasticity tensor

𝑷˙=𝖠:𝑭˙;𝖠:=ρ02W𝑭𝑭

where the fourth order tensor 𝖠 is call the first elasticity tensor. This tensor has major symmetries but not minor symmetries. In index notation with respect to an orthonormal basis

PΛ™iJ=𝖠iJkLFΛ™kL=ρ02WFiJFkLFΛ™kL

Proof:

We have

𝑷˙=ρ0t(W𝑭)=ρ0𝑭(Wt)=ρ0𝑭(W𝑭:𝑭˙)

Using the product rule, we have

𝑷˙=ρ02W𝑭𝑭:𝑭˙+ρ0W𝑭:2𝑭𝑭t=ρ02W𝑭𝑭:𝑭˙+ρ0W𝑭:t(𝑭𝑭)

Therefore,

𝑷˙=ρ02W𝑭𝑭:𝑭˙=𝖠:𝑭˙

Second elasticity tensor

Similarly, if we start off with the relation

𝑺=2ρ0Wπ‘ͺ

the material time derivative of 𝑺 can be expressed as

Second elasticity tensor

𝑺˙=𝖒:12π‘ͺΛ™;𝖒=4ρ02Wπ‘ͺπ‘ͺ

where the fourth order tensor 𝖒 is called the material elasticity tensor or the second elasticity tensor. Since this tensor relates symmetric second order tensors it has minor symmetries. It also has major symmetries because the two partial derivatives are with the same quantity and an interchange does not change things. In index notation with respect to an orthonormal basis

SΛ™IJ=12𝖒IJKLCΛ™KL=12(4ρ02WCIJCKL)CΛ™KL

Proof:

We have

𝑺˙=2ρ0t(Wπ‘ͺ)=2ρ0π‘ͺ(Wt)=2ρ0π‘ͺ(Wπ‘ͺ:π‘ͺΛ™)

Again using the product rule, we have

𝑺˙=2ρ02Wπ‘ͺπ‘ͺ:π‘ͺΛ™+2ρ0Wπ‘ͺ:2π‘ͺπ‘ͺt=2ρ02Wπ‘ͺπ‘ͺ:π‘ͺΛ™

Therefore,

𝑺˙=4ρ02Wπ‘ͺπ‘ͺ:12π‘ͺΛ™=𝖒:12π‘ͺΛ™

Relation between first and second elasticity tensors

The first and second elasticity tensors are related by

𝖠iJkL=δikSJL+FiN𝖒NJMLFkM

Proof:

Recall that the first and second Piola-Kirchhoff stresses are related by

𝑷=𝑭𝑺

Taking the material time derivative of both sides gives

𝑷˙=𝑭˙𝑺+𝑭𝑺˙

Using the expression for 𝑺˙ above, we get

𝑷˙=𝑭˙𝑺+12𝑭(𝖒:π‘ͺΛ™)

Now

π‘ͺ=𝑭T𝑭π‘ͺΛ™=𝑭T˙𝑭+𝑭T𝑭˙

Therefore,

𝑷˙=𝑭˙𝑺+12𝑭[𝖒:(𝑭T˙𝑭+𝑭T𝑭˙)]

Now

𝖒:(𝑭T˙𝑭)𝖒IJKLFΛ™mKFmLand𝖒:(𝑭T𝑭˙)=𝖒IJLKFmLFΛ™mK=𝖒IJKLFmLFΛ™mK

That means

𝖒:(𝑭T˙𝑭)=𝖒:(𝑭T𝑭˙)

which gives us

𝑷˙=𝑭˙𝑺+𝑭[𝖒:(𝑭T𝑭˙)]

In index notation,

PΛ™iJ=FΛ™iLSLJ+FiN[𝖒NJMLFkMFΛ™kL]=SLJδikFΛ™kL+𝖒NJMLFiNFkMFΛ™kL=(δikSJL+FiN𝖒NJMLFkM)FΛ™kL=𝖠iJkLFΛ™kL

Therefore,

𝖠iJkL=δikSJL+FiN𝖒NJMLFkM

Fourth elasticity tensor

Now we will compute the spatial elasticity tensor for the rate constitutive equation for a hyperelastic material. This tensor relates an objective rate of stress (Cauchy or Kirchhoff) to the rate of deformation tensor. We can show that

Fourth elasticity tensor for the Kirchhoff stress

β„’φ[τ]=𝖼:𝒅𝖼τ:𝒅

where

𝖼ijkl=FiKFjLFkMFlN𝖒KLMN=4ρ0FiKFjLFkMFlN2WCKLCMN.

The fourth order tensor 𝖼 is called the spatial elasticity tensor or the fourth elasticity tensor. Clearly, 𝖼ijkl cannot be derived from the store energy function W because of the dependence on the deformation gradient.

Proof:

Recall that the Lie derivative of the Kirchhoff stress is defined as

β„’φ[τ]=τ=𝑭𝑺˙𝑭T

We have found that

𝑺˙=𝖒:12π‘ͺΛ™

We also know from Continuum_mechanics/Time_derivatives_and_rates#Time_derivative_of_strain that

π‘ͺΛ™=2𝑭T𝒅𝑭

where 𝒅 is the spatial rate of deformation tensor. Therefore,

β„’φ[τ]=𝑭[𝖒:(𝑭T𝒅𝑭)]𝑭T

In index notation,

(β„’φ[τ])ij=FiK𝖒KLMNFkMdklFlNFjL=FiKFjLFkMFlN𝖒KLMNdkl=𝖼ijkldkl

or,

β„’φ[τ]=𝖼:𝒅

where

𝖼ijkl=FiKFjLFkMFlN𝖒KLMN

Alternatively, we may define 𝖼 in terms of the Cauchy stress σ, in which case the constitutive relation is written as

Fourth elasticity tensor for the Cauchy stress

β„’φ[σ]=𝖼:𝒅=J1𝖼τ:𝒅

where

𝖼ijkl=J1FiKFjLFkMFlN𝖒KLMN

The proof of this relation between the spatial and material elasticity tensors is very similar to that for the rate of Kirchhoff stress. Many authors define this quantity 𝖼ijkl as the spatial elasticity tensor. Note the factor of J1. This form of the spatial elasticity tensor is crucial for some of the calculations that follow.

Relation between first and fourth elasticity tensors

The first and fourth elasticity tensors are related by

𝖠iJkL=FJj1[𝖼ijkl+τjlδik]FLl1

In the above equation 𝖼ijkl is the elasticity tensor that relates the rate of Kirchhoff stress to the rate of deformation.

Instead, if we use the Cauchy stress and the spatial elasticity tensor 𝖼ijkl that relates the Cauchy stress to the rate of deformation), the above relation becomes

𝖠iJkL=JFJj1[𝖼ijkl+σjlδik]FLl1

Proof:

Recall that

𝑷˙=𝖠:𝑭˙

Therefore,

𝑷˙𝑭T=(𝖠:𝑭˙)𝑭T

Also recall that

𝖠iJkL=δikSJL+FiN𝖒NJMLFkM

Therefore, using index notation,

PΛ™iJFlJ=δikSJLFlJFΛ™kL+FiN𝖒NJMLFkMFlJ)FΛ™kL

Now,

𝒍=𝑭˙𝑭1or𝑭˙=𝒍𝑭

In index notation

FΛ™kL=lkmFmL

Using this we get

PΛ™iJFlJ=δikSJLFlJlkmFmL+FiN𝖒NJMLFkMFlJ)lkmFmL

or,

PΛ™iJFlJ=(δikFlJSJLFmL+FiNFlJFkMFmL𝖒NJML)lkm

Now,

τ=𝑭𝑺𝑭Torτlm=FlJSJLFmL

Therefore,

PΛ™iJFlJ=(δikτlm+FiNFlJFkMFmL𝖒NJML)lkm

Also,

𝖼ilkm=FiNFlJFkMFmL𝖒NJML

So we have

PΛ™iJFlJ=(δikτlm+𝖼ilkm)lkm=𝖺ilkmlkm

Note:

The fourth order tensor

𝖺ijkl=δikτjl+𝖼ijkl

which depends on the symmetry of τ is called the third elasticity tensor, i.e.,

𝑷˙𝑭T=𝖺:𝒍

Therefore, the relation between the first and third elasticity tensors is

𝑷˙𝑭T=(𝖠:𝑭˙)𝑭T=𝖺:𝒍=𝖺:(𝑭˙𝑭1)

or,

𝖠:𝑭˙=[𝖺:(𝑭˙𝑭1)]𝑭T

In index notation

𝖠iJkLFΛ™kL=𝖺ijklFΛ™kLFLl1FjJ1

Therefore,

𝖠iJkL=𝖺ijklFLl1FjJ1=(δikτjl+𝖼ijkl)FLl1FjJ1

An isotropic spatial elasticity tensor?

An isotropic spatial elasticity tensor cannot be derived from a stored energy function if the constitutive relation is of the form

β„’φ[τ]=𝖼:𝒅

where

𝖼=λ11+2μ𝖨s.

Since a significant number of finite element codes use such a constitutive equation, (also called the equation of a hypoelastic material of grade 0) it is worth examining why such a model is incompatible with elasticity.

Start with a constant and isotropic material elasticity tensor

Let us start of with an isotropic elastic material model in the reference configuration. The simplest such model is the St. Venant-Kirchhoff hyperelastic model

𝑺=λtr(𝑬)1+2μ𝑬

where 𝑺 is the second Piola-Kirchhoff stress, 𝑬 is the Lagrangian Green strain, and λ,μ are material constants. We can show that this equation can be derived from a stored energy function.

Taking the material time derivative of this equation, we get

𝑺˙=λtr(𝑬˙)1+2μ𝑬˙

Now,

𝑺˙=𝖒:12π‘ͺΛ™=𝖒:𝑬˙

where 𝖒 is the second (material) elasticity tensor.

Therefore,

𝖒:𝑬˙=λtr(𝑬˙)1+2μ𝑬˙

which implies that

𝖒=λ11+2μ𝖨s.

In index notation,

𝖒IJKL=λδIJδKL+μ(δIKδJL+δJKδIL)

Now, from the relations between the second elasticity tensor and the fourth (spatial) elasticity tensor, we have

𝖼ijkl=FiIFjJFkKFlL𝖒IJKL

Therefore, in this case,

𝖼ijkl=FiIFjJFkKFlL[λδIJδKL+μ(δIKδJL+δJKδIL)]

or,

𝖼ijkl=λbijbkl+μ(bikbjl+bilbjk)

where 𝒃=𝑭𝑭T. So we see that the spatial elasticity tensor 𝖼 cannot be a constant tensor unless 𝒃=1.

Alternatively, if we define

𝖼ijkl=J1FiIFjJFkKFlL𝖒IJKL

we get

𝖼ijkl=J1λbijbkl+J1μ(bikbjl+bilbjk)

Start with a constant and isotropic spatial elasticity tensor

Let us now look at the situation where we start off with a constant and isotropic spatial elasticity tensor, i.e.,

𝖼=λ11+2μ𝖨s

In index notation,

𝖼ijkl=λδijδkl+μ(δikδjl+δjkδil)

Since

𝖼ijkl=FiIFjJFkKFlL𝖒IJKL

multiplying both sides by FAi1FBj1FCk1FDl1 we have,

FAi1FBj1FCk1FDl1𝖼ijkl=𝖒ABCD

Therefore, substituting in the expression for a constant and isotropic 𝖼, we have

𝖒ABCD=FAi1FBj1FCk1FDl1[λδijδkl+μ(δikδjl+δjkδil)]

or,

𝖒ABCD=λFAj1FBj1FCl1FDl1+μ(FAk1FBl1FCk1FDl1+FAl1FBk1FCk1FDl1)

or,

𝖒ABCD=λ(𝑭1𝑭T)AB(𝑭1𝑭T)CD+μ(𝑭1𝑭T)AC(𝑭1𝑭T)BD+𝑭1𝑭T)AD(𝑭1𝑭T)BC)

Since π‘ͺ=𝑭T𝑭 which gives us π‘ͺ1=𝑭1𝑭T, we can write

𝖒ABCD=λπ‘ͺAB1π‘ͺCD1+μ(π‘ͺAC1π‘ͺBD1+π‘ͺAD1π‘ͺBC1)

Alternatively, if we define

𝖼ijkl=J1FiIFjJFkKFlL𝖒IJKL

we get

𝖒ABCD=JFAi1FBj1FCk1FDl1𝖼ijkl

and therefore,

𝖒ABCD=Jλπ‘ͺAB1π‘ͺCD1+Jμ(π‘ͺAC1π‘ͺBD1+π‘ͺAD1π‘ͺBC1)

Hypoelastic material of grade 0

A hypoelastic material of grade zero is one for which the stress-strain relation in rate form can be expressed as

β„’φ[σ]=𝖼:𝒅

where 𝖼 is constant. When the material is isotropic we have

𝖼=λ11+2μ𝖨s

We want to show that hypoelastic material models of grade 0 cannot be derived from a stored energy function. To do that, recall that

𝖒IJKL=4ρ02WCIJCKL

and

SIJ=2ρ0WCIJ

For a material elasticity tensor 𝖒 to be derivable from a stored energy function it has to satisfy the Bernstein integrability conditions. We have

SIJCKL=2ρ02WCIJCKL=12𝖒IJKL

Also, due to the interchangeability of derivatives,

CMN(SIJCKL)=CKL(SIJCMN)

Therefore,

𝖒IJKLCMN=𝖒IJMNCKL

These integrability conditions have to be satisfied by any material elasticity tensor.

At this stage we will use the relation

𝖒ABCD=Jλπ‘ͺAB1π‘ͺCD1+Jμ(π‘ͺAC1π‘ͺBD1+π‘ͺAD1π‘ͺBC1)

If we plug this into the integrability condition we will see that

(λ+μ)CKL1(CIM1CJN1+CIN1CJM1)=(λ+μ)CMN1(CIK1CJL1+CIL1CJK1)

If we multiply both sides by CMNCIKCJL we are left with

λ+μ=0

This is an unphysical situation and hence shows that a hypoelastic material of grade zero requires that λ=μ for it to be derivable from a stored energy function.

Proof:

Let us simplify the notation by writing 𝑨:=π‘ͺ1. Then,

𝖒IJKL=JλAIJAKL+Jμ(AIKAJL+AJKAIL)

Then,

𝖒IJKLCMN=λJCMNAIJAKL+Jλ[AIJCMNAKL+AIJAKLCMN]+μJCMN(AIKAJL+AJKAIL)+Jμ[AIKCMNAJL+AIKAJLCMN+AJKCMNAIL+AJKAILCMN]

and

𝖒IJMNCKL=λJCKLAIJAMN+Jλ[AIJCKLAMN+AIJAMNCKL]+μJCKL(AIMAJN+AJMAIN)+Jμ[AIMCKLAJN+AIMAJNCKL+AJMCKLAIN+AJMAINCKL]

As this stage we use the identities (see Nonlinear finite elements/Kinematics#Some_useful_results for proofs)

JCIJ=J2CIJ1=J2AIJ

and

CIJ1CKL=12(CIK1CJL1+CJK1CIL1)orAIJCKL=12(AIKAJL+AJKAIL)

Therefore we have

𝖒IJKLCMN=Jλ2AMNAIJAKLJλ2[(AIMAJN+AJMAIN)AKL+(AKMALN+ALMAKN)AIJ]+Jμ2AMN(AIKAJL+AJKAIL)Jμ2[(AIMAKN+AKMAIN)AJL+(AJMALN+ALMAJN)AIK+(AJMAKN+AKMAJN)AIL+(AIMALN+ALMAIN)AJK]

and

𝖒IJMNCKL=Jλ2AKLAIJAMNJλ2[(AIKAJL+AJKAIL)AMN+(AKMALN+ALMAKN)AIJ]+Jμ2AKL(AIMAJN+AJMAIN)Jμ2[(AIKALM+AKMAIL)AJN+(AJKALN+AJLAKN)AIM+(AJKALM+AJLAKM)AIN+(AIKALN+AILAKN)AJM]

Equating the two, we see that the terms that cancel out are

Jλ2AKLAIJAMN;Jλ2(AKMALN+ALMAKN)AIJ;

and

Jμ2[(AIMAKN+AKMAIN)AJL+(AJMALN+ALMAJN)AIK+(AJMAKN+AKMAJN)AIL+(AIMALN+ALMAIN)AJK]

Therefore,

𝖒IJKLCMN=𝖒IJMNCKL

implies that

Jλ2(AIMAJN+AJMAIN)AKL+Jμ2AMN(AIKAJL+AJKAIL)=Jλ2(AIKAJL+AJKAIL)AMN+Jμ2AKL(AIMAJN+AJMAIN)

or,

(λ+μ)AMN(AIKAJL+AJKAIL)=(λ+μ)AKL(AIMAJN+AJMAIN)

In other words,

(λ+μ)CKL1(CIM1CJN1+CIN1CJM1)=(λ+μ)CMN1(CIK1CJL1+CIL1CJK1)

Now, if we multiply both sides by CMN we get

(λ+μ)CKL1(δINCJN1+δIMCJM1)=(λ+μ)δMM(CIK1CJL1+CIL1CJK1)

or,

2(λ+μ)CKL1CJI1=3(λ+μ)(CIK1CJL1+CIL1CJK1)

Next, multiplying both sides by CIK gives

2(λ+μ)δILCJI1=3(λ+μ)(δKKCJL1+δKLCJK1)

or,

2(λ+μ)CJL1=12(λ+μ)CJL1

Finally, multiplying both sides by CJL gives

(λ+μ)δJJ=6(λ+μ)δJJ

Therefore,

λ+μ=0

Template:Subpage navbar