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 |