| 53 |
|
#include "brains/SimInfo.hpp" |
| 54 |
|
#include "math/Vector3.hpp" |
| 55 |
|
#include "primitives/Molecule.hpp" |
| 56 |
+ |
#include "primitives/StuntDouble.hpp" |
| 57 |
|
#include "UseTheForce/fCutoffPolicy.h" |
| 58 |
|
#include "UseTheForce/DarkSide/fElectrostaticSummationMethod.h" |
| 59 |
|
#include "UseTheForce/DarkSide/fElectrostaticScreeningMethod.h" |
| 85 |
|
|
| 86 |
|
SimInfo::SimInfo(ForceField* ff, Globals* simParams) : |
| 87 |
|
forceField_(ff), simParams_(simParams), |
| 88 |
< |
ndf_(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), |
| 88 |
> |
ndf_(0), fdf_local(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), |
| 89 |
|
nGlobalMols_(0), nGlobalAtoms_(0), nGlobalCutoffGroups_(0), |
| 90 |
|
nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), |
| 91 |
|
nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nRigidBodies_(0), |
| 92 |
|
nIntegrableObjects_(0), nCutoffGroups_(0), nConstraints_(0), |
| 93 |
< |
sman_(NULL), fortranInitialized_(false) { |
| 93 |
> |
sman_(NULL), fortranInitialized_(false), calcBoxDipole_(false) { |
| 94 |
|
|
| 95 |
|
MoleculeStamp* molStamp; |
| 96 |
|
int nMolWithSameStamp; |
| 291 |
|
|
| 292 |
|
} |
| 293 |
|
|
| 294 |
+ |
int SimInfo::getFdf() { |
| 295 |
+ |
#ifdef IS_MPI |
| 296 |
+ |
MPI_Allreduce(&fdf_local,&fdf_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); |
| 297 |
+ |
#else |
| 298 |
+ |
fdf_ = fdf_local; |
| 299 |
+ |
#endif |
| 300 |
+ |
return fdf_; |
| 301 |
+ |
} |
| 302 |
+ |
|
| 303 |
|
void SimInfo::calcNdfRaw() { |
| 304 |
|
int ndfRaw_local; |
| 305 |
|
|
| 603 |
|
|
| 604 |
|
setupElectrostaticSummationMethod( isError ); |
| 605 |
|
setupSwitchingFunction(); |
| 606 |
+ |
setupAccumulateBoxDipole(); |
| 607 |
|
|
| 608 |
|
if(isError){ |
| 609 |
|
sprintf( painCave.errMsg, |
| 663 |
|
int usePBC = simParams_->getUsePeriodicBoundaryConditions(); |
| 664 |
|
int useRF; |
| 665 |
|
int useSF; |
| 666 |
+ |
int useSP; |
| 667 |
+ |
int useBoxDipole; |
| 668 |
|
std::string myMethod; |
| 669 |
|
|
| 670 |
|
// set the useRF logical |
| 675 |
|
if (simParams_->haveElectrostaticSummationMethod()) { |
| 676 |
|
std::string myMethod = simParams_->getElectrostaticSummationMethod(); |
| 677 |
|
toUpper(myMethod); |
| 678 |
< |
if (myMethod == "REACTION_FIELD") { |
| 678 |
> |
if (myMethod == "REACTION_FIELD"){ |
| 679 |
|
useRF=1; |
| 680 |
< |
} else { |
| 681 |
< |
if (myMethod == "SHIFTED_FORCE") { |
| 682 |
< |
useSF = 1; |
| 683 |
< |
} |
| 680 |
> |
} else if (myMethod == "SHIFTED_FORCE"){ |
| 681 |
> |
useSF = 1; |
| 682 |
> |
} else if (myMethod == "SHIFTED_POTENTIAL"){ |
| 683 |
> |
useSP = 1; |
| 684 |
|
} |
| 685 |
|
} |
| 686 |
+ |
|
| 687 |
+ |
if (simParams_->haveAccumulateBoxDipole()) |
| 688 |
+ |
if (simParams_->getAccumulateBoxDipole()) |
| 689 |
+ |
useBoxDipole = 1; |
| 690 |
|
|
| 691 |
|
//loop over all of the atom types |
| 692 |
|
for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { |
| 757 |
|
MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
| 758 |
|
|
| 759 |
|
temp = useSF; |
| 760 |
< |
MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
| 760 |
> |
MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
| 761 |
|
|
| 762 |
+ |
temp = useSP; |
| 763 |
+ |
MPI_Allreduce(&temp, &useSP, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
| 764 |
+ |
|
| 765 |
+ |
temp = useBoxDipole; |
| 766 |
+ |
MPI_Allreduce(&temp, &useBoxDipole, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
| 767 |
+ |
|
| 768 |
|
#endif |
| 769 |
|
|
| 770 |
|
fInfo_.SIM_uses_PBC = usePBC; |
| 782 |
|
fInfo_.SIM_uses_FLARB = useFLARB; |
| 783 |
|
fInfo_.SIM_uses_RF = useRF; |
| 784 |
|
fInfo_.SIM_uses_SF = useSF; |
| 785 |
+ |
fInfo_.SIM_uses_SP = useSP; |
| 786 |
+ |
fInfo_.SIM_uses_BoxDipole = useBoxDipole; |
| 787 |
|
|
| 788 |
|
if( myMethod == "REACTION_FIELD") { |
| 789 |
|
|
| 815 |
|
} |
| 816 |
|
|
| 817 |
|
//calculate mass ratio of cutoff group |
| 818 |
< |
std::vector<double> mfact; |
| 818 |
> |
std::vector<RealType> mfact; |
| 819 |
|
SimInfo::MoleculeIterator mi; |
| 820 |
|
Molecule* mol; |
| 821 |
|
Molecule::CutoffGroupIterator ci; |
| 822 |
|
CutoffGroup* cg; |
| 823 |
|
Molecule::AtomIterator ai; |
| 824 |
|
Atom* atom; |
| 825 |
< |
double totalMass; |
| 825 |
> |
RealType totalMass; |
| 826 |
|
|
| 827 |
|
//to avoid memory reallocation, reserve enough space for mfact |
| 828 |
|
mfact.reserve(getNCutoffGroups()); |
| 982 |
|
notifyFortranCutoffPolicy(&cp); |
| 983 |
|
|
| 984 |
|
// Check the Skin Thickness for neighborlists |
| 985 |
< |
double skin; |
| 985 |
> |
RealType skin; |
| 986 |
|
if (simParams_->haveSkinThickness()) { |
| 987 |
|
skin = simParams_->getSkinThickness(); |
| 988 |
|
notifyFortranSkinThickness(&skin); |
| 1072 |
|
int errorOut; |
| 1073 |
|
int esm = NONE; |
| 1074 |
|
int sm = UNDAMPED; |
| 1075 |
< |
double alphaVal; |
| 1076 |
< |
double dielectric; |
| 1075 |
> |
RealType alphaVal; |
| 1076 |
> |
RealType dielectric; |
| 1077 |
|
|
| 1078 |
|
errorOut = isError; |
| 1079 |
|
alphaVal = simParams_->getDampingAlpha(); |
| 1174 |
|
|
| 1175 |
|
// send switching function notification to switcheroo |
| 1176 |
|
setFunctionType(&ft); |
| 1177 |
+ |
|
| 1178 |
+ |
} |
| 1179 |
|
|
| 1180 |
+ |
void SimInfo::setupAccumulateBoxDipole() { |
| 1181 |
+ |
|
| 1182 |
+ |
// we only call setAccumulateBoxDipole if the accumulateBoxDipole parameter is true |
| 1183 |
+ |
if ( simParams_->haveAccumulateBoxDipole() ) |
| 1184 |
+ |
if ( simParams_->getAccumulateBoxDipole() ) { |
| 1185 |
+ |
setAccumulateBoxDipole(); |
| 1186 |
+ |
calcBoxDipole_ = true; |
| 1187 |
+ |
} |
| 1188 |
+ |
|
| 1189 |
|
} |
| 1190 |
|
|
| 1191 |
|
void SimInfo::addProperty(GenericData* genData) { |
| 1244 |
|
Molecule* mol; |
| 1245 |
|
|
| 1246 |
|
Vector3d comVel(0.0); |
| 1247 |
< |
double totalMass = 0.0; |
| 1247 |
> |
RealType totalMass = 0.0; |
| 1248 |
|
|
| 1249 |
|
|
| 1250 |
|
for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { |
| 1251 |
< |
double mass = mol->getMass(); |
| 1251 |
> |
RealType mass = mol->getMass(); |
| 1252 |
|
totalMass += mass; |
| 1253 |
|
comVel += mass * mol->getComVel(); |
| 1254 |
|
} |
| 1255 |
|
|
| 1256 |
|
#ifdef IS_MPI |
| 1257 |
< |
double tmpMass = totalMass; |
| 1257 |
> |
RealType tmpMass = totalMass; |
| 1258 |
|
Vector3d tmpComVel(comVel); |
| 1259 |
< |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1260 |
< |
MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1259 |
> |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1260 |
> |
MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1261 |
|
#endif |
| 1262 |
|
|
| 1263 |
|
comVel /= totalMass; |
| 1270 |
|
Molecule* mol; |
| 1271 |
|
|
| 1272 |
|
Vector3d com(0.0); |
| 1273 |
< |
double totalMass = 0.0; |
| 1273 |
> |
RealType totalMass = 0.0; |
| 1274 |
|
|
| 1275 |
|
for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { |
| 1276 |
< |
double mass = mol->getMass(); |
| 1276 |
> |
RealType mass = mol->getMass(); |
| 1277 |
|
totalMass += mass; |
| 1278 |
|
com += mass * mol->getCom(); |
| 1279 |
|
} |
| 1280 |
|
|
| 1281 |
|
#ifdef IS_MPI |
| 1282 |
< |
double tmpMass = totalMass; |
| 1282 |
> |
RealType tmpMass = totalMass; |
| 1283 |
|
Vector3d tmpCom(com); |
| 1284 |
< |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1285 |
< |
MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1284 |
> |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1285 |
> |
MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1286 |
|
#endif |
| 1287 |
|
|
| 1288 |
|
com /= totalMass; |
| 1306 |
|
Molecule* mol; |
| 1307 |
|
|
| 1308 |
|
|
| 1309 |
< |
double totalMass = 0.0; |
| 1309 |
> |
RealType totalMass = 0.0; |
| 1310 |
|
|
| 1311 |
|
|
| 1312 |
|
for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { |
| 1313 |
< |
double mass = mol->getMass(); |
| 1313 |
> |
RealType mass = mol->getMass(); |
| 1314 |
|
totalMass += mass; |
| 1315 |
|
com += mass * mol->getCom(); |
| 1316 |
|
comVel += mass * mol->getComVel(); |
| 1317 |
|
} |
| 1318 |
|
|
| 1319 |
|
#ifdef IS_MPI |
| 1320 |
< |
double tmpMass = totalMass; |
| 1320 |
> |
RealType tmpMass = totalMass; |
| 1321 |
|
Vector3d tmpCom(com); |
| 1322 |
|
Vector3d tmpComVel(comVel); |
| 1323 |
< |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1324 |
< |
MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1325 |
< |
MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1323 |
> |
MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1324 |
> |
MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1325 |
> |
MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1326 |
|
#endif |
| 1327 |
|
|
| 1328 |
|
com /= totalMass; |
| 1341 |
|
void SimInfo::getInertiaTensor(Mat3x3d &inertiaTensor, Vector3d &angularMomentum){ |
| 1342 |
|
|
| 1343 |
|
|
| 1344 |
< |
double xx = 0.0; |
| 1345 |
< |
double yy = 0.0; |
| 1346 |
< |
double zz = 0.0; |
| 1347 |
< |
double xy = 0.0; |
| 1348 |
< |
double xz = 0.0; |
| 1349 |
< |
double yz = 0.0; |
| 1344 |
> |
RealType xx = 0.0; |
| 1345 |
> |
RealType yy = 0.0; |
| 1346 |
> |
RealType zz = 0.0; |
| 1347 |
> |
RealType xy = 0.0; |
| 1348 |
> |
RealType xz = 0.0; |
| 1349 |
> |
RealType yz = 0.0; |
| 1350 |
|
Vector3d com(0.0); |
| 1351 |
|
Vector3d comVel(0.0); |
| 1352 |
|
|
| 1358 |
|
Vector3d thisq(0.0); |
| 1359 |
|
Vector3d thisv(0.0); |
| 1360 |
|
|
| 1361 |
< |
double thisMass = 0.0; |
| 1361 |
> |
RealType thisMass = 0.0; |
| 1362 |
|
|
| 1363 |
|
|
| 1364 |
|
|
| 1396 |
|
#ifdef IS_MPI |
| 1397 |
|
Mat3x3d tmpI(inertiaTensor); |
| 1398 |
|
Vector3d tmpAngMom; |
| 1399 |
< |
MPI_Allreduce(tmpI.getArrayPointer(), inertiaTensor.getArrayPointer(),9,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1400 |
< |
MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1399 |
> |
MPI_Allreduce(tmpI.getArrayPointer(), inertiaTensor.getArrayPointer(),9,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1400 |
> |
MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1401 |
|
#endif |
| 1402 |
|
|
| 1403 |
|
return; |
| 1418 |
|
Vector3d thisr(0.0); |
| 1419 |
|
Vector3d thisp(0.0); |
| 1420 |
|
|
| 1421 |
< |
double thisMass; |
| 1421 |
> |
RealType thisMass; |
| 1422 |
|
|
| 1423 |
|
for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { |
| 1424 |
|
thisMass = mol->getMass(); |
| 1431 |
|
|
| 1432 |
|
#ifdef IS_MPI |
| 1433 |
|
Vector3d tmpAngMom; |
| 1434 |
< |
MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); |
| 1434 |
> |
MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); |
| 1435 |
|
#endif |
| 1436 |
|
|
| 1437 |
|
return angularMomentum; |
| 1438 |
|
} |
| 1439 |
|
|
| 1440 |
< |
|
| 1440 |
> |
StuntDouble* SimInfo::getIOIndexToIntegrableObject(int index) { |
| 1441 |
> |
return IOIndexToIntegrableObject.at(index); |
| 1442 |
> |
} |
| 1443 |
> |
|
| 1444 |
> |
void SimInfo::setIOIndexToIntegrableObject(const std::vector<StuntDouble*>& v) { |
| 1445 |
> |
IOIndexToIntegrableObject= v; |
| 1446 |
> |
} |
| 1447 |
> |
|
| 1448 |
> |
/* |
| 1449 |
> |
void SimInfo::setStuntDoubleFromGlobalIndex(std::vector<StuntDouble*> v) { |
| 1450 |
> |
assert( v.size() == nAtoms_ + nRigidBodies_); |
| 1451 |
> |
sdByGlobalIndex_ = v; |
| 1452 |
> |
} |
| 1453 |
> |
|
| 1454 |
> |
StuntDouble* SimInfo::getStuntDoubleFromGlobalIndex(int index) { |
| 1455 |
> |
//assert(index < nAtoms_ + nRigidBodies_); |
| 1456 |
> |
return sdByGlobalIndex_.at(index); |
| 1457 |
> |
} |
| 1458 |
> |
*/ |
| 1459 |
|
}//end namespace oopse |
| 1460 |
|
|