# Forces for the Embedded Atom Model

It is not easy to find the forces written down for embedded atom model.

The potential is always written as

$U_i = F_i(\rho) + \frac{1}{2}\sum_{j \neq i} \phi(r_{ij})$

with

$\rho = \sum_{j\neq i} f_j(r_{ij})$

Probably the reason for forces written down is that different fitting functions for $F$, $\rho$, $\phi$ are used depending on the researcher.  For this reason, LAMMPS uses splines for their implementation.

As usual, the force is calculated by

$\boldsymbol F_i = \nabla U_i$

$\varphi_{ij} = \frac{\partial U}{\partial r_{ij}} = -\frac{\partial F(\rho_i)}{\partial \rho}\frac{\partial f_i}{\partial r_{ij}} \nabla_i r_{ij} -\frac{\partial F(\rho_j)}{\partial \rho}\frac{\partial f_j}{\partial r_{ij}}\nabla_i r_{ij}- \frac{1}{2}\frac{\partial \phi}{\partial r_{ij}}\nabla_i r_{ij}$

Cai and Ye suggest the following forms will in the following forms (with derivatives added by me). Note that we have replaced $r_{ij}$ with a more compact notation given by $r$.

$f(r) = f_e \mathrm{exp}\left(-\chi (r- r_e) \right)$

$\frac{\partial f}{\partial r} = -f_e \chi \mathrm{exp}\left(-\chi (r - r_e) \right)$

$F(\rho) = -F_0\left[1- \mathrm{ln}\left(\frac{\rho}{\rho_e} \right)^n \right]\left(\frac{\rho}{\rho_e} \right)^n + F_1 \left(\frac{\rho}{\rho_e} \right)$

$\frac{\partial F}{\partial \rho} = -F_0n\left( \frac{\rho}{\rho_e} \right)^{n-1}\left(\frac{1}{\rho_e} - 1 \right) + F_0n^2 \frac{\rho^{n-1}}{\rho_e^n}\mathrm{ln}\left(\frac{\rho}{\rho_e} \right) + \frac{F_1}{\rho_e}$

$\phi(r) = -\alpha\left[1 + \beta\left(\frac{r}{r_a}-1 \right) \right]\mathrm{exp}\left(\beta\left(\frac{r}{r_a}-1 \right) \right)$

$\frac{\partial \phi}{\partial r} = - \alpha \beta \left[1 + \beta\left(\frac{r}{r_a} -1 \right)\right] \mathrm{exp}\left(\beta\left(\frac{r}{r_a} -1 \right) \right) + -\alpha \frac{\beta}{r_a}\mathrm{exp}\left(\beta\left(\frac{r}{r_a} -1 \right) \right)$

# A Brief Introduction to Phase Field Modeling of Solidification

