# | Line 99 | Line 99 | namespace OpenMD { | |
---|---|---|
99 | nGroups_ = info_->getNLocalCutoffGroups(); | |
100 | // gather the information for atomtype IDs (atids): | |
101 | idents = info_->getIdentArray(); | |
102 | + | regions = info_->getRegions(); |
103 | AtomLocalToGlobal = info_->getGlobalAtomIndices(); | |
104 | cgLocalToGlobal = info_->getGlobalGroupIndices(); | |
105 | vector<int> globalGroupMembership = info_->getGlobalGroupMembership(); | |
# | Line 163 | Line 164 | namespace OpenMD { | |
164 | ||
165 | AtomPlanIntRow->gather(idents, identsRow); | |
166 | AtomPlanIntColumn->gather(idents, identsCol); | |
167 | + | |
168 | + | regionsRow.resize(nAtomsInRow_); |
169 | + | regionsCol.resize(nAtomsInCol_); |
170 | + | |
171 | + | AtomPlanIntRow->gather(regions, regionsRow); |
172 | + | AtomPlanIntColumn->gather(regions, regionsCol); |
173 | ||
174 | // allocate memory for the parallel objects | |
175 | atypesRow.resize(nAtomsInRow_); | |
# | Line 308 | Line 315 | namespace OpenMD { | |
315 | ||
316 | void ForceMatrixDecomposition::createGtypeCutoffMap() { | |
317 | ||
318 | + | GrCut.clear(); |
319 | + | GrCutSq.clear(); |
320 | + | GrlistSq.clear(); |
321 | + | |
322 | RealType tol = 1e-6; | |
323 | largestRcut_ = 0.0; | |
324 | int atid; | |
# | Line 419 | Line 430 | namespace OpenMD { | |
430 | ||
431 | RealType tradRcut = groupMax; | |
432 | ||
433 | + | GrCut.resize( gTypeCutoffs.size() ); |
434 | + | GrCutSq.resize( gTypeCutoffs.size() ); |
435 | + | GrlistSq.resize( gTypeCutoffs.size() ); |
436 | + | |
437 | + | |
438 | for (unsigned int i = 0; i < gTypeCutoffs.size(); i++) { | |
439 | + | GrCut[i].resize( gTypeCutoffs.size() , 0.0); |
440 | + | GrCutSq[i].resize( gTypeCutoffs.size(), 0.0 ); |
441 | + | GrlistSq[i].resize( gTypeCutoffs.size(), 0.0 ); |
442 | + | |
443 | for (unsigned int j = 0; j < gTypeCutoffs.size(); j++) { | |
444 | RealType thisRcut; | |
445 | switch(cutoffPolicy_) { | |
# | Line 442 | Line 462 | namespace OpenMD { | |
462 | break; | |
463 | } | |
464 | ||
465 | < | pair<int,int> key = make_pair(i,j); |
446 | < | gTypeCutoffMap[key].first = thisRcut; |
465 | > | GrCut[i][j] = thisRcut; |
466 | if (thisRcut > largestRcut_) largestRcut_ = thisRcut; | |
467 | < | gTypeCutoffMap[key].second = thisRcut*thisRcut; |
468 | < | gTypeCutoffMap[key].third = pow(thisRcut + skinThickness_, 2); |
467 | > | GrCutSq[i][j] = thisRcut * thisRcut; |
468 | > | GrlistSq[i][j] = pow(thisRcut + skinThickness_, 2); |
469 | > | |
470 | > | // pair<int,int> key = make_pair(i,j); |
471 | > | // gTypeCutoffMap[key].first = thisRcut; |
472 | > | // gTypeCutoffMap[key].third = pow(thisRcut + skinThickness_, 2); |
473 | // sanity check | |
474 | ||
475 | if (userChoseCutoff_) { | |
476 | < | if (abs(gTypeCutoffMap[key].first - userCutoff_) > 0.0001) { |
476 | > | if (abs(GrCut[i][j] - userCutoff_) > 0.0001) { |
477 | sprintf(painCave.errMsg, | |
478 | "ForceMatrixDecomposition::createGtypeCutoffMap " | |
479 | "user-specified rCut (%lf) does not match computed group Cutoff\n", userCutoff_); | |
# | Line 463 | Line 486 | namespace OpenMD { | |
486 | } | |
487 | } | |
488 | ||
489 | < | groupCutoffs ForceMatrixDecomposition::getGroupCutoffs(int cg1, int cg2) { |
489 | > | void ForceMatrixDecomposition::getGroupCutoffs(int &cg1, int &cg2, RealType &rcut, RealType &rcutsq, RealType &rlistsq) { |
490 | int i, j; | |
491 | #ifdef IS_MPI | |
492 | i = groupRowToGtype[cg1]; | |
# | Line 472 | Line 495 | namespace OpenMD { | |
495 | i = groupToGtype[cg1]; | |
496 | j = groupToGtype[cg2]; | |
497 | #endif | |
498 | < | return gTypeCutoffMap[make_pair(i,j)]; |
498 | > | rcut = GrCut[i][j]; |
499 | > | rcutsq = GrCutSq[i][j]; |
500 | > | rlistsq = GrlistSq[i][j]; |
501 | > | return; |
502 | > | //return gTypeCutoffMap[make_pair(i,j)]; |
503 | } | |
504 | ||
505 | int ForceMatrixDecomposition::getTopologicalDistance(int atom1, int atom2) { | |
506 | for (unsigned int j = 0; j < toposForAtom[atom1].size(); j++) { | |
507 | if (toposForAtom[atom1][j] == atom2) | |
508 | return topoDist[atom1][j]; | |
509 | < | } |
509 | > | } |
510 | return 0; | |
511 | } | |
512 | ||
# | Line 939 | Line 966 | namespace OpenMD { | |
966 | ||
967 | ||
968 | ||
969 | < | int ForceMatrixDecomposition::getNAtomsInRow() { |
969 | > | int& ForceMatrixDecomposition::getNAtomsInRow() { |
970 | #ifdef IS_MPI | |
971 | return nAtomsInRow_; | |
972 | #else | |
# | Line 950 | Line 977 | namespace OpenMD { | |
977 | /** | |
978 | * returns the list of atoms belonging to this group. | |
979 | */ | |
980 | < | vector<int> ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){ |
980 | > | vector<int>& ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){ |
981 | #ifdef IS_MPI | |
982 | return groupListRow_[cg1]; | |
983 | #else | |
# | Line 958 | Line 985 | namespace OpenMD { | |
985 | #endif | |
986 | } | |
987 | ||
988 | < | vector<int> ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){ |
988 | > | vector<int>& ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){ |
989 | #ifdef IS_MPI | |
990 | return groupListCol_[cg2]; | |
991 | #else | |
# | Line 981 | Line 1008 | namespace OpenMD { | |
1008 | return d; | |
1009 | } | |
1010 | ||
1011 | < | Vector3d ForceMatrixDecomposition::getGroupVelocityColumn(int cg2){ |
1011 | > | Vector3d& ForceMatrixDecomposition::getGroupVelocityColumn(int cg2){ |
1012 | #ifdef IS_MPI | |
1013 | return cgColData.velocity[cg2]; | |
1014 | #else | |
# | Line 989 | Line 1016 | namespace OpenMD { | |
1016 | #endif | |
1017 | } | |
1018 | ||
1019 | < | Vector3d ForceMatrixDecomposition::getAtomVelocityColumn(int atom2){ |
1019 | > | Vector3d& ForceMatrixDecomposition::getAtomVelocityColumn(int atom2){ |
1020 | #ifdef IS_MPI | |
1021 | return atomColData.velocity[atom2]; | |
1022 | #else | |
# | Line 1027 | Line 1054 | namespace OpenMD { | |
1054 | return d; | |
1055 | } | |
1056 | ||
1057 | < | RealType ForceMatrixDecomposition::getMassFactorRow(int atom1) { |
1057 | > | RealType& ForceMatrixDecomposition::getMassFactorRow(int atom1) { |
1058 | #ifdef IS_MPI | |
1059 | return massFactorsRow[atom1]; | |
1060 | #else | |
# | Line 1035 | Line 1062 | namespace OpenMD { | |
1062 | #endif | |
1063 | } | |
1064 | ||
1065 | < | RealType ForceMatrixDecomposition::getMassFactorColumn(int atom2) { |
1065 | > | RealType& ForceMatrixDecomposition::getMassFactorColumn(int atom2) { |
1066 | #ifdef IS_MPI | |
1067 | return massFactorsCol[atom2]; | |
1068 | #else | |
# | Line 1058 | Line 1085 | namespace OpenMD { | |
1085 | return d; | |
1086 | } | |
1087 | ||
1088 | < | vector<int> ForceMatrixDecomposition::getExcludesForAtom(int atom1) { |
1088 | > | vector<int>& ForceMatrixDecomposition::getExcludesForAtom(int atom1) { |
1089 | return excludesForAtom[atom1]; | |
1090 | } | |
1091 | ||
# | Line 1148 | Line 1175 | namespace OpenMD { | |
1175 | idat.excluded = excludeAtomPair(atom1, atom2); | |
1176 | ||
1177 | #ifdef IS_MPI | |
1178 | < | idat.atypes = make_pair( atypesRow[atom1], atypesCol[atom2]); |
1179 | < | //idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), |
1180 | < | // ff_->getAtomType(identsCol[atom2]) ); |
1181 | < | |
1178 | > | //idat.atypes = make_pair( atypesRow[atom1], atypesCol[atom2]); |
1179 | > | idat.atid1 = identsRow[atom1]; |
1180 | > | idat.atid2 = identsCol[atom2]; |
1181 | > | |
1182 | > | if (regionsRow[atom1] >= 0 && regionsCol[atom2] >= 0) { |
1183 | > | idat.sameRegion = (regionsRow[atom1] == regionsCol[atom2]); |
1184 | > | } else { |
1185 | > | idat.sameRegion = false; |
1186 | > | } |
1187 | > | |
1188 | if (storageLayout_ & DataStorage::dslAmat) { | |
1189 | idat.A1 = &(atomRowData.aMat[atom1]); | |
1190 | idat.A2 = &(atomColData.aMat[atom2]); | |
# | Line 1204 | Line 1237 | namespace OpenMD { | |
1237 | ||
1238 | #else | |
1239 | ||
1240 | < | idat.atypes = make_pair( atypesLocal[atom1], atypesLocal[atom2]); |
1240 | > | //idat.atypes = make_pair( atypesLocal[atom1], atypesLocal[atom2]); |
1241 | > | idat.atid1 = idents[atom1]; |
1242 | > | idat.atid2 = idents[atom2]; |
1243 | ||
1244 | + | if (regions[atom1] >= 0 && regions[atom2] >= 0) { |
1245 | + | idat.sameRegion = (regions[atom1] == regions[atom2]); |
1246 | + | } else { |
1247 | + | idat.sameRegion = false; |
1248 | + | } |
1249 | + | |
1250 | if (storageLayout_ & DataStorage::dslAmat) { | |
1251 | idat.A1 = &(snap_->atomData.aMat[atom1]); | |
1252 | idat.A2 = &(snap_->atomData.aMat[atom2]); | |
# | Line 1316 | Line 1357 | namespace OpenMD { | |
1357 | * first element of pair is row-indexed CutoffGroup | |
1358 | * second element of pair is column-indexed CutoffGroup | |
1359 | */ | |
1360 | < | vector<pair<int, int> > ForceMatrixDecomposition::buildNeighborList() { |
1361 | < | |
1362 | < | vector<pair<int, int> > neighborList; |
1360 | > | void ForceMatrixDecomposition::buildNeighborList(vector<pair<int,int> >& neighborList) { |
1361 | > | |
1362 | > | neighborList.clear(); |
1363 | groupCutoffs cuts; | |
1364 | bool doAllPairs = false; | |
1365 | ||
1366 | RealType rList_ = (largestRcut_ + skinThickness_); | |
1367 | + | RealType rcut, rcutsq, rlistsq; |
1368 | Snapshot* snap_ = sman_->getCurrentSnapshot(); | |
1369 | Mat3x3d box; | |
1370 | Mat3x3d invBox; | |
# | Line 1506 | Line 1548 | namespace OpenMD { | |
1548 | if (usePeriodicBoundaryConditions_) { | |
1549 | snap_->wrapVector(dr); | |
1550 | } | |
1551 | < | cuts = getGroupCutoffs( (*j1), (*j2) ); |
1552 | < | if (dr.lengthSquare() < cuts.third) { |
1551 | > | getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1552 | > | if (dr.lengthSquare() < rlistsq) { |
1553 | neighborList.push_back(make_pair((*j1), (*j2))); | |
1554 | } | |
1555 | } | |
# | Line 1533 | Line 1575 | namespace OpenMD { | |
1575 | if (usePeriodicBoundaryConditions_) { | |
1576 | snap_->wrapVector(dr); | |
1577 | } | |
1578 | < | cuts = getGroupCutoffs( (*j1), (*j2) ); |
1579 | < | if (dr.lengthSquare() < cuts.third) { |
1578 | > | getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1579 | > | if (dr.lengthSquare() < rlistsq) { |
1580 | neighborList.push_back(make_pair((*j1), (*j2))); | |
1581 | } | |
1582 | } | |
# | Line 1554 | Line 1596 | namespace OpenMD { | |
1596 | if (usePeriodicBoundaryConditions_) { | |
1597 | snap_->wrapVector(dr); | |
1598 | } | |
1599 | < | cuts = getGroupCutoffs( j1, j2 ); |
1600 | < | if (dr.lengthSquare() < cuts.third) { |
1599 | > | getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq); |
1600 | > | if (dr.lengthSquare() < rlistsq) { |
1601 | neighborList.push_back(make_pair(j1, j2)); | |
1602 | } | |
1603 | } | |
# | Line 1569 | Line 1611 | namespace OpenMD { | |
1611 | if (usePeriodicBoundaryConditions_) { | |
1612 | snap_->wrapVector(dr); | |
1613 | } | |
1614 | < | cuts = getGroupCutoffs( j1, j2 ); |
1615 | < | if (dr.lengthSquare() < cuts.third) { |
1614 | > | getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq ); |
1615 | > | if (dr.lengthSquare() < rlistsq) { |
1616 | neighborList.push_back(make_pair(j1, j2)); | |
1617 | } | |
1618 | } | |
# | Line 1583 | Line 1625 | namespace OpenMD { | |
1625 | saved_CG_positions_.clear(); | |
1626 | for (int i = 0; i < nGroups_; i++) | |
1627 | saved_CG_positions_.push_back(snap_->cgData.position[i]); | |
1586 | – | |
1587 | – | return neighborList; |
1628 | } | |
1629 | } //end namespace OpenMD |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |