Nonlinear finite elements/Kinematics - polar decomposition

Polar decomposition

The w:Polar decomposition theorem states that any second order tensor whose determinant is positive can be decomposed uniquely into a symmetric part and an orthogonal part.

In continuum mechanics, the deformation gradient $𝑭$ is such a tensor because $\det (𝐅)>0$. Therefore we can write

$𝑭=𝑹\cdot 𝑼=𝑽\cdot 𝑹$

where $𝑹$ is an orthogonal tensor ($𝑹\cdot {𝑹}^T=1$) and $𝑼,𝑽$ are symmetric tensors ($𝑼={𝑼}^T$ and $𝑽={𝑽}^T$) called the right stretch tensor and the left stretch tensor, respectively. This decomposition is called the polar decomposition of $𝑭$.

Recall that the right Cauchy-Green deformation tensor is defined as

$𝑪={𝑭}^T\cdot 𝑭$

Clearly this is a symmetric tensor. From the polar decomposition of $𝑭$ we have

$𝑪={𝑼}^T\cdot {𝑹}^T\cdot 𝑹\cdot 𝑼=𝑼\cdot 𝑼={𝑼}^2$

If you know $𝑪$ then you can calculate $𝑼$ and hence $𝑹$ using $𝑹=𝑭\cdot {𝑼}^{-1}$.

How do you find the square root of a tensor?

If you want to find $𝑼$ given $𝑪$ you will need to take the square root of $𝑪$. How does one do that?

We use what is called the spectral decomposition or eigenprojection of $𝑪$. The spectral decomposition involves expressing $𝑪$ in terms of its eigenvalues and eigenvectors. The tensor product of the eigenvectors acts as a basis while the eigenvalues give the magnitude of the projection.

Thus,

$𝑪=\sum _{i=1}^3{\lambda }_i^2{𝑵}_i\otimes {𝑵}_i$

where ${\lambda }_i^2$ are the principal values (eigenvalues) of $𝑪$ and ${𝑵}_i$ are the principal directions (eigenvectors) of $𝑪$.

Therefore,

${𝑼}^2=\sum _{i=1}^3{\lambda }_i^2{𝑵}_i\otimes {𝑵}_i$

Since the basis does not change, we then have

$𝑼=\sum _{i=1}^3{\lambda }_i{𝑵}_i\otimes {𝑵}_i$

Therefore the ${\lambda }_i$ can be interpreted as principal stretches and the vectors ${𝑵}_i$ are the directions of the principal stretches.

Exercise:

If

$𝑼=\sum _{i=1}^3{\lambda }_i{𝑵}_i\otimes {𝑵}_i$

show that

${𝑼}^2=𝑼\cdot 𝑼=\sum _{i=1}^3{\lambda }_i^2{𝑵}_i\otimes {𝑵}_i.$

Example of polar decomposition

Let us assume that the motion is given by

$\begin{array}{cc}{x}_1& =\frac{1}{4}[4X_1+(9-3X_1-5X_2-X_1{X}_2)t]\\ x_2& =X_2+(4+2X_1)t\end{array}$

The adjacent figure shows how a unit square subjected to this motion evolves over time.

File:MotionExample.png An example of a motion.

Deformation gradient

The deformation gradient is given by

$𝑭=\frac{\partial 𝐱}{\partial 𝑿}⟹F_{ij}=\frac{\partial x_i}{\partial X_j}$

Therefore

$\begin{array}{cc}{F}_{11}& =\frac{\partial x_1}{\partial X_1}=\frac{1}{4}[4+(-3-X_2)t]\\ F_{12}& =\frac{\partial x_1}{\partial X_2}=\frac{1}{4}[(-5-X_1)t]\\ F_{21}& =\frac{\partial x_2}{\partial X_1}=2t\\ F_{22}& =\frac{\partial x_2}{\partial X_2}=1\end{array}$

At $t=1$ at the position $𝑿=(0,0)$ we have

$𝐅=\left[\begin{array}{cc}\frac{\partial x_1}{\partial X_1}& \frac{\partial x_1}{\partial X_2}\\ \frac{\partial x_2}{\partial X_1}& \frac{\partial x_2}{\partial X_2}\end{array}\right]=\frac{1}{4}\left[\begin{array}{cc}1& -5\\ 8& 4\end{array}\right]$

You can calculate the deformation gradient at other points in a similar manner.

Right Cauchy-Green deformation tensor

We have

$𝑪={𝑭}^T\cdot 𝑭$

Therefore,

$𝐂={𝐅}^T𝐅=\frac{1}{16}\left[\begin{array}{cc}65& 27\\ 27& 41\end{array}\right]$

To compute $𝑼$ we have to find the eigenvalues and eigenvectors of $𝑪$. The eigenvalue problem is

$(𝐂-{\lambda }^2𝐈)𝐍=\mathrm{𝟎}$

where

$𝐈=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right]$

To find the eigenvalues we solve the characteristic equation

