1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
|
#ifndef CBLJCUT_H
#define CBLJCUT_H
#include <iostream>
#include <cstdlib>
#include "CbPotential.h"
namespace ATC
{
/**
* @class CbLjCut
* @brief Class for computing Cauchy-Born quantities for a Lennard-Jones/cut material
* (A factor of one-half is already included to split the
* bond energy between atoms)
*/
class CbLjCut : public CbPotential
{
public:
//! Constructor - initializes coefficients and enables PAIRWISE term.
CbLjCut(double eps, double sig, double cutoff_radius)
: CbPotential(Interactions(PAIRWISE)),
A (4.0*eps*pow(sig, 12)),
B (4.0*eps*pow(sig, 6)),
rcut(cutoff_radius)
{ }
//! Returns the cutoff readius of the LJ potential.
double cutoff_radius() const { return rcut; }
//! Returns the LJ energy between two interacting atoms (6,12).
double phi(const double &r) const
{
const double r6i = 1.0/((r*r*r)*(r*r*r));
return r6i*(A*r6i - B);
}
//! Returns the first derivative of the LJ energy (7,13).
double phi_r(const double &r) const
{
const double ri = 1.0/r;
const double r6i = (ri*ri*ri)*(ri*ri*ri);
return r6i*ri*(6.0*B - 12.0*A*r6i);
}
//! Returns the second derivative of the LJ energy (8,14).
double phi_rr(const double &r) const
{
const double r2i = 1.0/(r*r);
const double r6i = r2i*r2i*r2i;
return r6i*r2i*(13.0*12.0*A*r6i - 7.0*6.0*B);
}
//! Returns the third derivative of the LJ bond energy (9,15).
double phi_rrr(const double &r) const
{
const double r3i = 1.0/(r*r*r);
const double r9i = r3i*r3i*r3i;
return r9i*(8.0*7.0*6.0*B - 14.0*13.0*12.0*A*r3i*r3i);
}
const double A, B; //!< phi = Ar^-12 + Br^-6
const double rcut; //!< cutoff force
};
}
#endif
|
