57 |
|
#include "primitives/Torsion.hpp" |
58 |
|
#include "primitives/Inversion.hpp" |
59 |
|
#include "nonbonded/NonBondedInteraction.hpp" |
60 |
< |
#include "perturbations/ElectricField.hpp" |
60 |
> |
#include "perturbations/UniformField.hpp" |
61 |
|
#include "parallel/ForceMatrixDecomposition.hpp" |
62 |
|
|
63 |
|
#include <cstdio> |
405 |
|
interactionMan_->setSimInfo(info_); |
406 |
|
interactionMan_->initialize(); |
407 |
|
|
408 |
< |
// We want to delay the cutoffs until after the interaction |
409 |
< |
// manager has set up the atom-atom interactions so that we can |
410 |
< |
// query them for suggested cutoff values |
408 |
> |
//! We want to delay the cutoffs until after the interaction |
409 |
> |
//! manager has set up the atom-atom interactions so that we can |
410 |
> |
//! query them for suggested cutoff values |
411 |
|
setupCutoffs(); |
412 |
|
|
413 |
|
info_->prepareTopology(); |
417 |
|
if (doHeatFlux_) doParticlePot_ = true; |
418 |
|
|
419 |
|
doElectricField_ = info_->getSimParams()->getOutputElectricField(); |
420 |
+ |
doSitePotential_ = info_->getSimParams()->getOutputSitePotential(); |
421 |
|
|
422 |
|
} |
423 |
|
|
424 |
|
ForceFieldOptions& fopts = forceField_->getForceFieldOptions(); |
425 |
|
|
426 |
< |
// Force fields can set options on how to scale van der Waals and |
427 |
< |
// electrostatic interactions for atoms connected via bonds, bends |
428 |
< |
// and torsions in this case the topological distance between |
429 |
< |
// atoms is: |
430 |
< |
// 0 = topologically unconnected |
431 |
< |
// 1 = bonded together |
432 |
< |
// 2 = connected via a bend |
433 |
< |
// 3 = connected via a torsion |
426 |
> |
//! Force fields can set options on how to scale van der Waals and |
427 |
> |
//! electrostatic interactions for atoms connected via bonds, bends |
428 |
> |
//! and torsions in this case the topological distance between |
429 |
> |
//! atoms is: |
430 |
> |
//! 0 = topologically unconnected |
431 |
> |
//! 1 = bonded together |
432 |
> |
//! 2 = connected via a bend |
433 |
> |
//! 3 = connected via a torsion |
434 |
|
|
435 |
|
vdwScale_.reserve(4); |
436 |
|
fill(vdwScale_.begin(), vdwScale_.end(), 0.0); |
448 |
|
electrostaticScale_[2] = fopts.getelectrostatic13scale(); |
449 |
|
electrostaticScale_[3] = fopts.getelectrostatic14scale(); |
450 |
|
|
451 |
< |
if (info_->getSimParams()->haveElectricField()) { |
452 |
< |
ElectricField* eField = new ElectricField(info_); |
451 |
> |
if (info_->getSimParams()->haveUniformField()) { |
452 |
> |
UniformField* eField = new UniformField(info_); |
453 |
|
perturbations_.push_back(eField); |
454 |
|
} |
455 |
|
|
638 |
|
// Collect from all nodes. This should eventually be moved into a |
639 |
|
// SystemDecomposition, but this is a better place than in |
640 |
|
// Thermo to do the collection. |
641 |
< |
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &bondPotential, 1, MPI::REALTYPE, |
642 |
< |
MPI::SUM); |
643 |
< |
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &bendPotential, 1, MPI::REALTYPE, |
644 |
< |
MPI::SUM); |
645 |
< |
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &torsionPotential, 1, |
646 |
< |
MPI::REALTYPE, MPI::SUM); |
647 |
< |
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &inversionPotential, 1, |
648 |
< |
MPI::REALTYPE, MPI::SUM); |
641 |
> |
|
642 |
> |
MPI_Allreduce(MPI_IN_PLACE, &bondPotential, 1, MPI_REALTYPE, |
643 |
> |
MPI_SUM, MPI_COMM_WORLD); |
644 |
> |
MPI_Allreduce(MPI_IN_PLACE, &bendPotential, 1, MPI_REALTYPE, |
645 |
> |
MPI_SUM, MPI_COMM_WORLD); |
646 |
> |
MPI_Allreduce(MPI_IN_PLACE, &torsionPotential, 1, |
647 |
> |
MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD); |
648 |
> |
MPI_Allreduce(MPI_IN_PLACE, &inversionPotential, 1, |
649 |
> |
MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD); |
650 |
|
#endif |
651 |
|
|
652 |
|
Snapshot* curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot(); |
717 |
|
potVec exPot(0.0); |
718 |
|
Vector3d eField1(0.0); |
719 |
|
Vector3d eField2(0.0); |
720 |
+ |
RealType sPot1(0.0); |
721 |
+ |
RealType sPot2(0.0); |
722 |
+ |
|
723 |
|
vector<int>::iterator ia, jb; |
724 |
|
|
725 |
|
int loopStart, loopEnd; |
735 |
|
idat.dVdFQ1 = &dVdFQ1; |
736 |
|
idat.dVdFQ2 = &dVdFQ2; |
737 |
|
idat.eField1 = &eField1; |
738 |
< |
idat.eField2 = &eField2; |
738 |
> |
idat.eField2 = &eField2; |
739 |
> |
idat.sPot1 = &sPot1; |
740 |
> |
idat.sPot2 = &sPot2; |
741 |
|
idat.f1 = &f1; |
742 |
|
idat.sw = &sw; |
743 |
|
idat.shiftedPot = (cutoffMethod_ == SHIFTED_POTENTIAL) ? true : false; |
744 |
|
idat.shiftedForce = (cutoffMethod_ == SHIFTED_FORCE || cutoffMethod_ == TAYLOR_SHIFTED) ? true : false; |
745 |
|
idat.doParticlePot = doParticlePot_; |
746 |
|
idat.doElectricField = doElectricField_; |
747 |
+ |
idat.doSitePotential = doSitePotential_; |
748 |
|
sdat.doParticlePot = doParticlePot_; |
749 |
|
|
750 |
|
loopEnd = PAIR_LOOP; |
784 |
|
fij.zero(); |
785 |
|
eField1.zero(); |
786 |
|
eField2.zero(); |
787 |
+ |
sPot1 = 0.0; |
788 |
+ |
sPot2 = 0.0; |
789 |
|
} |
790 |
|
|
791 |
|
in_switching_region = switcher_->getSwitch(rgrpsq, sw, dswdr, |
987 |
|
} |
988 |
|
|
989 |
|
#ifdef IS_MPI |
990 |
< |
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, stressTensor.getArrayPointer(), 9, |
991 |
< |
MPI::REALTYPE, MPI::SUM); |
990 |
> |
MPI_Allreduce(MPI_IN_PLACE, stressTensor.getArrayPointer(), 9, |
991 |
> |
MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD); |
992 |
|
#endif |
993 |
|
curSnapshot->setStressTensor(stressTensor); |
994 |
|
|