1 |
gezelter |
1743 |
#ifndef TYPES_MORSEBONDTYPE_HPP |
2 |
|
|
#define TYPES_MORSEBONDTYPE_HPP |
3 |
|
|
|
4 |
|
|
#include "types/BondType.hpp" |
5 |
|
|
|
6 |
|
|
namespace oopse { |
7 |
|
|
/** |
8 |
|
|
* @class MorseBondType |
9 |
|
|
* |
10 |
|
|
* MorseBondType is a more realistic bond potential. |
11 |
|
|
* |
12 |
|
|
* The functional form is given by: \f$V(r) = D_e (1 - e^{-\beta (r |
13 |
|
|
* - r_0)})^2\f$ where \f$D_e\f$ is the bond dissociation energy (in |
14 |
|
|
* kcal / mol), \f$\beta$\f is an inverse distance parameter related |
15 |
|
|
* to the force constant. \f$\beta = \sqrt{\frac{k}{2 D_e}}\f$, and |
16 |
|
|
* \f$r_0\f$ is the equilibrium bond length. |
17 |
|
|
*/ |
18 |
|
|
class MorseBondType : public BondType { |
19 |
|
|
|
20 |
|
|
public: |
21 |
|
|
|
22 |
tim |
1746 |
MorseBondType( double myR0, double myD, double myBeta) |
23 |
|
|
: BondType(myR0), De(myD), beta(myBeta) { |
24 |
gezelter |
1743 |
} |
25 |
|
|
|
26 |
|
|
void setWellDepth(double myD) { De = myD;} |
27 |
tim |
1746 |
|
28 |
gezelter |
1743 |
void setBeta(double myBeta) { beta = myBeta; } |
29 |
tim |
1746 |
|
30 |
gezelter |
1743 |
void setWellDepthAndForceConstant(double myD, double myK) { |
31 |
|
|
De = myD; |
32 |
|
|
beta = sqrt(myK/(2.0*De)); |
33 |
|
|
} |
34 |
|
|
|
35 |
|
|
double getWellDepth() {return De;} |
36 |
tim |
1746 |
|
37 |
gezelter |
1743 |
double getBeta() {return beta;} |
38 |
tim |
1746 |
|
39 |
gezelter |
1743 |
double getForceConstant() {return 2.0*De*beta*beta;} |
40 |
|
|
|
41 |
tim |
1748 |
virtual void calcForce(double r, double& V, double& dVdr) { |
42 |
gezelter |
1743 |
double dr, eterm, eterm2; |
43 |
|
|
|
44 |
|
|
dr = r - r0; |
45 |
|
|
eterm = exp(-beta*dr); |
46 |
|
|
eterm2 = eterm*eterm; |
47 |
|
|
|
48 |
|
|
V = De*(1 - 2.0*eterm + eterm2); |
49 |
|
|
dVdr = 2.0 * De * beta * (eterm - eterm2); |
50 |
|
|
} |
51 |
|
|
|
52 |
|
|
private: |
53 |
|
|
|
54 |
|
|
double De; |
55 |
|
|
double beta; |
56 |
|
|
|
57 |
|
|
}; |
58 |
|
|
} |
59 |
|
|
#endif |