| 84 |  |  | 
| 85 |  | SimInfo::SimInfo(ForceField* ff, Globals* simParams) : | 
| 86 |  | forceField_(ff), simParams_(simParams), | 
| 87 | < | ndf_(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), | 
| 87 | > | ndf_(0), fdf_local(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), | 
| 88 |  | nGlobalMols_(0), nGlobalAtoms_(0), nGlobalCutoffGroups_(0), | 
| 89 |  | nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), | 
| 90 |  | nAtoms_(0), nBonds_(0),  nBends_(0), nTorsions_(0), nRigidBodies_(0), | 
| 91 |  | nIntegrableObjects_(0),  nCutoffGroups_(0), nConstraints_(0), | 
| 92 | < | sman_(NULL), fortranInitialized_(false) { | 
| 92 | > | sman_(NULL), fortranInitialized_(false), calcBoxDipole_(false) { | 
| 93 |  |  | 
| 94 |  | MoleculeStamp* molStamp; | 
| 95 |  | int nMolWithSameStamp; | 
| 290 |  |  | 
| 291 |  | } | 
| 292 |  |  | 
| 293 | + | int SimInfo::getFdf() { | 
| 294 | + | #ifdef IS_MPI | 
| 295 | + | MPI_Allreduce(&fdf_local,&fdf_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); | 
| 296 | + | #else | 
| 297 | + | fdf_ = fdf_local; | 
| 298 | + | #endif | 
| 299 | + | return fdf_; | 
| 300 | + | } | 
| 301 | + |  | 
| 302 |  | void SimInfo::calcNdfRaw() { | 
| 303 |  | int ndfRaw_local; | 
| 304 |  |  | 
| 602 |  |  | 
| 603 |  | setupElectrostaticSummationMethod( isError ); | 
| 604 |  | setupSwitchingFunction(); | 
| 605 | + | setupAccumulateBoxDipole(); | 
| 606 |  |  | 
| 607 |  | if(isError){ | 
| 608 |  | sprintf( painCave.errMsg, | 
| 662 |  | int usePBC = simParams_->getUsePeriodicBoundaryConditions(); | 
| 663 |  | int useRF; | 
| 664 |  | int useSF; | 
| 665 | + | int useSP; | 
| 666 | + | int useBoxDipole; | 
| 667 |  | std::string myMethod; | 
| 668 |  |  | 
| 669 |  | // set the useRF logical | 
| 674 |  | if (simParams_->haveElectrostaticSummationMethod()) { | 
| 675 |  | std::string myMethod = simParams_->getElectrostaticSummationMethod(); | 
| 676 |  | toUpper(myMethod); | 
| 677 | < | if (myMethod == "REACTION_FIELD") { | 
| 677 | > | if (myMethod == "REACTION_FIELD"){ | 
| 678 |  | useRF=1; | 
| 679 | < | } else { | 
| 680 | < | if (myMethod == "SHIFTED_FORCE") { | 
| 681 | < | useSF = 1; | 
| 682 | < | } | 
| 679 | > | } else if (myMethod == "SHIFTED_FORCE"){ | 
| 680 | > | useSF = 1; | 
| 681 | > | } else if (myMethod == "SHIFTED_POTENTIAL"){ | 
| 682 | > | useSP = 1; | 
| 683 |  | } | 
| 684 |  | } | 
| 685 | + |  | 
| 686 | + | if (simParams_->haveAccumulateBoxDipole()) | 
| 687 | + | if (simParams_->getAccumulateBoxDipole()) | 
| 688 | + | useBoxDipole = 1; | 
| 689 |  |  | 
| 690 |  | //loop over all of the atom types | 
| 691 |  | for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { | 
| 756 |  | MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 757 |  |  | 
| 758 |  | temp = useSF; | 
| 759 | < | MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 759 | > | MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 760 |  |  | 
| 761 | + | temp = useSP; | 
| 762 | + | MPI_Allreduce(&temp, &useSP, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 763 | + |  | 
| 764 | + | temp = useBoxDipole; | 
| 765 | + | MPI_Allreduce(&temp, &useBoxDipole, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 766 | + |  | 
| 767 |  | #endif | 
| 768 |  |  | 
| 769 |  | fInfo_.SIM_uses_PBC = usePBC; | 
| 781 |  | fInfo_.SIM_uses_FLARB = useFLARB; | 
| 782 |  | fInfo_.SIM_uses_RF = useRF; | 
| 783 |  | fInfo_.SIM_uses_SF = useSF; | 
| 784 | + | fInfo_.SIM_uses_SP = useSP; | 
| 785 | + | fInfo_.SIM_uses_BoxDipole = useBoxDipole; | 
| 786 |  |  | 
| 787 |  | if( myMethod == "REACTION_FIELD") { | 
| 788 |  |  | 
| 814 |  | } | 
| 815 |  |  | 
| 816 |  | //calculate mass ratio of cutoff group | 
| 817 | < | std::vector<double> mfact; | 
| 817 | > | std::vector<RealType> mfact; | 
| 818 |  | SimInfo::MoleculeIterator mi; | 
| 819 |  | Molecule* mol; | 
| 820 |  | Molecule::CutoffGroupIterator ci; | 
| 821 |  | CutoffGroup* cg; | 
| 822 |  | Molecule::AtomIterator ai; | 
| 823 |  | Atom* atom; | 
| 824 | < | double totalMass; | 
| 824 | > | RealType totalMass; | 
| 825 |  |  | 
| 826 |  | //to avoid memory reallocation, reserve enough space for mfact | 
| 827 |  | mfact.reserve(getNCutoffGroups()); | 
| 981 |  | notifyFortranCutoffPolicy(&cp); | 
| 982 |  |  | 
| 983 |  | // Check the Skin Thickness for neighborlists | 
| 984 | < | double skin; | 
| 984 | > | RealType skin; | 
| 985 |  | if (simParams_->haveSkinThickness()) { | 
| 986 |  | skin = simParams_->getSkinThickness(); | 
| 987 |  | notifyFortranSkinThickness(&skin); | 
| 1000 |  | rsw_ = 0.85 * rcut_; | 
| 1001 |  | sprintf(painCave.errMsg, | 
| 1002 |  | "SimCreator Warning: No value was set for the switchingRadius.\n" | 
| 1003 | < | "\tOOPSE will use a default value of 85\% of the cutoffRadius" | 
| 1003 | > | "\tOOPSE will use a default value of 85 percent of the cutoffRadius.\n" | 
| 1004 |  | "\tswitchingRadius = %f. for this simulation\n", rsw_); | 
| 1005 |  | painCave.isFatal = 0; | 
| 1006 |  | simError(); | 
| 1013 |  | painCave.isFatal = 0; | 
| 1014 |  | simError(); | 
| 1015 |  | } | 
| 1016 | < |  | 
| 1016 | > | } | 
| 1017 | > |  | 
| 1018 |  | notifyFortranCutoffs(&rcut_, &rsw_); | 
| 1019 |  |  | 
| 1020 |  | } else { | 
| 1071 |  | int errorOut; | 
| 1072 |  | int esm =  NONE; | 
| 1073 |  | int sm = UNDAMPED; | 
| 1074 | < | double alphaVal; | 
| 1075 | < | double dielectric; | 
| 1074 | > | RealType alphaVal; | 
| 1075 | > | RealType dielectric; | 
| 1076 |  |  | 
| 1077 |  | errorOut = isError; | 
| 1078 |  | alphaVal = simParams_->getDampingAlpha(); | 
| 1173 |  |  | 
| 1174 |  | // send switching function notification to switcheroo | 
| 1175 |  | setFunctionType(&ft); | 
| 1176 | + |  | 
| 1177 | + | } | 
| 1178 | + |  | 
| 1179 | + | void SimInfo::setupAccumulateBoxDipole() { | 
| 1180 | + |  | 
| 1181 | + | // we only call setAccumulateBoxDipole if the accumulateBoxDipole parameter is true | 
| 1182 | + | if ( simParams_->haveAccumulateBoxDipole() ) | 
| 1183 | + | if ( simParams_->getAccumulateBoxDipole() ) { | 
| 1184 | + | setAccumulateBoxDipole(); | 
| 1185 | + | calcBoxDipole_ = true; | 
| 1186 | + | } | 
| 1187 |  |  | 
| 1188 |  | } | 
| 1189 |  |  | 
| 1243 |  | Molecule* mol; | 
| 1244 |  |  | 
| 1245 |  | Vector3d comVel(0.0); | 
| 1246 | < | double totalMass = 0.0; | 
| 1246 | > | RealType totalMass = 0.0; | 
| 1247 |  |  | 
| 1248 |  |  | 
| 1249 |  | for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { | 
| 1250 | < | double mass = mol->getMass(); | 
| 1250 | > | RealType mass = mol->getMass(); | 
| 1251 |  | totalMass += mass; | 
| 1252 |  | comVel += mass * mol->getComVel(); | 
| 1253 |  | } | 
| 1254 |  |  | 
| 1255 |  | #ifdef IS_MPI | 
| 1256 | < | double tmpMass = totalMass; | 
| 1256 | > | RealType tmpMass = totalMass; | 
| 1257 |  | Vector3d tmpComVel(comVel); | 
| 1258 | < | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1259 | < | MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1258 | > | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1259 | > | MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1260 |  | #endif | 
| 1261 |  |  | 
| 1262 |  | comVel /= totalMass; | 
| 1269 |  | Molecule* mol; | 
| 1270 |  |  | 
| 1271 |  | Vector3d com(0.0); | 
| 1272 | < | double totalMass = 0.0; | 
| 1272 | > | RealType totalMass = 0.0; | 
| 1273 |  |  | 
| 1274 |  | for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { | 
| 1275 | < | double mass = mol->getMass(); | 
| 1275 | > | RealType mass = mol->getMass(); | 
| 1276 |  | totalMass += mass; | 
| 1277 |  | com += mass * mol->getCom(); | 
| 1278 |  | } | 
| 1279 |  |  | 
| 1280 |  | #ifdef IS_MPI | 
| 1281 | < | double tmpMass = totalMass; | 
| 1281 | > | RealType tmpMass = totalMass; | 
| 1282 |  | Vector3d tmpCom(com); | 
| 1283 | < | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1284 | < | MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1283 | > | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1284 | > | MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1285 |  | #endif | 
| 1286 |  |  | 
| 1287 |  | com /= totalMass; | 
| 1305 |  | Molecule* mol; | 
| 1306 |  |  | 
| 1307 |  |  | 
| 1308 | < | double totalMass = 0.0; | 
| 1308 | > | RealType totalMass = 0.0; | 
| 1309 |  |  | 
| 1310 |  |  | 
| 1311 |  | for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { | 
| 1312 | < | double mass = mol->getMass(); | 
| 1312 | > | RealType mass = mol->getMass(); | 
| 1313 |  | totalMass += mass; | 
| 1314 |  | com += mass * mol->getCom(); | 
| 1315 |  | comVel += mass * mol->getComVel(); | 
| 1316 |  | } | 
| 1317 |  |  | 
| 1318 |  | #ifdef IS_MPI | 
| 1319 | < | double tmpMass = totalMass; | 
| 1319 | > | RealType tmpMass = totalMass; | 
| 1320 |  | Vector3d tmpCom(com); | 
| 1321 |  | Vector3d tmpComVel(comVel); | 
| 1322 | < | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1323 | < | MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1324 | < | MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1322 | > | MPI_Allreduce(&tmpMass,&totalMass,1,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1323 | > | MPI_Allreduce(tmpCom.getArrayPointer(), com.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1324 | > | MPI_Allreduce(tmpComVel.getArrayPointer(), comVel.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1325 |  | #endif | 
| 1326 |  |  | 
| 1327 |  | com /= totalMass; | 
| 1340 |  | void SimInfo::getInertiaTensor(Mat3x3d &inertiaTensor, Vector3d &angularMomentum){ | 
| 1341 |  |  | 
| 1342 |  |  | 
| 1343 | < | double xx = 0.0; | 
| 1344 | < | double yy = 0.0; | 
| 1345 | < | double zz = 0.0; | 
| 1346 | < | double xy = 0.0; | 
| 1347 | < | double xz = 0.0; | 
| 1348 | < | double yz = 0.0; | 
| 1343 | > | RealType xx = 0.0; | 
| 1344 | > | RealType yy = 0.0; | 
| 1345 | > | RealType zz = 0.0; | 
| 1346 | > | RealType xy = 0.0; | 
| 1347 | > | RealType xz = 0.0; | 
| 1348 | > | RealType yz = 0.0; | 
| 1349 |  | Vector3d com(0.0); | 
| 1350 |  | Vector3d comVel(0.0); | 
| 1351 |  |  | 
| 1357 |  | Vector3d thisq(0.0); | 
| 1358 |  | Vector3d thisv(0.0); | 
| 1359 |  |  | 
| 1360 | < | double thisMass = 0.0; | 
| 1360 | > | RealType thisMass = 0.0; | 
| 1361 |  |  | 
| 1362 |  |  | 
| 1363 |  |  | 
| 1395 |  | #ifdef IS_MPI | 
| 1396 |  | Mat3x3d tmpI(inertiaTensor); | 
| 1397 |  | Vector3d tmpAngMom; | 
| 1398 | < | MPI_Allreduce(tmpI.getArrayPointer(), inertiaTensor.getArrayPointer(),9,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1399 | < | MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1398 | > | MPI_Allreduce(tmpI.getArrayPointer(), inertiaTensor.getArrayPointer(),9,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1399 | > | MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1400 |  | #endif | 
| 1401 |  |  | 
| 1402 |  | return; | 
| 1417 |  | Vector3d thisr(0.0); | 
| 1418 |  | Vector3d thisp(0.0); | 
| 1419 |  |  | 
| 1420 | < | double thisMass; | 
| 1420 | > | RealType thisMass; | 
| 1421 |  |  | 
| 1422 |  | for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { | 
| 1423 |  | thisMass = mol->getMass(); | 
| 1430 |  |  | 
| 1431 |  | #ifdef IS_MPI | 
| 1432 |  | Vector3d tmpAngMom; | 
| 1433 | < | MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_DOUBLE,MPI_SUM, MPI_COMM_WORLD); | 
| 1433 | > | MPI_Allreduce(tmpAngMom.getArrayPointer(), angularMomentum.getArrayPointer(),3,MPI_REALTYPE,MPI_SUM, MPI_COMM_WORLD); | 
| 1434 |  | #endif | 
| 1435 |  |  | 
| 1436 |  | return angularMomentum; |