# | Line 35 | Line 35 | |
---|---|---|
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 | */ | |
# | Line 308 | Line 308 | namespace OpenMD { | |
308 | ||
309 | void ForceMatrixDecomposition::createGtypeCutoffMap() { | |
310 | ||
311 | + | GrCut.clear(); |
312 | + | GrCutSq.clear(); |
313 | + | GrlistSq.clear(); |
314 | + | |
315 | RealType tol = 1e-6; | |
316 | largestRcut_ = 0.0; | |
317 | int atid; | |
# | Line 418 | Line 422 | namespace OpenMD { | |
422 | #endif | |
423 | ||
424 | RealType tradRcut = groupMax; | |
425 | + | |
426 | + | GrCut.resize( gTypeCutoffs.size() ); |
427 | + | GrCutSq.resize( gTypeCutoffs.size() ); |
428 | + | GrlistSq.resize( gTypeCutoffs.size() ); |
429 | ||
430 | + | |
431 | for (unsigned int i = 0; i < gTypeCutoffs.size(); i++) { | |
432 | + | GrCut[i].resize( gTypeCutoffs.size() , 0.0); |
433 | + | GrCutSq[i].resize( gTypeCutoffs.size(), 0.0 ); |
434 | + | GrlistSq[i].resize( gTypeCutoffs.size(), 0.0 ); |
435 | + | |
436 | for (unsigned int j = 0; j < gTypeCutoffs.size(); j++) { | |
437 | RealType thisRcut; | |
438 | switch(cutoffPolicy_) { | |
# | Line 442 | Line 455 | namespace OpenMD { | |
455 | break; | |
456 | } | |
457 | ||
458 | < | pair<int,int> key = make_pair(i,j); |
446 | < | gTypeCutoffMap[key].first = thisRcut; |
458 | > | GrCut[i][j] = thisRcut; |
459 | if (thisRcut > largestRcut_) largestRcut_ = thisRcut; | |
460 | < | gTypeCutoffMap[key].second = thisRcut*thisRcut; |
461 | < | gTypeCutoffMap[key].third = pow(thisRcut + skinThickness_, 2); |
460 | > | GrCutSq[i][j] = thisRcut * thisRcut; |
461 | > | GrlistSq[i][j] = pow(thisRcut + skinThickness_, 2); |
462 | > | |
463 | > | // pair<int,int> key = make_pair(i,j); |
464 | > | // gTypeCutoffMap[key].first = thisRcut; |
465 | > | // gTypeCutoffMap[key].third = pow(thisRcut + skinThickness_, 2); |
466 | // sanity check | |
467 | ||
468 | if (userChoseCutoff_) { | |
469 | < | if (abs(gTypeCutoffMap[key].first - userCutoff_) > 0.0001) { |
469 | > | if (abs(GrCut[i][j] - userCutoff_) > 0.0001) { |
470 | sprintf(painCave.errMsg, | |
471 | "ForceMatrixDecomposition::createGtypeCutoffMap " | |
472 | "user-specified rCut (%lf) does not match computed group Cutoff\n", userCutoff_); | |
# | Line 463 | Line 479 | namespace OpenMD { | |
479 | } | |
480 | } | |
481 | ||
482 | < | groupCutoffs ForceMatrixDecomposition::getGroupCutoffs(int cg1, int cg2) { |
482 | > | void ForceMatrixDecomposition::getGroupCutoffs(int &cg1, int &cg2, RealType &rcut, RealType &rcutsq, RealType &rlistsq) { |
483 | int i, j; | |
484 | #ifdef IS_MPI | |
485 | i = groupRowToGtype[cg1]; | |
# | Line 472 | Line 488 | namespace OpenMD { | |
488 | i = groupToGtype[cg1]; | |
489 | j = groupToGtype[cg2]; | |
490 | #endif | |
491 | < | return gTypeCutoffMap[make_pair(i,j)]; |
491 | > | rcut = GrCut[i][j]; |
492 | > | rcutsq = GrCutSq[i][j]; |
493 | > | rlistsq = GrlistSq[i][j]; |
494 | > | return; |
495 | > | //return gTypeCutoffMap[make_pair(i,j)]; |
496 | } | |
497 | ||
498 | int ForceMatrixDecomposition::getTopologicalDistance(int atom1, int atom2) { | |
499 | for (unsigned int j = 0; j < toposForAtom[atom1].size(); j++) { | |
500 | if (toposForAtom[atom1][j] == atom2) | |
501 | return topoDist[atom1][j]; | |
502 | < | } |
502 | > | } |
503 | return 0; | |
504 | } | |
505 | ||
# | Line 559 | Line 579 | namespace OpenMD { | |
579 | atomColData.electricField.end(), V3Zero); | |
580 | } | |
581 | ||
562 | – | if (storageLayout_ & DataStorage::dslFlucQForce) { |
563 | – | fill(atomRowData.flucQFrc.begin(), atomRowData.flucQFrc.end(), |
564 | – | 0.0); |
565 | – | fill(atomColData.flucQFrc.begin(), atomColData.flucQFrc.end(), |
566 | – | 0.0); |
567 | – | } |
568 | – | |
582 | #endif | |
583 | // even in parallel, we need to zero out the local arrays: | |
584 | ||
# | Line 639 | Line 652 | namespace OpenMD { | |
652 | AtomPlanMatrixColumn->gather(snap_->atomData.aMat, | |
653 | atomColData.aMat); | |
654 | } | |
655 | < | |
656 | < | // if needed, gather the atomic eletrostatic frames |
657 | < | if (storageLayout_ & DataStorage::dslElectroFrame) { |
658 | < | AtomPlanMatrixRow->gather(snap_->atomData.electroFrame, |
659 | < | atomRowData.electroFrame); |
660 | < | AtomPlanMatrixColumn->gather(snap_->atomData.electroFrame, |
661 | < | atomColData.electroFrame); |
655 | > | |
656 | > | // if needed, gather the atomic eletrostatic information |
657 | > | if (storageLayout_ & DataStorage::dslDipole) { |
658 | > | AtomPlanVectorRow->gather(snap_->atomData.dipole, |
659 | > | atomRowData.dipole); |
660 | > | AtomPlanVectorColumn->gather(snap_->atomData.dipole, |
661 | > | atomColData.dipole); |
662 | } | |
663 | ||
664 | + | if (storageLayout_ & DataStorage::dslQuadrupole) { |
665 | + | AtomPlanMatrixRow->gather(snap_->atomData.quadrupole, |
666 | + | atomRowData.quadrupole); |
667 | + | AtomPlanMatrixColumn->gather(snap_->atomData.quadrupole, |
668 | + | atomColData.quadrupole); |
669 | + | } |
670 | + | |
671 | // if needed, gather the atomic fluctuating charge values | |
672 | if (storageLayout_ & DataStorage::dslFlucQPosition) { | |
673 | AtomPlanRealRow->gather(snap_->atomData.flucQPos, | |
# | Line 679 | Line 699 | namespace OpenMD { | |
699 | snap_->atomData.density[i] += rho_tmp[i]; | |
700 | } | |
701 | ||
702 | + | // this isn't necessary if we don't have polarizable atoms, but |
703 | + | // we'll leave it here for now. |
704 | if (storageLayout_ & DataStorage::dslElectricField) { | |
705 | ||
706 | AtomPlanVectorRow->scatter(atomRowData.electricField, | |
# | Line 686 | Line 708 | namespace OpenMD { | |
708 | ||
709 | int n = snap_->atomData.electricField.size(); | |
710 | vector<Vector3d> field_tmp(n, V3Zero); | |
711 | < | AtomPlanVectorColumn->scatter(atomColData.electricField, field_tmp); |
711 | > | AtomPlanVectorColumn->scatter(atomColData.electricField, |
712 | > | field_tmp); |
713 | for (int i = 0; i < n; i++) | |
714 | snap_->atomData.electricField[i] += field_tmp[i]; | |
715 | } | |
# | Line 784 | Line 807 | namespace OpenMD { | |
807 | for (int i = 0; i < nq; i++) | |
808 | snap_->atomData.flucQFrc[i] += fqfrc_tmp[i]; | |
809 | ||
810 | + | } |
811 | + | |
812 | + | if (storageLayout_ & DataStorage::dslElectricField) { |
813 | + | |
814 | + | int nef = snap_->atomData.electricField.size(); |
815 | + | vector<Vector3d> efield_tmp(nef, V3Zero); |
816 | + | |
817 | + | AtomPlanVectorRow->scatter(atomRowData.electricField, efield_tmp); |
818 | + | for (int i = 0; i < nef; i++) { |
819 | + | snap_->atomData.electricField[i] += efield_tmp[i]; |
820 | + | efield_tmp[i] = 0.0; |
821 | + | } |
822 | + | |
823 | + | AtomPlanVectorColumn->scatter(atomColData.electricField, efield_tmp); |
824 | + | for (int i = 0; i < nef; i++) |
825 | + | snap_->atomData.electricField[i] += efield_tmp[i]; |
826 | } | |
827 | ||
828 | + | |
829 | nLocal_ = snap_->getNumberOfAtoms(); | |
830 | ||
831 | vector<potVec> pot_temp(nLocal_, | |
# | Line 919 | Line 959 | namespace OpenMD { | |
959 | ||
960 | ||
961 | ||
962 | < | int ForceMatrixDecomposition::getNAtomsInRow() { |
962 | > | int& ForceMatrixDecomposition::getNAtomsInRow() { |
963 | #ifdef IS_MPI | |
964 | return nAtomsInRow_; | |
965 | #else | |
# | Line 930 | Line 970 | namespace OpenMD { | |
970 | /** | |
971 | * returns the list of atoms belonging to this group. | |
972 | */ | |
973 | < | vector<int> ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){ |
973 | > | vector<int>& ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){ |
974 | #ifdef IS_MPI | |
975 | return groupListRow_[cg1]; | |
976 | #else | |
# | Line 938 | Line 978 | namespace OpenMD { | |
978 | #endif | |
979 | } | |
980 | ||
981 | < | vector<int> ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){ |
981 | > | vector<int>& ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){ |
982 | #ifdef IS_MPI | |
983 | return groupListCol_[cg2]; | |
984 | #else | |
# | Line 955 | Line 995 | namespace OpenMD { | |
995 | d = snap_->cgData.position[cg2] - snap_->cgData.position[cg1]; | |
996 | #endif | |
997 | ||
998 | < | snap_->wrapVector(d); |
998 | > | if (usePeriodicBoundaryConditions_) { |
999 | > | snap_->wrapVector(d); |
1000 | > | } |
1001 | return d; | |
1002 | } | |
1003 | ||
1004 | < | Vector3d ForceMatrixDecomposition::getGroupVelocityColumn(int cg2){ |
1004 | > | Vector3d& ForceMatrixDecomposition::getGroupVelocityColumn(int cg2){ |
1005 | #ifdef IS_MPI | |
1006 | return cgColData.velocity[cg2]; | |
1007 | #else | |
# | Line 967 | Line 1009 | namespace OpenMD { | |
1009 | #endif | |
1010 | } | |
1011 | ||
1012 | < | Vector3d ForceMatrixDecomposition::getAtomVelocityColumn(int atom2){ |
1012 | > | Vector3d& ForceMatrixDecomposition::getAtomVelocityColumn(int atom2){ |
1013 | #ifdef IS_MPI | |
1014 | return atomColData.velocity[atom2]; | |
1015 | #else | |
# | Line 985 | Line 1027 | namespace OpenMD { | |
1027 | #else | |
1028 | d = snap_->cgData.position[cg1] - snap_->atomData.position[atom1]; | |
1029 | #endif | |
1030 | < | |
1031 | < | snap_->wrapVector(d); |
1030 | > | if (usePeriodicBoundaryConditions_) { |
1031 | > | snap_->wrapVector(d); |
1032 | > | } |
1033 | return d; | |
1034 | } | |
1035 | ||
# | Line 998 | Line 1041 | namespace OpenMD { | |
1041 | #else | |
1042 | d = snap_->cgData.position[cg2] - snap_->atomData.position[atom2]; | |
1043 | #endif | |
1044 | < | |
1045 | < | snap_->wrapVector(d); |
1044 | > | if (usePeriodicBoundaryConditions_) { |
1045 | > | snap_->wrapVector(d); |
1046 | > | } |
1047 | return d; | |
1048 | } | |
1049 | ||
1050 | < | RealType ForceMatrixDecomposition::getMassFactorRow(int atom1) { |
1050 | > | RealType& ForceMatrixDecomposition::getMassFactorRow(int atom1) { |
1051 | #ifdef IS_MPI | |
1052 | return massFactorsRow[atom1]; | |
1053 | #else | |
# | Line 1011 | Line 1055 | namespace OpenMD { | |
1055 | #endif | |
1056 | } | |
1057 | ||
1058 | < | RealType ForceMatrixDecomposition::getMassFactorColumn(int atom2) { |
1058 | > | RealType& ForceMatrixDecomposition::getMassFactorColumn(int atom2) { |
1059 | #ifdef IS_MPI | |
1060 | return massFactorsCol[atom2]; | |
1061 | #else | |
# | Line 1028 | Line 1072 | namespace OpenMD { | |
1072 | #else | |
1073 | d = snap_->atomData.position[atom2] - snap_->atomData.position[atom1]; | |
1074 | #endif | |
1075 | < | |
1076 | < | snap_->wrapVector(d); |
1075 | > | if (usePeriodicBoundaryConditions_) { |
1076 | > | snap_->wrapVector(d); |
1077 | > | } |
1078 | return d; | |
1079 | } | |
1080 | ||
1081 | < | vector<int> ForceMatrixDecomposition::getExcludesForAtom(int atom1) { |
1081 | > | vector<int>& ForceMatrixDecomposition::getExcludesForAtom(int atom1) { |
1082 | return excludesForAtom[atom1]; | |
1083 | } | |
1084 | ||
# | Line 1124 | Line 1169 | namespace OpenMD { | |
1169 | ||
1170 | #ifdef IS_MPI | |
1171 | idat.atypes = make_pair( atypesRow[atom1], atypesCol[atom2]); | |
1172 | + | idat.atid1 = identsRow[atom1]; |
1173 | + | idat.atid2 = identsCol[atom2]; |
1174 | //idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), | |
1175 | // ff_->getAtomType(identsCol[atom2]) ); | |
1176 | ||
# | Line 1132 | Line 1179 | namespace OpenMD { | |
1179 | idat.A2 = &(atomColData.aMat[atom2]); | |
1180 | } | |
1181 | ||
1135 | – | if (storageLayout_ & DataStorage::dslElectroFrame) { |
1136 | – | idat.eFrame1 = &(atomRowData.electroFrame[atom1]); |
1137 | – | idat.eFrame2 = &(atomColData.electroFrame[atom2]); |
1138 | – | } |
1139 | – | |
1182 | if (storageLayout_ & DataStorage::dslTorque) { | |
1183 | idat.t1 = &(atomRowData.torque[atom1]); | |
1184 | idat.t2 = &(atomColData.torque[atom2]); | |
1185 | } | |
1186 | ||
1187 | + | if (storageLayout_ & DataStorage::dslDipole) { |
1188 | + | idat.dipole1 = &(atomRowData.dipole[atom1]); |
1189 | + | idat.dipole2 = &(atomColData.dipole[atom2]); |
1190 | + | } |
1191 | + | |
1192 | + | if (storageLayout_ & DataStorage::dslQuadrupole) { |
1193 | + | idat.quadrupole1 = &(atomRowData.quadrupole[atom1]); |
1194 | + | idat.quadrupole2 = &(atomColData.quadrupole[atom2]); |
1195 | + | } |
1196 | + | |
1197 | if (storageLayout_ & DataStorage::dslDensity) { | |
1198 | idat.rho1 = &(atomRowData.density[atom1]); | |
1199 | idat.rho2 = &(atomColData.density[atom2]); | |
# | Line 1175 | Line 1227 | namespace OpenMD { | |
1227 | #else | |
1228 | ||
1229 | idat.atypes = make_pair( atypesLocal[atom1], atypesLocal[atom2]); | |
1230 | + | idat.atid1 = idents[atom1]; |
1231 | + | idat.atid2 = idents[atom2]; |
1232 | ||
1233 | if (storageLayout_ & DataStorage::dslAmat) { | |
1234 | idat.A1 = &(snap_->atomData.aMat[atom1]); | |
1235 | idat.A2 = &(snap_->atomData.aMat[atom2]); | |
1236 | } | |
1237 | ||
1184 | – | if (storageLayout_ & DataStorage::dslElectroFrame) { |
1185 | – | idat.eFrame1 = &(snap_->atomData.electroFrame[atom1]); |
1186 | – | idat.eFrame2 = &(snap_->atomData.electroFrame[atom2]); |
1187 | – | } |
1188 | – | |
1238 | if (storageLayout_ & DataStorage::dslTorque) { | |
1239 | idat.t1 = &(snap_->atomData.torque[atom1]); | |
1240 | idat.t2 = &(snap_->atomData.torque[atom2]); | |
1241 | } | |
1242 | ||
1243 | + | if (storageLayout_ & DataStorage::dslDipole) { |
1244 | + | idat.dipole1 = &(snap_->atomData.dipole[atom1]); |
1245 | + | idat.dipole2 = &(snap_->atomData.dipole[atom2]); |
1246 | + | } |
1247 | + | |
1248 | + | if (storageLayout_ & DataStorage::dslQuadrupole) { |
1249 | + | idat.quadrupole1 = &(snap_->atomData.quadrupole[atom1]); |
1250 | + | idat.quadrupole2 = &(snap_->atomData.quadrupole[atom2]); |
1251 | + | } |
1252 | + | |
1253 | if (storageLayout_ & DataStorage::dslDensity) { | |
1254 | idat.rho1 = &(snap_->atomData.density[atom1]); | |
1255 | idat.rho2 = &(snap_->atomData.density[atom2]); | |
# | Line 1281 | Line 1340 | namespace OpenMD { | |
1340 | * first element of pair is row-indexed CutoffGroup | |
1341 | * second element of pair is column-indexed CutoffGroup | |
1342 | */ | |
1343 | < | vector<pair<int, int> > ForceMatrixDecomposition::buildNeighborList() { |
1344 | < | |
1345 | < | vector<pair<int, int> > neighborList; |
1343 | > | void ForceMatrixDecomposition::buildNeighborList(vector<pair<int,int> >& neighborList) { |
1344 | > | |
1345 | > | neighborList.clear(); |
1346 | groupCutoffs cuts; | |
1347 | bool doAllPairs = false; | |
1348 | ||
1349 | + | RealType rList_ = (largestRcut_ + skinThickness_); |
1350 | + | RealType rcut, rcutsq, rlistsq; |
1351 | + | Snapshot* snap_ = sman_->getCurrentSnapshot(); |
1352 | + | Mat3x3d box; |
1353 | + | Mat3x3d invBox; |
1354 | + | |
1355 | + | Vector3d rs, scaled, dr; |
1356 | + | Vector3i whichCell; |
1357 | + | int cellIndex; |
1358 | + | |
1359 | #ifdef IS_MPI | |
1360 | cellListRow_.clear(); | |
1361 | cellListCol_.clear(); | |
1362 | #else | |
1363 | cellList_.clear(); | |
1364 | #endif | |
1365 | < | |
1366 | < | RealType rList_ = (largestRcut_ + skinThickness_); |
1367 | < | Snapshot* snap_ = sman_->getCurrentSnapshot(); |
1368 | < | Mat3x3d Hmat = snap_->getHmat(); |
1369 | < | Vector3d Hx = Hmat.getColumn(0); |
1370 | < | Vector3d Hy = Hmat.getColumn(1); |
1371 | < | Vector3d Hz = Hmat.getColumn(2); |
1372 | < | |
1373 | < | nCells_.x() = (int) ( Hx.length() )/ rList_; |
1374 | < | nCells_.y() = (int) ( Hy.length() )/ rList_; |
1375 | < | nCells_.z() = (int) ( Hz.length() )/ rList_; |
1376 | < | |
1365 | > | |
1366 | > | if (!usePeriodicBoundaryConditions_) { |
1367 | > | box = snap_->getBoundingBox(); |
1368 | > | invBox = snap_->getInvBoundingBox(); |
1369 | > | } else { |
1370 | > | box = snap_->getHmat(); |
1371 | > | invBox = snap_->getInvHmat(); |
1372 | > | } |
1373 | > | |
1374 | > | Vector3d boxX = box.getColumn(0); |
1375 | > | Vector3d boxY = box.getColumn(1); |
1376 | > | Vector3d boxZ = box.getColumn(2); |
1377 | > | |
1378 | > | nCells_.x() = (int) ( boxX.length() )/ rList_; |
1379 | > | nCells_.y() = (int) ( boxY.length() )/ rList_; |
1380 | > | nCells_.z() = (int) ( boxZ.length() )/ rList_; |
1381 | > | |
1382 | // handle small boxes where the cell offsets can end up repeating cells | |
1383 | ||
1384 | if (nCells_.x() < 3) doAllPairs = true; | |
1385 | if (nCells_.y() < 3) doAllPairs = true; | |
1386 | if (nCells_.z() < 3) doAllPairs = true; | |
1387 | < | |
1314 | < | Mat3x3d invHmat = snap_->getInvHmat(); |
1315 | < | Vector3d rs, scaled, dr; |
1316 | < | Vector3i whichCell; |
1317 | < | int cellIndex; |
1387 | > | |
1388 | int nCtot = nCells_.x() * nCells_.y() * nCells_.z(); | |
1389 | < | |
1389 | > | |
1390 | #ifdef IS_MPI | |
1391 | cellListRow_.resize(nCtot); | |
1392 | cellListCol_.resize(nCtot); | |
1393 | #else | |
1394 | cellList_.resize(nCtot); | |
1395 | #endif | |
1396 | < | |
1396 | > | |
1397 | if (!doAllPairs) { | |
1398 | #ifdef IS_MPI | |
1399 | < | |
1399 | > | |
1400 | for (int i = 0; i < nGroupsInRow_; i++) { | |
1401 | rs = cgRowData.position[i]; | |
1402 | ||
1403 | // scaled positions relative to the box vectors | |
1404 | < | scaled = invHmat * rs; |
1404 | > | scaled = invBox * rs; |
1405 | ||
1406 | // wrap the vector back into the unit box by subtracting integer box | |
1407 | // numbers | |
# | Line 1359 | Line 1429 | namespace OpenMD { | |
1429 | rs = cgColData.position[i]; | |
1430 | ||
1431 | // scaled positions relative to the box vectors | |
1432 | < | scaled = invHmat * rs; |
1432 | > | scaled = invBox * rs; |
1433 | ||
1434 | // wrap the vector back into the unit box by subtracting integer box | |
1435 | // numbers | |
# | Line 1383 | Line 1453 | namespace OpenMD { | |
1453 | // add this cutoff group to the list of groups in this cell; | |
1454 | cellListCol_[cellIndex].push_back(i); | |
1455 | } | |
1456 | < | |
1456 | > | |
1457 | #else | |
1458 | for (int i = 0; i < nGroups_; i++) { | |
1459 | rs = snap_->cgData.position[i]; | |
1460 | ||
1461 | // scaled positions relative to the box vectors | |
1462 | < | scaled = invHmat * rs; |
1462 | > | scaled = invBox * rs; |
1463 | ||
1464 | // wrap the vector back into the unit box by subtracting integer box | |
1465 | // numbers | |
# | Line 1458 | Line 1528 | namespace OpenMD { | |
1528 | // & column indicies and will divide labor in the | |
1529 | // force evaluation later. | |
1530 | dr = cgColData.position[(*j2)] - cgRowData.position[(*j1)]; | |
1531 | < | snap_->wrapVector(dr); |
1532 | < | cuts = getGroupCutoffs( (*j1), (*j2) ); |
1533 | < | if (dr.lengthSquare() < cuts.third) { |
1531 | > | if (usePeriodicBoundaryConditions_) { |
1532 | > | snap_->wrapVector(dr); |
1533 | > | } |
1534 | > | getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1535 | > | if (dr.lengthSquare() < rlistsq) { |
1536 | neighborList.push_back(make_pair((*j1), (*j2))); | |
1537 | } | |
1538 | } | |
# | Line 1480 | Line 1552 | namespace OpenMD { | |
1552 | // allows atoms within a single cutoff group to | |
1553 | // interact with each other. | |
1554 | ||
1483 | – | |
1484 | – | |
1555 | if (m2 != m1 || (*j2) >= (*j1) ) { | |
1556 | ||
1557 | dr = snap_->cgData.position[(*j2)] - snap_->cgData.position[(*j1)]; | |
1558 | < | snap_->wrapVector(dr); |
1559 | < | cuts = getGroupCutoffs( (*j1), (*j2) ); |
1560 | < | if (dr.lengthSquare() < cuts.third) { |
1558 | > | if (usePeriodicBoundaryConditions_) { |
1559 | > | snap_->wrapVector(dr); |
1560 | > | } |
1561 | > | getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1562 | > | if (dr.lengthSquare() < rlistsq) { |
1563 | neighborList.push_back(make_pair((*j1), (*j2))); | |
1564 | } | |
1565 | } | |
# | Line 1504 | Line 1576 | namespace OpenMD { | |
1576 | for (int j1 = 0; j1 < nGroupsInRow_; j1++) { | |
1577 | for (int j2 = 0; j2 < nGroupsInCol_; j2++) { | |
1578 | dr = cgColData.position[j2] - cgRowData.position[j1]; | |
1579 | < | snap_->wrapVector(dr); |
1580 | < | cuts = getGroupCutoffs( j1, j2 ); |
1581 | < | if (dr.lengthSquare() < cuts.third) { |
1579 | > | if (usePeriodicBoundaryConditions_) { |
1580 | > | snap_->wrapVector(dr); |
1581 | > | } |
1582 | > | getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq); |
1583 | > | if (dr.lengthSquare() < rlistsq) { |
1584 | neighborList.push_back(make_pair(j1, j2)); | |
1585 | } | |
1586 | } | |
# | Line 1517 | Line 1591 | namespace OpenMD { | |
1591 | // include self group interactions j2 == j1 | |
1592 | for (int j2 = j1; j2 < nGroups_; j2++) { | |
1593 | dr = snap_->cgData.position[j2] - snap_->cgData.position[j1]; | |
1594 | < | snap_->wrapVector(dr); |
1595 | < | cuts = getGroupCutoffs( j1, j2 ); |
1596 | < | if (dr.lengthSquare() < cuts.third) { |
1594 | > | if (usePeriodicBoundaryConditions_) { |
1595 | > | snap_->wrapVector(dr); |
1596 | > | } |
1597 | > | getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq ); |
1598 | > | if (dr.lengthSquare() < rlistsq) { |
1599 | neighborList.push_back(make_pair(j1, j2)); | |
1600 | } | |
1601 | } | |
# | Line 1532 | Line 1608 | namespace OpenMD { | |
1608 | saved_CG_positions_.clear(); | |
1609 | for (int i = 0; i < nGroups_; i++) | |
1610 | saved_CG_positions_.push_back(snap_->cgData.position[i]); | |
1535 | – | |
1536 | – | return neighborList; |
1611 | } | |
1612 | } //end namespace OpenMD |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |