| 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). |
| 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 |
|
*/ |
| 67 |
|
using namespace std; |
| 68 |
|
namespace OpenMD { |
| 69 |
|
|
| 70 |
< |
ForceManager::ForceManager(SimInfo * info) : info_(info) { |
| 70 |
> |
ForceManager::ForceManager(SimInfo * info) : info_(info), switcher_(NULL), |
| 71 |
> |
initialized_(false) { |
| 72 |
|
forceField_ = info_->getForceField(); |
| 73 |
|
interactionMan_ = new InteractionManager(); |
| 74 |
|
fDecomp_ = new ForceMatrixDecomposition(info_, interactionMan_); |
| 75 |
+ |
thermo = new Thermo(info_); |
| 76 |
|
} |
| 77 |
|
|
| 78 |
+ |
ForceManager::~ForceManager() { |
| 79 |
+ |
perturbations_.clear(); |
| 80 |
+ |
|
| 81 |
+ |
delete switcher_; |
| 82 |
+ |
delete interactionMan_; |
| 83 |
+ |
delete fDecomp_; |
| 84 |
+ |
delete thermo; |
| 85 |
+ |
} |
| 86 |
+ |
|
| 87 |
|
/** |
| 88 |
|
* setupCutoffs |
| 89 |
|
* |
| 98 |
|
* simulation for suggested cutoff values (e.g. 2.5 * sigma). |
| 99 |
|
* Use the maximum suggested value that was found. |
| 100 |
|
* |
| 101 |
< |
* cutoffMethod : (one of HARD, SWITCHED, SHIFTED_FORCE, |
| 101 |
> |
* cutoffMethod : (one of HARD, SWITCHED, SHIFTED_FORCE, TAYLOR_SHIFTED, |
| 102 |
|
* or SHIFTED_POTENTIAL) |
| 103 |
|
* If cutoffMethod was explicitly set, use that choice. |
| 104 |
|
* If cutoffMethod was not explicitly set, use SHIFTED_FORCE |
| 171 |
|
stringToCutoffMethod["SWITCHED"] = SWITCHED; |
| 172 |
|
stringToCutoffMethod["SHIFTED_POTENTIAL"] = SHIFTED_POTENTIAL; |
| 173 |
|
stringToCutoffMethod["SHIFTED_FORCE"] = SHIFTED_FORCE; |
| 174 |
+ |
stringToCutoffMethod["TAYLOR_SHIFTED"] = TAYLOR_SHIFTED; |
| 175 |
|
|
| 176 |
|
if (simParams_->haveCutoffMethod()) { |
| 177 |
|
string cutMeth = toUpperCopy(simParams_->getCutoffMethod()); |
| 181 |
|
sprintf(painCave.errMsg, |
| 182 |
|
"ForceManager::setupCutoffs: Could not find chosen cutoffMethod %s\n" |
| 183 |
|
"\tShould be one of: " |
| 184 |
< |
"HARD, SWITCHED, SHIFTED_POTENTIAL, or SHIFTED_FORCE\n", |
| 184 |
> |
"HARD, SWITCHED, SHIFTED_POTENTIAL, TAYLOR_SHIFTED,\n" |
| 185 |
> |
"\tor SHIFTED_FORCE\n", |
| 186 |
|
cutMeth.c_str()); |
| 187 |
|
painCave.isFatal = 1; |
| 188 |
|
painCave.severity = OPENMD_ERROR; |
| 226 |
|
cutoffMethod_ = SHIFTED_POTENTIAL; |
| 227 |
|
} else if (myMethod == "SHIFTED_FORCE") { |
| 228 |
|
cutoffMethod_ = SHIFTED_FORCE; |
| 229 |
+ |
} else if (myMethod == "TAYLOR_SHIFTED") { |
| 230 |
+ |
cutoffMethod_ = TAYLOR_SHIFTED; |
| 231 |
|
} |
| 232 |
|
|
| 233 |
|
if (simParams_->haveSwitchingRadius()) |
| 234 |
|
rSwitch_ = simParams_->getSwitchingRadius(); |
| 235 |
|
|
| 236 |
< |
if (myMethod == "SHIFTED_POTENTIAL" || myMethod == "SHIFTED_FORCE") { |
| 236 |
> |
if (myMethod == "SHIFTED_POTENTIAL" || myMethod == "SHIFTED_FORCE" || |
| 237 |
> |
myMethod == "TAYLOR_SHIFTED") { |
| 238 |
|
if (simParams_->haveSwitchingRadius()){ |
| 239 |
|
sprintf(painCave.errMsg, |
| 240 |
|
"ForceManager::setupCutoffs : DEPRECATED ERROR MESSAGE\n" |
| 449 |
|
perturbations_.push_back(eField); |
| 450 |
|
} |
| 451 |
|
|
| 452 |
+ |
usePeriodicBoundaryConditions_ = info_->getSimParams()->getUsePeriodicBoundaryConditions(); |
| 453 |
+ |
|
| 454 |
|
fDecomp_->distributeInitialData(); |
| 455 |
< |
|
| 455 |
> |
|
| 456 |
|
initialized_ = true; |
| 457 |
< |
|
| 457 |
> |
|
| 458 |
|
} |
| 459 |
< |
|
| 459 |
> |
|
| 460 |
|
void ForceManager::calcForces() { |
| 461 |
|
|
| 462 |
|
if (!initialized_) initialize(); |
| 463 |
< |
|
| 463 |
> |
|
| 464 |
|
preCalculation(); |
| 465 |
|
shortRangeInteractions(); |
| 466 |
|
longRangeInteractions(); |
| 659 |
|
|
| 660 |
|
void ForceManager::longRangeInteractions() { |
| 661 |
|
|
| 644 |
– |
|
| 662 |
|
Snapshot* curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot(); |
| 663 |
|
DataStorage* config = &(curSnapshot->atomData); |
| 664 |
|
DataStorage* cgConfig = &(curSnapshot->cgData); |
| 698 |
|
RealType rCutSq; |
| 699 |
|
bool in_switching_region; |
| 700 |
|
RealType sw, dswdr, swderiv; |
| 701 |
< |
vector<int> atomListColumn, atomListRow, atomListLocal; |
| 701 |
> |
vector<int> atomListColumn, atomListRow; |
| 702 |
|
InteractionData idat; |
| 703 |
|
SelfData sdat; |
| 704 |
|
RealType mf; |
| 728 |
|
idat.f1 = &f1; |
| 729 |
|
idat.sw = &sw; |
| 730 |
|
idat.shiftedPot = (cutoffMethod_ == SHIFTED_POTENTIAL) ? true : false; |
| 731 |
< |
idat.shiftedForce = (cutoffMethod_ == SHIFTED_FORCE) ? true : false; |
| 731 |
> |
idat.shiftedForce = (cutoffMethod_ == SHIFTED_FORCE || cutoffMethod_ == TAYLOR_SHIFTED) ? true : false; |
| 732 |
|
idat.doParticlePot = doParticlePot_; |
| 733 |
|
idat.doElectricField = doElectricField_; |
| 734 |
|
sdat.doParticlePot = doParticlePot_; |
| 743 |
|
|
| 744 |
|
if (iLoop == loopStart) { |
| 745 |
|
bool update_nlist = fDecomp_->checkNeighborList(); |
| 746 |
< |
if (update_nlist) |
| 746 |
> |
if (update_nlist) { |
| 747 |
> |
if (!usePeriodicBoundaryConditions_) |
| 748 |
> |
Mat3x3d bbox = thermo->getBoundingBox(); |
| 749 |
|
neighborList = fDecomp_->buildNeighborList(); |
| 750 |
< |
} |
| 750 |
> |
} |
| 751 |
> |
} |
| 752 |
|
|
| 753 |
|
for (vector<pair<int, int> >::iterator it = neighborList.begin(); |
| 754 |
|
it != neighborList.end(); ++it) { |
| 768 |
|
idat.rcut = &cuts.first; |
| 769 |
|
if (iLoop == PAIR_LOOP) { |
| 770 |
|
vij = 0.0; |
| 771 |
< |
fij = V3Zero; |
| 772 |
< |
eField1 = V3Zero; |
| 773 |
< |
eField2 = V3Zero; |
| 771 |
> |
fij.zero(); |
| 772 |
> |
eField1.zero(); |
| 773 |
> |
eField2.zero(); |
| 774 |
|
} |
| 775 |
|
|
| 776 |
|
in_switching_region = switcher_->getSwitch(rgrpsq, sw, dswdr, |
| 795 |
|
vpair = 0.0; |
| 796 |
|
workPot = 0.0; |
| 797 |
|
exPot = 0.0; |
| 798 |
< |
f1 = V3Zero; |
| 798 |
> |
f1.zero(); |
| 799 |
|
dVdFQ1 = 0.0; |
| 800 |
|
dVdFQ2 = 0.0; |
| 801 |
|
|