| 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(); |
| 68 |
|
massFactors = info_->getMassFactors(); |
| 69 |
+ |
|
| 70 |
|
PairList excludes = info_->getExcludedInteractions(); |
| 71 |
|
PairList oneTwo = info_->getOneTwoInteractions(); |
| 72 |
|
PairList oneThree = info_->getOneThreeInteractions(); |
| 109 |
|
identsRow.resize(nAtomsInRow_); |
| 110 |
|
identsCol.resize(nAtomsInCol_); |
| 111 |
|
|
| 112 |
< |
AtomCommIntRow->gather(identsLocal, identsRow); |
| 113 |
< |
AtomCommIntColumn->gather(identsLocal, identsCol); |
| 112 |
> |
AtomCommIntRow->gather(idents, identsRow); |
| 113 |
> |
AtomCommIntColumn->gather(idents, identsCol); |
| 114 |
|
|
| 115 |
|
AtomCommIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal); |
| 116 |
|
AtomCommIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal); |
| 233 |
|
set<AtomType*> atypes = info_->getSimulatedAtomTypes(); |
| 234 |
|
vector<RealType> atypeCutoff; |
| 235 |
|
atypeCutoff.resize( atypes.size() ); |
| 236 |
< |
|
| 236 |
> |
|
| 237 |
|
for (set<AtomType*>::iterator at = atypes.begin(); |
| 238 |
|
at != atypes.end(); ++at){ |
| 238 |
– |
rc = interactionMan_->getSuggestedCutoffRadius(*at); |
| 239 |
|
atid = (*at)->getIdent(); |
| 240 |
< |
atypeCutoff[atid] = rc; |
| 240 |
> |
|
| 241 |
> |
if (userChoseCutoff_) |
| 242 |
> |
atypeCutoff[atid] = userCutoff_; |
| 243 |
> |
else |
| 244 |
> |
atypeCutoff[atid] = interactionMan_->getSuggestedCutoffRadius(*at); |
| 245 |
|
} |
| 246 |
|
|
| 247 |
|
vector<RealType> gTypeCutoffs; |
| 313 |
|
for (vector<int>::iterator ia = atomList.begin(); |
| 314 |
|
ia != atomList.end(); ++ia) { |
| 315 |
|
int atom1 = (*ia); |
| 316 |
< |
atid = identsLocal[atom1]; |
| 316 |
> |
atid = idents[atom1]; |
| 317 |
|
if (atypeCutoff[atid] > groupCutoff[cg1]) { |
| 318 |
|
groupCutoff[cg1] = atypeCutoff[atid]; |
| 319 |
|
} |
| 382 |
|
if (abs(gTypeCutoffMap[key].first - userCutoff_) > 0.0001) { |
| 383 |
|
sprintf(painCave.errMsg, |
| 384 |
|
"ForceMatrixDecomposition::createGtypeCutoffMap " |
| 385 |
< |
"user-specified rCut does not match computed group Cutoff\n"); |
| 385 |
> |
"user-specified rCut (%lf) does not match computed group Cutoff\n", userCutoff_); |
| 386 |
|
painCave.severity = OPENMD_ERROR; |
| 387 |
|
painCave.isFatal = 1; |
| 388 |
|
simError(); |
| 414 |
|
} |
| 415 |
|
|
| 416 |
|
void ForceMatrixDecomposition::zeroWorkArrays() { |
| 417 |
< |
|
| 418 |
< |
for (int j = 0; j < N_INTERACTION_FAMILIES; j++) { |
| 415 |
< |
longRangePot_[j] = 0.0; |
| 416 |
< |
} |
| 417 |
> |
pairwisePot = 0.0; |
| 418 |
> |
embeddingPot = 0.0; |
| 419 |
|
|
| 420 |
|
#ifdef IS_MPI |
| 421 |
|
if (storageLayout_ & DataStorage::dslForce) { |
| 432 |
|
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 433 |
|
|
| 434 |
|
fill(pot_col.begin(), pot_col.end(), |
| 435 |
< |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 434 |
< |
|
| 435 |
< |
pot_local = Vector<RealType, N_INTERACTION_FAMILIES>(0.0); |
| 435 |
> |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 436 |
|
|
| 437 |
|
if (storageLayout_ & DataStorage::dslParticlePot) { |
| 438 |
|
fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), 0.0); |
| 606 |
|
AtomCommPotRow->scatter(pot_row, pot_temp); |
| 607 |
|
|
| 608 |
|
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
| 609 |
< |
pot_local += pot_temp[ii]; |
| 609 |
> |
pairwisePot += pot_temp[ii]; |
| 610 |
|
|
| 611 |
|
fill(pot_temp.begin(), pot_temp.end(), |
| 612 |
|
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
| 614 |
|
AtomCommPotColumn->scatter(pot_col, pot_temp); |
| 615 |
|
|
| 616 |
|
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
| 617 |
< |
pot_local += pot_temp[ii]; |
| 618 |
< |
|
| 617 |
> |
pairwisePot += pot_temp[ii]; |
| 618 |
|
#endif |
| 619 |
+ |
|
| 620 |
|
} |
| 621 |
|
|
| 622 |
|
int ForceMatrixDecomposition::getNAtomsInRow() { |
| 691 |
|
#ifdef IS_MPI |
| 692 |
|
return massFactorsRow[atom1]; |
| 693 |
|
#else |
| 694 |
+ |
cerr << "mfs = " << massFactors.size() << " atom1 = " << atom1 << "\n"; |
| 695 |
|
return massFactors[atom1]; |
| 696 |
|
#endif |
| 697 |
|
} |
| 755 |
|
for (vector<int>::iterator i = skipsForAtom[atom1].begin(); |
| 756 |
|
i != skipsForAtom[atom1].end(); ++i) { |
| 757 |
|
if ( (*i) == unique_id_2 ) return true; |
| 758 |
< |
} |
| 758 |
> |
} |
| 759 |
|
|
| 760 |
+ |
return false; |
| 761 |
|
} |
| 762 |
|
|
| 763 |
|
|
| 778 |
|
} |
| 779 |
|
|
| 780 |
|
// filling interaction blocks with pointers |
| 781 |
< |
void ForceMatrixDecomposition::fillInteractionData(InteractionData idat, |
| 781 |
> |
void ForceMatrixDecomposition::fillInteractionData(InteractionData &idat, |
| 782 |
|
int atom1, int atom2) { |
| 783 |
|
#ifdef IS_MPI |
| 784 |
|
|
| 822 |
|
|
| 823 |
|
#else |
| 824 |
|
|
| 825 |
< |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
| 826 |
< |
ff_->getAtomType(identsLocal[atom2]) ); |
| 825 |
> |
idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), |
| 826 |
> |
ff_->getAtomType(idents[atom2]) ); |
| 827 |
|
|
| 828 |
|
if (storageLayout_ & DataStorage::dslAmat) { |
| 829 |
|
idat.A1 = &(snap_->atomData.aMat[atom1]); |
| 840 |
|
idat.t2 = &(snap_->atomData.torque[atom2]); |
| 841 |
|
} |
| 842 |
|
|
| 843 |
< |
if (storageLayout_ & DataStorage::dslDensity) { |
| 843 |
> |
if (storageLayout_ & DataStorage::dslDensity) { |
| 844 |
|
idat.rho1 = &(snap_->atomData.density[atom1]); |
| 845 |
|
idat.rho2 = &(snap_->atomData.density[atom2]); |
| 846 |
|
} |
| 864 |
|
} |
| 865 |
|
|
| 866 |
|
|
| 867 |
< |
void ForceMatrixDecomposition::unpackInteractionData(InteractionData idat, int atom1, int atom2) { |
| 867 |
> |
void ForceMatrixDecomposition::unpackInteractionData(InteractionData &idat, int atom1, int atom2) { |
| 868 |
|
#ifdef IS_MPI |
| 869 |
|
pot_row[atom1] += 0.5 * *(idat.pot); |
| 870 |
|
pot_col[atom2] += 0.5 * *(idat.pot); |
| 872 |
|
atomRowData.force[atom1] += *(idat.f1); |
| 873 |
|
atomColData.force[atom2] -= *(idat.f1); |
| 874 |
|
#else |
| 875 |
< |
longRangePot_ += *(idat.pot); |
| 876 |
< |
|
| 875 |
> |
pairwisePot += *(idat.pot); |
| 876 |
> |
|
| 877 |
|
snap_->atomData.force[atom1] += *(idat.f1); |
| 878 |
|
snap_->atomData.force[atom2] -= *(idat.f1); |
| 879 |
|
#endif |
| 881 |
|
} |
| 882 |
|
|
| 883 |
|
|
| 884 |
< |
void ForceMatrixDecomposition::fillSkipData(InteractionData idat, |
| 884 |
> |
void ForceMatrixDecomposition::fillSkipData(InteractionData &idat, |
| 885 |
|
int atom1, int atom2) { |
| 886 |
|
#ifdef IS_MPI |
| 887 |
|
idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), |
| 891 |
|
idat.eFrame1 = &(atomRowData.electroFrame[atom1]); |
| 892 |
|
idat.eFrame2 = &(atomColData.electroFrame[atom2]); |
| 893 |
|
} |
| 894 |
+ |
|
| 895 |
|
if (storageLayout_ & DataStorage::dslTorque) { |
| 896 |
|
idat.t1 = &(atomRowData.torque[atom1]); |
| 897 |
|
idat.t2 = &(atomColData.torque[atom2]); |
| 898 |
|
} |
| 899 |
+ |
|
| 900 |
+ |
if (storageLayout_ & DataStorage::dslSkippedCharge) { |
| 901 |
+ |
idat.skippedCharge1 = &(atomRowData.skippedCharge[atom1]); |
| 902 |
+ |
idat.skippedCharge2 = &(atomColData.skippedCharge[atom2]); |
| 903 |
+ |
} |
| 904 |
|
#else |
| 905 |
< |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
| 906 |
< |
ff_->getAtomType(identsLocal[atom2]) ); |
| 905 |
> |
idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), |
| 906 |
> |
ff_->getAtomType(idents[atom2]) ); |
| 907 |
|
|
| 908 |
|
if (storageLayout_ & DataStorage::dslElectroFrame) { |
| 909 |
|
idat.eFrame1 = &(snap_->atomData.electroFrame[atom1]); |
| 910 |
|
idat.eFrame2 = &(snap_->atomData.electroFrame[atom2]); |
| 911 |
|
} |
| 912 |
+ |
|
| 913 |
|
if (storageLayout_ & DataStorage::dslTorque) { |
| 914 |
|
idat.t1 = &(snap_->atomData.torque[atom1]); |
| 915 |
|
idat.t2 = &(snap_->atomData.torque[atom2]); |
| 916 |
|
} |
| 917 |
+ |
|
| 918 |
+ |
if (storageLayout_ & DataStorage::dslSkippedCharge) { |
| 919 |
+ |
idat.skippedCharge1 = &(snap_->atomData.skippedCharge[atom1]); |
| 920 |
+ |
idat.skippedCharge2 = &(snap_->atomData.skippedCharge[atom2]); |
| 921 |
+ |
} |
| 922 |
|
#endif |
| 923 |
|
} |
| 924 |
|
|
| 925 |
+ |
|
| 926 |
+ |
void ForceMatrixDecomposition::unpackSkipData(InteractionData &idat, int atom1, int atom2) { |
| 927 |
+ |
#ifdef IS_MPI |
| 928 |
+ |
pot_row[atom1] += 0.5 * *(idat.pot); |
| 929 |
+ |
pot_col[atom2] += 0.5 * *(idat.pot); |
| 930 |
+ |
#else |
| 931 |
+ |
pairwisePot += *(idat.pot); |
| 932 |
+ |
#endif |
| 933 |
+ |
|
| 934 |
+ |
} |
| 935 |
+ |
|
| 936 |
+ |
|
| 937 |
|
/* |
| 938 |
|
* buildNeighborList |
| 939 |
|
* |
| 969 |
|
int cellIndex; |
| 970 |
|
int nCtot = nCells_.x() * nCells_.y() * nCells_.z(); |
| 971 |
|
|
| 946 |
– |
cerr << "flag1\n"; |
| 972 |
|
#ifdef IS_MPI |
| 973 |
|
cellListRow_.resize(nCtot); |
| 974 |
|
cellListCol_.resize(nCtot); |
| 975 |
|
#else |
| 976 |
|
cellList_.resize(nCtot); |
| 977 |
|
#endif |
| 978 |
< |
cerr << "flag2\n"; |
| 978 |
> |
|
| 979 |
|
#ifdef IS_MPI |
| 980 |
|
for (int i = 0; i < nGroupsInRow_; i++) { |
| 981 |
|
rs = cgRowData.position[i]; |
| 1136 |
|
saved_CG_positions_.clear(); |
| 1137 |
|
for (int i = 0; i < nGroups_; i++) |
| 1138 |
|
saved_CG_positions_.push_back(snap_->cgData.position[i]); |
| 1139 |
< |
|
| 1139 |
> |
|
| 1140 |
|
return neighborList; |
| 1141 |
|
} |
| 1142 |
|
} //end namespace OpenMD |