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; |
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_) { |
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_); |
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]; |
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) { |
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 |
|
|
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]); |
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; |
1531 |
|
if (usePeriodicBoundaryConditions_) { |
1532 |
|
snap_->wrapVector(dr); |
1533 |
|
} |
1534 |
< |
cuts = getGroupCutoffs( (*j1), (*j2) ); |
1535 |
< |
if (dr.lengthSquare() < cuts.third) { |
1534 |
> |
getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1535 |
> |
if (dr.lengthSquare() < rlistsq) { |
1536 |
|
neighborList.push_back(make_pair((*j1), (*j2))); |
1537 |
|
} |
1538 |
|
} |
1558 |
|
if (usePeriodicBoundaryConditions_) { |
1559 |
|
snap_->wrapVector(dr); |
1560 |
|
} |
1561 |
< |
cuts = getGroupCutoffs( (*j1), (*j2) ); |
1562 |
< |
if (dr.lengthSquare() < cuts.third) { |
1561 |
> |
getGroupCutoffs( (*j1), (*j2), rcut, rcutsq, rlistsq ); |
1562 |
> |
if (dr.lengthSquare() < rlistsq) { |
1563 |
|
neighborList.push_back(make_pair((*j1), (*j2))); |
1564 |
|
} |
1565 |
|
} |
1579 |
|
if (usePeriodicBoundaryConditions_) { |
1580 |
|
snap_->wrapVector(dr); |
1581 |
|
} |
1582 |
< |
cuts = getGroupCutoffs( j1, j2 ); |
1583 |
< |
if (dr.lengthSquare() < cuts.third) { |
1582 |
> |
getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq); |
1583 |
> |
if (dr.lengthSquare() < rlistsq) { |
1584 |
|
neighborList.push_back(make_pair(j1, j2)); |
1585 |
|
} |
1586 |
|
} |
1594 |
|
if (usePeriodicBoundaryConditions_) { |
1595 |
|
snap_->wrapVector(dr); |
1596 |
|
} |
1597 |
< |
cuts = getGroupCutoffs( j1, j2 ); |
1598 |
< |
if (dr.lengthSquare() < cuts.third) { |
1597 |
> |
getGroupCutoffs( j1, j2, rcut, rcutsq, rlistsq ); |
1598 |
> |
if (dr.lengthSquare() < rlistsq) { |
1599 |
|
neighborList.push_back(make_pair(j1, j2)); |
1600 |
|
} |
1601 |
|
} |
1608 |
|
saved_CG_positions_.clear(); |
1609 |
|
for (int i = 0; i < nGroups_; i++) |
1610 |
|
saved_CG_positions_.push_back(snap_->cgData.position[i]); |
1586 |
– |
|
1587 |
– |
return neighborList; |
1611 |
|
} |
1612 |
|
} //end namespace OpenMD |