libmultiscale/test/single_scale/AKANTU/test_2D_Impulse_Al-111-LJ54d9314b39ddfeatures/coupler_compon…
libmultiscale/test/single_scale/AKANTU/test_2D_Impulse_Al-111-LJ
54d9314b39ddfeatures/coupler_compon…
test_2D_Impulse_Al-111-LJ
test_2D_Impulse_Al-111-LJ
readme.tex
readme.tex
\documentclass[a4paper,11pt]{article}
\usepackage[utf8]{inputenc}
\usepackage[french]{babel}
\usepackage[T1]{fontenc} % accents codés dans la fonte
%\usepackage{layout}
\usepackage{a4wide}
\usepackage{amsmath}
\usepackage[usenames]{color}
\usepackage{url}
\usepackage{hyperref}
%\usepackage{comment}
%\usepackage{textcomp}
\usepackage{parskip}
\usepackage{stmaryrd}
\usepackage{tikz}
\usetikzlibrary{calc,patterns,snakes}
\setlength{\hoffset}{0pt}
\setlength{\oddsidemargin}{-1cm} % Marge gauche sur pages impaires
\setlength{\evensidemargin}{-1cm} % Marge gauche sur pages paires
\setlength{\marginparwidth}{0cm} % Largeur de note dans la marge
\setlength{\textwidth}{18cm} % Largeur de la zone de texte (17cm)
\setlength{\voffset}{0pt} % Bon pour DOS
\setlength{\marginparsep}{0pt} % Séparation de la marge
\setlength{\topmargin}{0cm} % Pas de marge en haut
\setlength{\headheight}{0cm} % Haut de page
\setlength{\headsep}{0cm} % Entre le haut de page et le texte
\setlength{\footskip}{1cm} % Bas de page + séparation
\setlength{\textheight}{25.5cm} % Hauteur de la zone de texte (25cm)
\begin{document}
\section{Getting linearized LJ}
We want to fit the parameters from an atomic lennard jones potential:
\
\begin{equation}
\phi(r) = 4 \epsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right]
\end{equation}
But we want to model it in a thin sheet of 111 aluminum plane. Thus the pattern is
an hexagonal lattice with the inter-atomic distance $r_0 = \frac{\sqrt{2}}{2} a$ where $a$ is the
lattice parameter of the FCC structure of aluminum (i.e. $a = 4.0320$ \AA).
The atomic mass of the aluminum is $m_{al} = 26.98153860 g/mol$.
Also we want to take
into account the first neighbors only, so that the energy per atom at ground state is:
\begin{equation}
E = 6 \cdot \phi(r_0) = 24 \epsilon \left[ \left( \frac{\sigma}{r_0} \right)^{12} - \left( \frac{\sigma}{r_0} \right)^{6} \right]
\end{equation}
From this we can compute the length parameter $r_0$ so that we are at the minimum of energy:
\
\begin{align}
\frac{\partial E}{\partial r_0} &= 24 \frac{ \epsilon} {r_0} \left[ -12 \left( \frac{\sigma}{r_0} \right)^{12} + 6 \left( \frac{\sigma}{r_0} \right)^{6} \right] = 0 \\
& \Leftrightarrow 12 \left( \frac{\sigma}{r_0} \right)^{12} = 6 \left( \frac{\sigma}{r_0} \right)^{6} \\
& \Leftrightarrow 2 = \left( \frac{r_0}{\sigma} \right)^{6} \Leftrightarrow 2^{1/6} \sigma = r_0
\end{align}
So that this last relation between $r_0$ and $\sigma$ ensures the ground state equilibrium.\\
If we linearize this LJ potential we have:
\begin{equation}
f(r) = \phi''(r_0)(r-r_0) = \frac{72 \epsilon}{r_0^2} (r - r_0) = k (r - r_0)
\end{equation}
\section{The elastic tensor}
It can be demonstrated \cite{tadmor} that the elastic fourth order tensor can be derived to be:
\
\begin{align}
C_{ijkl}(r_0) &= \frac{1}{\Omega} \sum_{p=1}^6 \left[ \frac{\phi''(r_0)}{r_0^2} - \frac{\phi'(r_0)}{r_0^3} \right] X_i^p X_j^p X_k^p X_l^p \\
&= \frac{1}{\Omega} \sum_{p=1}^6 \left[
\frac{4 \epsilon} {r_0^2} \frac{\left[ 156 \left( \frac{\sigma}{r_0} \right)^{12}
- 42 \left( \frac{\sigma}{r_0} \right)^{6} \right]}{r_0^2} -
\frac{4 \epsilon} {r_0} \frac{\left[ -12 \left( \frac{\sigma}{r_0} \right)^{12}
+ 6 \left( \frac{\sigma}{r_0} \right)^{6} \right]}{r_0^3} \right]
X_i^p X_j^p X_k^p X_l^p \\
&= \frac{1}{\Omega} \sum_{p=1}^6 \frac{4 \epsilon} {r_0^4} \left[
\left[ 156 \left( \frac{\sigma}{r_0} \right)^{12}
- 42 \left( \frac{\sigma}{r_0} \right)^{6} \right] -
\left[ -12 \left( \frac{\sigma}{r_0} \right)^{12}
+ 6 \left( \frac{\sigma}{r_0} \right)^{6} \right] \right]
X_i^p X_j^p X_k^p X_l^p \\
&= \frac{1}{\Omega} \frac{4 \epsilon} {r_0^4} \left[
168 \left( \frac{\sigma}{r_0} \right)^{12}
- 48 \left( \frac{\sigma}{r_0} \right)^{6} \right]
\sum_{p=1}^6 X_i^p X_j^p X_k^p X_l^p \\
&= \frac{1}{\Omega} \xi(r_0) \sum_{p=1}^6 X_i^p X_j^p X_k^p X_l^p
\end{align}
With $X^p$ the position of atom $p$ within the 6 neighbors and $\Omega$ is the lattice volume.
\
\begin{align}
\xi(r_0) &= \frac{4 \epsilon} {r_0^4} \left( \frac{\sigma}{r_0} \right)^{6} \left[ 168 \left( \frac{\sigma}{r_0} \right)^{6} - 48 \right] \\
&= \frac{4 \epsilon} {r_0^4} \left( 2^{-1/6} \right)^{6} \left[ 168 \left( 2^{-1/6} \right)^{6} - 48 \right] \\
&= \frac{2 \epsilon} {r_0^4} \left[ 84 - 48 \right] = \frac{72 \epsilon} {r_0^4}
\end{align}
where we used the relation between $r_0$ and $\sigma$.
\begin{align}
C_{1111} &= \frac{1}{\Omega} \frac{72 \epsilon} {r_0^4} \sum_{p=1}^6 (X_1^p)^4 \\
C_{2222} &= \frac{1}{\Omega} \frac{72 \epsilon} {r_0^4} \sum_{p=1}^6 (X_2^p)^4 \\
C_{1212} &= \frac{1}{\Omega} C_{1122} = \frac{1}{\Omega} \frac{72 \epsilon} {r_0^4} \sum_{p=1}^6 (X_1^p)^2 (X_2^p)^2 \\
\end{align}
Let us then get interested in the sums $\sum_{p=1}^6 (X_1^p)^4$
$\sum_{p=1}^6 (X_2^p)^4$ and $\sum_{p=1}^6 (X_1^p)^2 (X_2^p)^2$.
For this let us assume that the orientation is that the X axis lattice is of length $r_0$
and the Y axis is of length $\sqrt{3}/2 r_0$ forming the lattice presented in figure
\ref{toto}.
Then the positions of the 6 atoms are:
$(r_0,0,0)$, $(-r_0,0,0)$,
$\frac{1}{2}(r_0,\sqrt{3} r_0,0)$ , $\frac{1}{2}(r_0,-\sqrt{3} r_0,0)$,
$\frac{1}{2}(-r_0,\sqrt{3} r_0,0)$ and $\frac{1}{2}(-r_0,-\sqrt{3} r_0,0)$.
the first sum is:
\begin{equation}
\sum_{p=1}^6 (X_1^p)^4 = 2 r_0^4 + 4 \frac{1}{2^4} r_0^4 = r_0^4 (2 + \frac{1}{4}) = r_0^4 \left(\frac{3}{2}\right)^2
\end{equation}
the second sum is:
\begin{equation}
\sum_{p=1}^6 (X_2^p)^4 = 2 \cdot 0 + 4 \frac{3^{4/2}}{2^4} r_0^4 = r_0^4 \left(\frac{3}{2} \right)^2
\end{equation}
And the third is:
\begin{equation}
\sum_{p=1}^6 (X_1^p)^2 (X_2^p)^2 = 2 \cdot r_0^2 \cdot 0^2 + 4 \cdot \frac{1}{2^4} r_0^4 3 = r_0^4 \frac{3}{4}
\end{equation}
Finally we have that:
\begin{align}
C_{1111} &= C_{2222} = \frac{1}{\Omega} \frac{72 \epsilon} {r_0^4} r_0^4 \left(\frac{3}{2}\right)^2 = \frac{1}{\Omega} 162 \epsilon \\
C_{1212} &= C_{1122} = \frac{1}{\Omega} \frac{72 \epsilon} {r_0^4} r_0^4 \frac{3}{4} = \frac{1}{\Omega} 54 \epsilon \\
\end{align}
\section{fitting the orthotropic law}
let us recall the orthotropic Hooke's law:
\
\begin{align}
\left(
\begin{array}{c}
\sigma_{11}\\
\sigma_{22}\\
\sigma_{33}\\
\sigma_{23}\\
\sigma_{31}\\
\sigma_{12}\\
\end{array}
\right) = \Gamma
\left(
\begin{array}{cccccc}
E_1 (1-\nu_{23}\nu_{32}) & E_1 (\nu_{21} + \nu_{31}\nu_{23}) & E_1 (\nu_{31} + \nu_{21}\nu_{32}) & 0 & 0 & 0 \\
E_1 (\nu_{21} + \nu_{31}\nu_{23}) & E_2 (1 - \nu_{13}\nu_{31}) & E_2 (\nu_{32} + \nu_{12}\nu_{31}) & 0 & 0 & 0 \\
E_3 (\nu_{31} + \nu_{21}\nu_{32}) & E_3 (\nu_{32} + \nu_{12}\nu_{31}) & E_3 (1 - \nu_{12}\nu_{21}) & 0 & 0 & 0 \\
0 & 0 & 0 & \mu_{23}/\Gamma & 0 & 0 \\
0 & 0 & 0 & 0 & \mu_{13}/\Gamma & 0 \\
0 & 0 & 0 & 0 & 0 & \mu_{12}/\Gamma \\
\end{array}
\right)
\left(
\begin{array}{c}
\epsilon_{11}\\
\epsilon_{22}\\
\epsilon_{33}\\
\epsilon_{23}\\
\epsilon_{31}\\
\epsilon_{12}\\
\end{array}
\right)
\end{align}
with $\Gamma = \frac{1}{1 - \nu_{12}\nu_{21} - \nu_{23}\nu_{32} - \nu_{31}\nu_{13} - 2 \nu_{21}\nu_{32}\nu_{13}}$
In 2D with plane strain it reduces to:
\begin{align}
\left(
\begin{array}{c}
\sigma_{11}\\
\sigma_{22}\\
\sigma_{12}\\
\end{array}
\right) = \Gamma
\left(
\begin{array}{cccccc}
E_1 & E_1 \nu_{21} & 0 \\
E_2 \nu_{12} & E_2 & 0 \\
0 & 0 & \mu/\Gamma \\
\end{array}
\right)
\left(
\begin{array}{c}
\epsilon_{11}\\
\epsilon_{22}\\
\epsilon_{12}\\
\end{array}
\right)
\end{align}
with $\Gamma = \frac{1}{1 - \nu_{12}\nu_{21}}$
So that we have the relations:
\
\begin{align}
C_{1111} &= \Gamma E_1 = \frac{1}{\Omega} 162 \epsilon \\
C_{2222} &= \Gamma E_2 = \frac{1}{\Omega} 162 \epsilon
\end{align}
which implies that $E_1 = E_2 = E$. We have consequently:
\begin{align}
C_{1122} &= \Gamma E \nu_{21} = \frac{E \nu_{21}}{1 - \nu_{12}\nu_{21}} = \frac{1}{\Omega} 54 \epsilon \\
&= C_{2211} = \Gamma E \nu_{12} = \frac{E \nu_{12}}{1 - \nu_{12}\nu_{21}} = \frac{1}{\Omega} 54 \epsilon \\
\end{align}
Then we have $\nu_{12} = \nu_{21} = \nu$. Let us summarize:
\begin{align}
C_{1111} &= C_{2222} = \frac{E}{1 - \nu^2} = \frac{1}{\Omega} 162 \epsilon \\
C_{1122} &= C_{2211} = \frac{E \nu}{1 - \nu^2} = \frac{1}{\Omega} 54 \epsilon
\end{align}
So that:
\begin{align}
\nu &= \frac{C_{1122}}{C_{1111}} = \frac{54}{162} = \frac{1}{3} \\
E &= \frac{1}{\Omega} \frac{8}{3} \cdot 54 \epsilon = \frac{1}{\Omega} 144 \epsilon \\
\mu &= \frac{1}{\Omega} 54 \epsilon = \frac{3}{8} E
\end{align}
Let us have $\Omega = \frac{r_0^3\sqrt{3}}{2} [m^3]$ with the z dimension taken as $r_0$.
So that let say we want a young's modulus of 68 GPa we have:
\begin{align}
E = 68 \cdot 10^9 N.m^{-2} \\
\epsilon = r_0^3\sqrt{3} \cdot 2.361111111 \cdot 10^8 Joules \\
\mu = 2.55 \cdot 10^{10} N.m^{-2} \\
\nu = 1/3
\end{align}
By choosing the lattice parameter $a = 4.0320$ we obtain $r_0 = \frac{\sqrt{2}}{2} a$ with
$a = 4.0320 \cdot 10^{-10}$ \AA, we obtain in real units
\begin{align}
a &= 4.0320 \cdot 10^{-10} \\
r_0 & = 2.8510545417441596184 \cdot 10^{-10} \\
E &= 68 \cdot 10^9 N.m^{-2} \\
\epsilon &= 9.4774964211648290153 \cdot 10^{-21} Joules \\
\mu &= 2.55 \cdot 10^{10} N.m^{-2} \\
\nu &= 1/3
\end{align}
By converting the distances to \AA~and the Joules to kcal/mol we have finally
\begin{align}
a &= 4.0320 [A] \\
r_0 & = 2.8510545417441596184 [A] \\
\sigma &= 2.5400008365880643002 \\
E &= 9.7874197348334701168 [kcal/mol/A^3] \\
\epsilon &= 1.3641211104312093062 [kcal/mol] \\
\mu &= 3.6702824005625512938 [kcal/mol/A^3] \\
\nu &= 1/3 \\
k &= \frac{72 \epsilon}{r_0^2} = 12.082988860820660663 [kcal/mol/A^2] \\
\rho &= \frac{m_{al}}{\frac{\sqrt{3}}{2} r_0^2 }
\end{align}
\end{document}
c4science · Help