This is copy and pasted from my Master’s thesis, and is a very short introduction into phase field modeling of solidification. It is my notes on Boettingers 2002 review article Phase-Field Simulation of Solidification
A free energy functional governing solid-liquid phase transition is defined as follows
$F = \int_V \left[ f(\phi, c, T) - \frac{\epsilon_c^2}{2}\vert \nabla c\vert^2 - \frac{\epsilon_\phi^2}{2}\vert \nabla \phi\vert^2 \right] \mathrm{dV}$
where $f$, $c$ and $\phi$ are the free energy density, concentration and phase field, respectively with $\epsilon_c$ and $\epsilon_\phi$ are the associated gradient entropy coefficients. The free energy density $f$ contains at least one minima for each phase (a two phase system is typically a double well). The time dependent changes can be derived from variational techniques in combination of mass continuity \cite{Boettinger2002} and are given by
$\frac{\partial \phi}{\partial t} = -M_{\phi}\left(\frac{\partial f}{\partial \phi} - \epsilon_\phi^2 \nabla^2 \phi \right)$
$\frac{\partial c}{\partial t} = \nabla \cdot \left(M_{c} c(1-c) \nabla \left(\frac{\partial f}{\partial c} - \epsilon_c^2 \nabla^2 c \right)\right) \label{eq:phase_conserved}$
where $M_{\phi}$ and $M_{c}$ are positive mobilities related to the material kinetics and the diffusion process. The difference between the two equations is due to the fact that $c$ (in the bottom equation) is conserved while $\phi$ (the top equation) is not. If the solid and liquid free energies are available the free energy density can be constructed as
$f(\phi, c, T) = \left[(1-c)W_A + cW_B\right] g(\phi) + (1-p(\phi))f^S(c,t) + p(\phi) f^L(c,T)$
where $g(\phi)$ is a double well, $W_A$ is the energy barrier between material A liquid phase and solid phase, $W_B$ is the energy barrier between material B liquid and solid phase, $p(\phi)$ is the interpolation function–typically $\phi$ or $C_1 \int g(\phi) \mathrm{d}\phi$ where $C_1$ is a constant. $f^S(c,t)$ is the solid phase free energy and $f^L(c,t)$ is the liquid phase free energy. For pure elements, the following simplification is frequently made. Take the solid state as the standard state $f_S^A(T) = 0$ and expand the difference between the liquid and solid free energies around the melting point
$f_L^A(T) - f_S^A(T) = \frac{L(T_M^A - T)}{T_M^A}$
where $T_M^A$ is the melting temperature of the material. This gives the free energy
$f_A(\phi,T) = W_A g(\phi) + L_A \frac{T-M^A -T}{T_M^A} p(\phi)$
The free enthalpy which is used to derive the heat equation is described as
$h = h_0 + C_p T + L \phi$
This gives our heat and phase field equations for a pure element as follows
$C_p \frac{\partial T}{\partial t} + L \frac{\partial \phi}{\partial t} = \nabla \cdot (k\nabla T)$
$\frac{\partial \phi}{\partial t} = M_\phi \epsilon_\phi\left(\nabla^2 \phi - \frac{2 W_A}{\epsilon_\phi} g'(\phi)\right) - \frac{M_\phi L}{T_M}(T_M - T)$
The reader is directed to the review paper from Boettinger et al. for further detail and examples.

I can never remember the Laplacians, gradients and divergence formulas in spherical and cylindrical coordinates. Usually I use Griffiths Introduction to Electrodynamics, however I don’t always carry the book with me. As such I am leaving them here as a useful place for me to reference. While there are many places to reference they are usually not front in center.

I hope you find them useful too.

##### Cylindrical Polar Coordinates

$\nabla U = \frac{\partial U}{\partial r}\boldsymbol{e}_r + \frac{1}{r}\frac{\partial U}{\partial \phi} \boldsymbol{e}_\phi + \frac{\partial U}{\partial z} \boldsymbol{e}_z$

$\nabla \cdot \boldsymbol{F} = \frac{1}{r}\left[\frac{\partial(rF_r)}{\partial r} + \frac{\partial(rF_\phi)}{\partial \phi} + \frac{\partial(rF_z)}{\partial z} \right]$

$\nabla^2 U = \frac{1}{r}\frac{\partial}{\partial r}r \frac{\partial U}{\partial r} +\frac{1}{r^2} \frac{\partial^2 U}{\partial \phi^2} + \frac{\partial^2U}{\partial z^2}$

##### Spherical Polar Coordinates

$\nabla U = \frac{\partial U}{\partial r}\boldsymbol{e}_r + \frac{1}{r}\frac{\partial U}{\partial \theta} \boldsymbol{e}_\theta + \frac{1}{r \mathrm{sin}\theta}\frac{\partial U}{\partial \phi} \boldsymbol{e}_\phi$

$\nabla \cdot \boldsymbol{F} = \frac{1}{r^2 \mathrm{sin}\theta}\left[\frac{\partial (r^2 \mathrm{sin} \theta F_r)}{\partial r} + \frac{\partial (\mathrm{sin}\theta F_\theta)}{\partial \theta} + \frac{\partial(rF_\phi)}{\partial \phi} \right]$

$\nabla^2 U = \frac{1}{r^2}\frac{\partial}{\partial r}r^2 \frac{\partial(U)}{\partial r} + \frac{1}{\mathrm{sin}\theta} \frac{\partial}{\partial \theta}\left(\mathrm{sin}\theta \frac{\partial}{\partial \theta}\right) +\frac{1}{\mathrm{sin}^2\theta} \frac{\partial^2 U}{\partial \phi^2}$