ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/primitives/Inversion.cpp
(Generate patch)

Comparing trunk/src/primitives/Inversion.cpp (file contents):
Revision 1390 by gezelter, Wed Nov 25 20:02:06 2009 UTC vs.
Revision 1890 by gezelter, Tue Jun 18 21:06:25 2013 UTC

# Line 35 | Line 35
35   *                                                                      
36   * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).            
37   * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
38 < * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
39 < * [4]  Vardeman & Gezelter, in progress (2009).                        
38 > * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
39 > * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 > * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41   */
42  
43 + #include "config.h"
44 + #include <cmath>
45 +
46   #include "primitives/Inversion.hpp"
47  
48   namespace OpenMD {
# Line 46 | Line 50 | namespace OpenMD {
50    Inversion::Inversion(Atom *atom1, Atom *atom2, Atom *atom3,
51                         Atom *atom4, InversionType *it) :
52      atom1_(atom1), atom2_(atom2), atom3_(atom3), atom4_(atom4),
53 <    inversionType_(it) { }
53 >    inversionType_(it) {
54 >    inversionKey_ = inversionType_->getKey();
55 >  }
56    
57 <  void Inversion::calcForce(RealType& angle) {
57 >  void Inversion::calcForce(RealType& angle, bool doParticlePot) {
58      
59      // In OpenMD's version of an inversion, the central atom
60      // comes first.  However, to get the planarity in a typical cosine
# Line 82 | Line 88 | namespace OpenMD {
88      if (cos_phi < -1.0) cos_phi = -1.0;
89  
90      RealType dVdcosPhi;
91 <    inversionType_->calcForce(cos_phi, potential_, dVdcosPhi);
91 >    switch (inversionKey_) {
92 >    case itCosAngle:
93 >      inversionType_->calcForce(cos_phi, potential_, dVdcosPhi);
94 >      break;
95 >    case itAngle:
96 >      RealType phi = acos(cos_phi);
97 >      RealType dVdPhi;
98 >      inversionType_->calcForce(phi, potential_, dVdPhi);
99 >      RealType sin_phi = sqrt(1.0 - cos_phi * cos_phi);  
100 >      if (fabs(sin_phi) < 1.0E-6) {
101 >        sin_phi = 1.0E-6;
102 >      }
103 >      dVdcosPhi = dVdPhi / sin_phi;
104 >      break;
105 >    }
106 >    
107      Vector3d f1 ;
108      Vector3d f2 ;
109      Vector3d f3 ;
# Line 109 | Line 130 | namespace OpenMD {
130      atom4_->addFrc(-f2);
131      atom3_->addFrc(-f3);
132  
133 <    atom1_->addParticlePot(potential_);
134 <    atom2_->addParticlePot(potential_);
135 <    atom3_->addParticlePot(potential_);
136 <    atom4_->addParticlePot(potential_);
137 <
133 >    if (doParticlePot) {
134 >      atom1_->addParticlePot(potential_);
135 >      atom2_->addParticlePot(potential_);
136 >      atom3_->addParticlePot(potential_);
137 >      atom4_->addParticlePot(potential_);
138 >    }
139 >    
140      angle = acos(cos_phi) /M_PI * 180.0;
141    }
142  

Comparing trunk/src/primitives/Inversion.cpp (property svn:keywords):
Revision 1390 by gezelter, Wed Nov 25 20:02:06 2009 UTC vs.
Revision 1890 by gezelter, Tue Jun 18 21:06:25 2013 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines