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

Comparing trunk/src/restraints/ThermoIntegrationForceManager.cpp (file contents):
Revision 417 by chrisfen, Thu Mar 10 15:10:24 2005 UTC vs.
Revision 665 by tim, Thu Oct 13 22:26:47 2005 UTC

# Line 1 | Line 1
1 < /*
1 > /*
2   * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   *
4   * The University of Notre Dame grants you ("Licensee") a
# Line 42 | Line 42
42   #include <cmath>
43   #include "restraints/ThermoIntegrationForceManager.hpp"
44   #include "integrators/Integrator.hpp"
45 + #include "math/SquareMatrix3.hpp"
46   #include "primitives/Molecule.hpp"
47   #include "utils/simError.h"
48   #include "utils/OOPSEConstant.hpp"
# Line 50 | Line 51 | namespace oopse {
51   namespace oopse {
52    
53    ThermoIntegrationForceManager::ThermoIntegrationForceManager(SimInfo* info):
54 <  ForceManager(info){
55 <    currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
56 <    simParam = info_->getSimParams();
54 >    ForceManager(info){
55 >      currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
56 >      simParam = info_->getSimParams();
57      
58 <    if (simParam->haveThermIntLambda()){
59 <      tIntLambda_ = simParam->getThermIntLambda();
60 <    }
61 <    else{
62 <      tIntLambda_ = 1.0;
63 <      sprintf(painCave.errMsg,
64 <              "ThermoIntegration error: the transformation parameter (lambda) was\n"
65 <              "\tnot specified. OOPSE will use a default value of %f. To set\n"
66 <              "\tlambda, use the thermodynamicIntegrationLambda variable.\n",
67 <              tIntLambda_);
68 <      painCave.isFatal = 0;
69 <      simError();
70 <    }
58 >      if (simParam->haveThermodynamicIntegrationLambda()){
59 >        tIntLambda_ = simParam->getThermodynamicIntegrationLambda();
60 >      }
61 >      else{
62 >        tIntLambda_ = 1.0;
63 >        sprintf(painCave.errMsg,
64 >                "ThermoIntegration error: the transformation parameter (lambda) was\n"
65 >                "\tnot specified. OOPSE will use a default value of %f. To set\n"
66 >                "\tlambda, use the thermodynamicIntegrationLambda variable.\n",
67 >                tIntLambda_);
68 >        painCave.isFatal = 0;
69 >        simError();
70 >      }
71      
72 <    if (simParam->haveThermIntK()){
73 <      tIntK_ = simParam->getThermIntK();
74 <    }
75 <    else{
76 <      tIntK_ = 1.0;
77 <      sprintf(painCave.errMsg,
78 <              "ThermoIntegration Warning: the tranformation parameter exponent\n"
79 <              "\t(k) was not specified. OOPSE will use a default value of %f.\n"
80 <              "\tTo set k, use the thermodynamicIntegrationK variable.\n",
81 <              tIntK_);
82 <      painCave.isFatal = 0;
83 <      simError();      
84 <    }
72 >      if (simParam->haveThermodynamicIntegrationK()){
73 >        tIntK_ = simParam->getThermodynamicIntegrationK();
74 >      }
75 >      else{
76 >        tIntK_ = 1.0;
77 >        sprintf(painCave.errMsg,
78 >                "ThermoIntegration Warning: the tranformation parameter exponent\n"
79 >                "\t(k) was not specified. OOPSE will use a default value of %f.\n"
80 >                "\tTo set k, use the thermodynamicIntegrationK variable.\n",
81 >                tIntK_);
82 >        painCave.isFatal = 0;
83 >        simError();      
84 >      }
85      
86 <    if (simParam->getUseSolidThermInt()) {
87 <      // build a restraint object
88 <      restraint_ =  new Restraints(info_, tIntLambda_, tIntK_);
86 >      if (simParam->getUseSolidThermInt()) {
87 >        // build a restraint object
88 >        restraint_ =  new Restraints(info_, tIntLambda_, tIntK_);
89        
90 <    }
90 >      }
91      
92 <    // build the scaling factor used to modulate the forces and torques
93 <    factor_ = pow(tIntLambda_, tIntK_);
94 <    
95 <  }
92 >      // build the scaling factor used to modulate the forces and torques
93 >      factor_ = pow(tIntLambda_, tIntK_);
94 >
95 >    }
96    
97    ThermoIntegrationForceManager::~ThermoIntegrationForceManager(){
98    }
# Line 105 | Line 106 | namespace oopse {
106      StuntDouble* integrableObject;
107      Vector3d frc;
108      Vector3d trq;
109 +    Mat3x3d tempTau;
110      
109 –    curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot();
110 –    
111      // perform the standard calcForces first
112      ForceManager::calcForces(needPotential, needStress);
113      
114 +    curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot();
115 +
116      // now scale forces and torques of all the integrableObjects
117        
118      for (mol = info_->beginMolecule(mi); mol != NULL;
# Line 128 | Line 130 | namespace oopse {
130            integrableObject->setTrq(trq);
131          }
132        }
131 –      
132 –      // set vraw to be the unmodulated potential
133 –      lrPot_ = curSnapshot->statData[Stats::POTENTIAL_ENERGY];
134 –      curSnapshot->statData[Stats::VRAW] = lrPot_;
135 –      
136 –      // modulate the potential and update the snapshot
137 –      lrPot_ *= factor_;
138 –      curSnapshot->statData[Stats::POTENTIAL_ENERGY] = lrPot_;
133      }
134 +  
135 +    // set vraw to be the unmodulated potential
136 +    lrPot_ = curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL];
137 +    curSnapshot->statData[Stats::VRAW] = lrPot_;
138      
139 +    // modulate the potential and update the snapshot
140 +    lrPot_ *= factor_;
141 +    curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] = lrPot_;
142 +    
143 +    // scale the pressure tensor
144 +    tempTau = curSnapshot->statData.getTau();
145 +    tempTau *= factor_;
146 +    curSnapshot->statData.setTau(tempTau);
147 +
148      // do crystal restraint forces for thermodynamic integration
149      if (simParam->getUseSolidThermInt()) {
150        
151        lrPot_ += restraint_->Calc_Restraint_Forces();
152 <      curSnapshot->statData[Stats::POTENTIAL_ENERGY] = lrPot_;
152 >      curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] = lrPot_;
153        
154        vHarm_ = restraint_->getVharm();
155        curSnapshot->statData[Stats::VHARM] = vHarm_;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines