Molecular dynamics implementation cheat sheet

A cheat-sheet for implementing a molecular dynamics simulation
Source: "Nanosystems: Molecular Machinery, Manufacturing, and Computation"

The MM2 model (3.3.2.)

[math]\displaystyle{ \mathcal{V}_{s} = \frac{1}{2}k_{s}(r - r_{0})^{2}\lbrack 1 - k_{cubic}(r - r_{0})\rbrack }[/math] … stretching (3.4)
[math]\displaystyle{ \mathcal{V}_{\theta} = \frac{1}{2}k_{\theta}(\theta - \theta_{0})^{2}\lbrack 1 + k_{sextic}(\theta - \theta_{0})^{4}\rbrack }[/math] … bending (3.5)
[math]\displaystyle{ k_{s\bot} = k_{\theta}/r_{0}^{2} }[/math] … per length angular stiffness [math]\displaystyle{ k_{s\bot} \approx k_{s}/20 }[/math] at [math]\displaystyle{ (r_{0},\theta_{0}) }[/math] (3.6)
[math]\displaystyle{ \mathcal{V}_{\omega} = \frac{1}{2}\lbrack V_{1}(1 - \cos(1\omega)) + V_{2}(1 - \cos(2\omega)) + V_{3}(1 + \cos(3\omega))\rbrack }[/math] … torsion (3.7)
[math]\displaystyle{ \mathcal{V}_{vdw} = \epsilon_{vdw}\lbrack 2.48 \times 10^{5}\exp( - 12.5\frac{r}{r_{vdw0}}) - 1.924(\frac{r}{r_{vdw0}})^{- 6}\rbrack }[/math] … vdW attraction & non-bonded Pauli repulsion (3.8)
[math]\displaystyle{ \mathcal{V}_{s\theta} = k_{s\theta}(\theta - \theta_{0})\lbrack(r_{A} - r_{A0}) + (r_{B} - r_{B0})\rbrack }[/math] … stretch-bend interaction (3.9)

---

[math]\displaystyle{ k_{s,C - C} = 440N/m }[/math] — [math]\displaystyle{ r_{0,C - C} = 111.3pm }[/math] … (Table 3.2.)
[math]\displaystyle{ k_{\theta,C - C - C} = 450zJ/rad^{2} }[/math] — [math]\displaystyle{ \theta_{0,C - C - C} = 1.911rad = 109.47{^\circ} }[/math] … (Table 3.3.)
[math]\displaystyle{ \epsilon_{vdw,Csp^{3}} = 357yJ }[/math] — [math]\displaystyle{ r_{vdw,Csp^{3}} = 190pm }[/math] … well depth [math]\displaystyle{ \epsilon_{vdw} }[/math] varies widely for other atoms (~10x) (Table 3.1.)
[math]\displaystyle{ V_{1,C - C - C - C} = 1.39zJ,\ V_{2,C - C - C - C} = 1.88zJ,\ V_{3,C - C - C - C} = 0.65zJ }[/math] … (Table 3.5.)
[math]\displaystyle{ k_{s\theta,C - C - C} = 1.2nN/rad }[/math] … (Table 3.6.)
[math]\displaystyle{ k_{cubic} = ?m^{- 1} }[/math] — [math]\displaystyle{ k_{sextic} = ?rad^{- 4} }[/math]

[math]\displaystyle{ n … nano = 10^{-9} }[/math] — [math]\displaystyle{ p … pico = 10^{-12} }[/math]
[math]\displaystyle{ z … zepto = 10^{-21} }[/math] — [math]\displaystyle{ y … yocto = 10^{-24} }[/math]

Bonds under large loads (beyond MM2) (3.3.3.)