$\det (𝐂-{\lambda }^2𝐈)=0$

Plugging in the numbers, we get

$\det \left[\begin{array}{cc}\frac{65}{16}-{\lambda }^2& \frac{27}{16}\\ \frac{27}{16}& \frac{41}{16}-{\lambda }^2\end{array}\right]=0$

or

${\lambda }^4-\frac{53}{8}{\lambda }^2+\frac{121}{16}=0$

This equation has two solutions

$\begin{array}{cc}{\lambda }_1^2& =\frac{53}{16}+\frac{3}{16}\sqrt{97}=5.159\\ {\lambda }_2^2& =\frac{53}{16}-\frac{3}{16}\sqrt{97}=1.466\end{array}$

Taking the square roots we get the values of the principal stretches

${\lambda }_1=2.2714{\lambda }_2=1.2107$

To compute the eigenvectors we plug into the eigenvalues into the eigenvalue problem to get

$\{\left[\begin{array}{cc}65& 27\\ 27& 41\end{array}\right]-{\lambda }_1^2\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right]\}\left[\begin{array}{c}{N}_1^{(1)}\\ N_2^{(1)}\end{array}\right]=\left[\begin{array}{c}0\\ 0\end{array}\right]$

Because this system of equations is not linearly independent, we need another equation to solve this system of equations for $N_1^{(1)}$ and $N_2^{(1)}$. This problem is eliminated by using the following equation (which implies that $𝐍$ is a unit vector)

$N_2^{(1)}=\sqrt{1-(N_1^{(1)}{)}^2}$

Solving, we get

${𝐍}_1=\left[\begin{array}{c}{N}_1^{(1)}\\ N_2^{(1)}\end{array}\right]=\left[\begin{array}{c}0.8385\\ 0.5449\end{array}\right]$

We can do the same thing for the other eigenvector ${𝐍}_2$ to get

${𝐍}_2=\left[\begin{array}{c}{N}_1^{(2)}\\ N_2^{(2)}\end{array}\right]=\left[\begin{array}{c}-0.5449\\ 0.8385\end{array}\right]$

Therefore,

${𝑵}_1\otimes {𝑵}_1={𝐍}_1{𝐍}_1^T=\left[\begin{array}{c}0.8385\\ 0.5449\end{array}\right]\left[\begin{array}{cc}0.8385& 0.5449\end{array}\right]=\left[\begin{array}{cc}0.7031& 0.4569\\ 0.4569& 0.2969\end{array}\right]$

and

${𝑵}_2\otimes {𝑵}_2={𝐍}_2{𝐍}_2^T=\left[\begin{array}{c}-0.5449\\ 0.8385\end{array}\right]\left[\begin{array}{cc}-0.5449& 0.8385\end{array}\right]=\left[\begin{array}{cc}0.2969& -0.4569\\ -0.4569& 0.7031\end{array}\right]$

Therefore,

$𝑪={\lambda }_1^2{𝑵}_1\otimes {𝑵}_1+{\lambda }_2^2{𝑵}_2\otimes {𝑵}_2⟹𝐂=5.159\left[\begin{array}{cc}0.7031& 0.4569\\ 0.4569& 0.2969\end{array}\right]+1.466\left[\begin{array}{cc}0.2969& -0.4569\\ -0.4569& 0.7031\end{array}\right]$

We usually don't see any problem to calculate $𝑪$ at this point and go straight to the right stretch tensor.

Right stretch

The right stretch tensor $𝑼$ is given by

$𝑼={\lambda }_1{𝑵}_1\otimes {𝑵}_1+{\lambda }_2{𝑵}_2\otimes {𝑵}_2⟹𝐔=2.2714\left[\begin{array}{cc}0.7031& 0.4569\\ 0.4569& 0.2969\end{array}\right]+1.2107\left[\begin{array}{cc}0.2969& -0.4569\\ -0.4569& 0.7031\end{array}\right]$

or

$𝐔=\left[\begin{array}{cc}1.9565& 0.4846\\ 0.4846& 1.5256\end{array}\right]$

We can invert this matrix to get

${𝐔}^{-1}=\left[\begin{array}{cc}0.5548& -0.1762\\ -0.1762& 0.7114\end{array}\right]$

Rotation

We can now find the rotation matrix by using th relation

$𝑹=𝑭\cdot {𝑼}^{-1}$

In matrix form,

$𝐑=\frac{1}{4}\left[\begin{array}{cc}1& -5\\ 8& 4\end{array}\right]\left[\begin{array}{cc}0.5548& -0.1762\\ -0.1762& 0.7114\end{array}\right]=\left[\begin{array}{cc}0.3590& -0.9334\\ 0.9334& 0.3590\end{array}\right]$

You can check whether this matrix is orthogonal by seeing whether $𝐑{𝐑}^T={𝐑}^T𝐑=𝐈$.

You thus get the polar decomposition of $𝑭$. In an actual calculation you have to be careful about floating point errors. Otherwise you might not get a matrix that is orthogonal.