| 61 |
|
nGroups_ = info_->getNLocalCutoffGroups(); |
| 62 |
|
cerr << "in dId, nGroups = " << nGroups_ << "\n"; |
| 63 |
|
// gather the information for atomtype IDs (atids): |
| 64 |
< |
identsLocal = info_->getIdentArray(); |
| 64 |
> |
idents = info_->getIdentArray(); |
| 65 |
|
AtomLocalToGlobal = info_->getGlobalAtomIndices(); |
| 66 |
|
cgLocalToGlobal = info_->getGlobalGroupIndices(); |
| 67 |
|
vector<int> globalGroupMembership = info_->getGlobalGroupMembership(); |
| 108 |
|
identsRow.resize(nAtomsInRow_); |
| 109 |
|
identsCol.resize(nAtomsInCol_); |
| 110 |
|
|
| 111 |
< |
AtomCommIntRow->gather(identsLocal, identsRow); |
| 112 |
< |
AtomCommIntColumn->gather(identsLocal, identsCol); |
| 111 |
> |
AtomCommIntRow->gather(idents, identsRow); |
| 112 |
> |
AtomCommIntColumn->gather(idents, identsCol); |
| 113 |
|
|
| 114 |
|
AtomCommIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal); |
| 115 |
|
AtomCommIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal); |
| 232 |
|
set<AtomType*> atypes = info_->getSimulatedAtomTypes(); |
| 233 |
|
vector<RealType> atypeCutoff; |
| 234 |
|
atypeCutoff.resize( atypes.size() ); |
| 235 |
< |
|
| 235 |
> |
|
| 236 |
|
for (set<AtomType*>::iterator at = atypes.begin(); |
| 237 |
|
at != atypes.end(); ++at){ |
| 238 |
– |
rc = interactionMan_->getSuggestedCutoffRadius(*at); |
| 238 |
|
atid = (*at)->getIdent(); |
| 239 |
< |
atypeCutoff[atid] = rc; |
| 239 |
> |
|
| 240 |
> |
if (userChoseCutoff_) |
| 241 |
> |
atypeCutoff[atid] = userCutoff_; |
| 242 |
> |
else |
| 243 |
> |
atypeCutoff[atid] = interactionMan_->getSuggestedCutoffRadius(*at); |
| 244 |
|
} |
| 245 |
|
|
| 246 |
|
vector<RealType> gTypeCutoffs; |
| 312 |
|
for (vector<int>::iterator ia = atomList.begin(); |
| 313 |
|
ia != atomList.end(); ++ia) { |
| 314 |
|
int atom1 = (*ia); |
| 315 |
< |
atid = identsLocal[atom1]; |
| 315 |
> |
atid = idents[atom1]; |
| 316 |
|
if (atypeCutoff[atid] > groupCutoff[cg1]) { |
| 317 |
|
groupCutoff[cg1] = atypeCutoff[atid]; |
| 318 |
|
} |
| 381 |
|
if (abs(gTypeCutoffMap[key].first - userCutoff_) > 0.0001) { |
| 382 |
|
sprintf(painCave.errMsg, |
| 383 |
|
"ForceMatrixDecomposition::createGtypeCutoffMap " |
| 384 |
< |
"user-specified rCut does not match computed group Cutoff\n"); |
| 384 |
> |
"user-specified rCut (%lf) does not match computed group Cutoff\n", userCutoff_); |
| 385 |
|
painCave.severity = OPENMD_ERROR; |
| 386 |
|
painCave.isFatal = 1; |
| 387 |
|
simError(); |
| 413 |
|
} |
| 414 |
|
|
| 415 |
|
void ForceMatrixDecomposition::zeroWorkArrays() { |
| 416 |
< |
|
| 417 |
< |
for (int j = 0; j < N_INTERACTION_FAMILIES; j++) { |
| 415 |
< |
longRangePot_[j] = 0.0; |
| 416 |
< |
} |
| 416 |
> |
pairwisePot = 0.0; |
| 417 |
> |
embeddingPot = 0.0; |
| 418 |
|
|
| 419 |
|
#ifdef IS_MPI |
| 420 |
|
if (storageLayout_ & DataStorage::dslForce) { |
| 431 |
|
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 432 |
|
|
| 433 |
|
fill(pot_col.begin(), pot_col.end(), |
| 434 |
< |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 434 |
< |
|
| 435 |
< |
pot_local = Vector<RealType, N_INTERACTION_FAMILIES>(0.0); |
| 434 |
> |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 435 |
|
|
| 436 |
|
if (storageLayout_ & DataStorage::dslParticlePot) { |
| 437 |
|
fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), 0.0); |
| 605 |
|
AtomCommPotRow->scatter(pot_row, pot_temp); |
| 606 |
|
|
| 607 |
|
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
| 608 |
< |
pot_local += pot_temp[ii]; |
| 608 |
> |
pairwisePot += pot_temp[ii]; |
| 609 |
|
|
| 610 |
|
fill(pot_temp.begin(), pot_temp.end(), |
| 611 |
|
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 613 |
|
AtomCommPotColumn->scatter(pot_col, pot_temp); |
| 614 |
|
|
| 615 |
|
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
| 616 |
< |
pot_local += pot_temp[ii]; |
| 618 |
< |
|
| 616 |
> |
pairwisePot += pot_temp[ii]; |
| 617 |
|
#endif |
| 618 |
+ |
|
| 619 |
|
} |
| 620 |
|
|
| 621 |
|
int ForceMatrixDecomposition::getNAtomsInRow() { |
| 753 |
|
for (vector<int>::iterator i = skipsForAtom[atom1].begin(); |
| 754 |
|
i != skipsForAtom[atom1].end(); ++i) { |
| 755 |
|
if ( (*i) == unique_id_2 ) return true; |
| 756 |
< |
} |
| 756 |
> |
} |
| 757 |
|
|
| 758 |
+ |
return false; |
| 759 |
|
} |
| 760 |
|
|
| 761 |
|
|
| 776 |
|
} |
| 777 |
|
|
| 778 |
|
// filling interaction blocks with pointers |
| 779 |
< |
void ForceMatrixDecomposition::fillInteractionData(InteractionData idat, |
| 779 |
> |
void ForceMatrixDecomposition::fillInteractionData(InteractionData &idat, |
| 780 |
|
int atom1, int atom2) { |
| 781 |
|
#ifdef IS_MPI |
| 782 |
|
|
| 820 |
|
|
| 821 |
|
#else |
| 822 |
|
|
| 823 |
< |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
| 824 |
< |
ff_->getAtomType(identsLocal[atom2]) ); |
| 823 |
> |
idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), |
| 824 |
> |
ff_->getAtomType(idents[atom2]) ); |
| 825 |
|
|
| 826 |
|
if (storageLayout_ & DataStorage::dslAmat) { |
| 827 |
|
idat.A1 = &(snap_->atomData.aMat[atom1]); |
| 838 |
|
idat.t2 = &(snap_->atomData.torque[atom2]); |
| 839 |
|
} |
| 840 |
|
|
| 841 |
< |
if (storageLayout_ & DataStorage::dslDensity) { |
| 841 |
> |
if (storageLayout_ & DataStorage::dslDensity) { |
| 842 |
|
idat.rho1 = &(snap_->atomData.density[atom1]); |
| 843 |
|
idat.rho2 = &(snap_->atomData.density[atom2]); |
| 844 |
|
} |
| 862 |
|
} |
| 863 |
|
|
| 864 |
|
|
| 865 |
< |
void ForceMatrixDecomposition::unpackInteractionData(InteractionData idat, int atom1, int atom2) { |
| 865 |
> |
void ForceMatrixDecomposition::unpackInteractionData(InteractionData &idat, int atom1, int atom2) { |
| 866 |
|
#ifdef IS_MPI |
| 867 |
|
pot_row[atom1] += 0.5 * *(idat.pot); |
| 868 |
|
pot_col[atom2] += 0.5 * *(idat.pot); |
| 870 |
|
atomRowData.force[atom1] += *(idat.f1); |
| 871 |
|
atomColData.force[atom2] -= *(idat.f1); |
| 872 |
|
#else |
| 873 |
< |
longRangePot_ += *(idat.pot); |
| 874 |
< |
|
| 873 |
> |
pairwisePot += *(idat.pot); |
| 874 |
> |
|
| 875 |
|
snap_->atomData.force[atom1] += *(idat.f1); |
| 876 |
|
snap_->atomData.force[atom2] -= *(idat.f1); |
| 877 |
|
#endif |
| 879 |
|
} |
| 880 |
|
|
| 881 |
|
|
| 882 |
< |
void ForceMatrixDecomposition::fillSkipData(InteractionData idat, |
| 882 |
> |
void ForceMatrixDecomposition::fillSkipData(InteractionData &idat, |
| 883 |
|
int atom1, int atom2) { |
| 884 |
+ |
// Still Missing:: skippedCharge fill must be added to DataStorage |
| 885 |
|
#ifdef IS_MPI |
| 886 |
|
idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), |
| 887 |
|
ff_->getAtomType(identsCol[atom2]) ); |
| 895 |
|
idat.t2 = &(atomColData.torque[atom2]); |
| 896 |
|
} |
| 897 |
|
#else |
| 898 |
< |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
| 899 |
< |
ff_->getAtomType(identsLocal[atom2]) ); |
| 898 |
> |
idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), |
| 899 |
> |
ff_->getAtomType(idents[atom2]) ); |
| 900 |
|
|
| 901 |
|
if (storageLayout_ & DataStorage::dslElectroFrame) { |
| 902 |
|
idat.eFrame1 = &(snap_->atomData.electroFrame[atom1]); |
| 909 |
|
#endif |
| 910 |
|
} |
| 911 |
|
|
| 912 |
+ |
|
| 913 |
+ |
void ForceMatrixDecomposition::unpackSkipData(InteractionData &idat, int atom1, int atom2) { |
| 914 |
+ |
#ifdef IS_MPI |
| 915 |
+ |
pot_row[atom1] += 0.5 * *(idat.pot); |
| 916 |
+ |
pot_col[atom2] += 0.5 * *(idat.pot); |
| 917 |
+ |
#else |
| 918 |
+ |
pairwisePot += *(idat.pot); |
| 919 |
+ |
#endif |
| 920 |
+ |
|
| 921 |
+ |
} |
| 922 |
+ |
|
| 923 |
+ |
|
| 924 |
|
/* |
| 925 |
|
* buildNeighborList |
| 926 |
|
* |
| 956 |
|
int cellIndex; |
| 957 |
|
int nCtot = nCells_.x() * nCells_.y() * nCells_.z(); |
| 958 |
|
|
| 946 |
– |
cerr << "flag1\n"; |
| 959 |
|
#ifdef IS_MPI |
| 960 |
|
cellListRow_.resize(nCtot); |
| 961 |
|
cellListCol_.resize(nCtot); |
| 962 |
|
#else |
| 963 |
|
cellList_.resize(nCtot); |
| 964 |
|
#endif |
| 965 |
< |
cerr << "flag2\n"; |
| 965 |
> |
|
| 966 |
|
#ifdef IS_MPI |
| 967 |
|
for (int i = 0; i < nGroupsInRow_; i++) { |
| 968 |
|
rs = cgRowData.position[i]; |
| 1123 |
|
saved_CG_positions_.clear(); |
| 1124 |
|
for (int i = 0; i < nGroups_; i++) |
| 1125 |
|
saved_CG_positions_.push_back(snap_->cgData.position[i]); |
| 1126 |
< |
|
| 1126 |
> |
|
| 1127 |
|
return neighborList; |
| 1128 |
|
} |
| 1129 |
|
} //end namespace OpenMD |