[math]\displaystyle{ \mathcal{V}_{morse} = D_{e}({1 - \exp\lbrack - \beta(r - r_{0})\rbrack}^{2} - 1) }[/math] … tensile (3.10)
[math]\displaystyle{ \mathcal{V}_{lippincott} = D_{e}\lbrack 1 - \exp( - \frac{k_{s}r_{0}(r - r_{0})^{2}}{2D_{e}r})\rbrack,\ \ r \geq r_{0} }[/math] … tensile – more accurate for larger distances (3.15)
For compressive loads [math]\displaystyle{ \mathcal{V}_{vdw} }[/math] see (3.8) above [math]\displaystyle{ F_{vdw} = - \frac{\partial}{\partial r}\mathcal{V}_{vdw} }[/math] … (3.16)
Note: Popular Lennard-Jones 6-12 potential is too steep in repulsive regime!
MM2 potential within 10% of experiment for [math]\displaystyle{ 0.5r_{vdw} }[/math] ([math]\displaystyle{ > 100yJ }[/math] repulsion) 🙂

[math]\displaystyle{ \beta = \sqrt{k_{s}/(2D_{e})} }[/math] … (3.13)
[math]\displaystyle{ D_{e} \approx D_{0} + \frac{\hslash}{2}\sqrt{k_{s}/\mu} = D_{0} + \frac{\hslash\omega}{2};\ \ \mu = \frac{m_{1}m_{2}}{m_{1} + m_{2}} }[/math] … (3.14)
[math]\displaystyle{ D_{0} }[/math] … potential well depth
[math]\displaystyle{ D_{e,C - C} = 556yJ }[/math] — [math]\displaystyle{ k_{s,C - C} = 440N/m }[/math] — [math]\displaystyle{ r_{0,C - C} = 152.3pm }[/math] … (Table 3.8)

Verlet Integration (not in Nanosystems)

[math]\displaystyle{ \overset{\rightarrow}{x}(t + \Delta t) = 2\overset{\rightarrow}{x}(t) - \overset{\rightarrow}{x}(t - \Delta t) + \overset{\rightarrow}{a}(t)\Delta t^{2} + O(\Delta t^{4}) }[/math]
[math]\displaystyle{ \overset{\rightarrow}{a}(t) = m\nabla_{\overset{\rightarrow}{x}}\mathcal{V}(t) }[/math]

MM2 vs MM3 (Chapter 3.3.2.g.)

★ MM2 prioritizes accuracy in energy and geometry
★ MM3 prioritizes accuracy of vibration frequencies
★ MM3 predicts greater angle-bending stiffness by ~1.5x or more
★ MM3 has a more complex functional form – e.g. it adds:
– stretch-torsion interaction, cubic bending, quartic stretching
★ MM3 has lower energies forces and stiffness in deep repulsive regime by ~10%
★ MM2 overall seems like a more conservative (safe side wrong) choice for diamondoid nanomachinery

Bond cleavage & radical coupling (3.4.2)

The Morse potential [math]\displaystyle{ \mathcal{V}_{morse} }[/math] can approximate homolytic reactions
[math]\displaystyle{ \mathcal{V}_{anti - morse} = \frac{1}{2}D_{e}({1 + \exp\lbrack - \beta(r - r_{0})\rbrack}^{2} - 1) }[/math] … unpaired spins (3.21)
But unpaired spins being the limiting factor in reaction rate should
typically be avoidable in piezomechanosynthesis.
See Nanosystems Chapter 8.4.3.b. Radical coupling and intersystem crossing.

Nonbonded & large compression limits (3.3.3.b.)

[math]\displaystyle{ k_{s,vdw} \approx \frac{12.5}{r_{vdw0}}F_{vdw} \approx 3.5 \times 10^{10}m^{- 1} \cdot F_{vdw} }[/math] … (3.18)
[math]\displaystyle{ \mathcal{V}_{vdw} \approx 0.08r_{vdw0}F_{vdw} \approx 2.9 \times 10^{- 11}m \cdot F_{vdw} }[/math] …(3.19)

Abstraction reactions (Chapter 3.4.3.)

[math]\displaystyle{ \mathcal{V}_{LEPS} }[/math] extended London-Eyring-Polyano-Sato potential – details omitted here

Continuum models of van der Waals attraction (Chapter 3.5.1.)

Hamaker constant – details omitted here

Related

• Snapback – sideward bond bending stiffness is only 1/20th of radial bond stretching/compressing stiffness
• Atomic orbitals