53 |
|
// surrounding cells (not just the 14 upper triangular blocks that |
54 |
|
// are used when the processor can see all pairs) |
55 |
|
#ifdef IS_MPI |
56 |
< |
cellOffsets_.push_back( Vector3i(-1, 0, 0) ); |
57 |
< |
cellOffsets_.push_back( Vector3i(-1,-1, 0) ); |
58 |
< |
cellOffsets_.push_back( Vector3i( 0,-1, 0) ); |
59 |
< |
cellOffsets_.push_back( Vector3i( 1,-1, 0) ); |
60 |
< |
cellOffsets_.push_back( Vector3i( 0, 0,-1) ); |
61 |
< |
cellOffsets_.push_back( Vector3i(-1, 0, 1) ); |
56 |
> |
cellOffsets_.clear(); |
57 |
|
cellOffsets_.push_back( Vector3i(-1,-1,-1) ); |
58 |
|
cellOffsets_.push_back( Vector3i( 0,-1,-1) ); |
59 |
< |
cellOffsets_.push_back( Vector3i( 1,-1,-1) ); |
59 |
> |
cellOffsets_.push_back( Vector3i( 1,-1,-1) ); |
60 |
> |
cellOffsets_.push_back( Vector3i(-1, 0,-1) ); |
61 |
> |
cellOffsets_.push_back( Vector3i( 0, 0,-1) ); |
62 |
|
cellOffsets_.push_back( Vector3i( 1, 0,-1) ); |
66 |
– |
cellOffsets_.push_back( Vector3i( 1, 1,-1) ); |
67 |
– |
cellOffsets_.push_back( Vector3i( 0, 1,-1) ); |
63 |
|
cellOffsets_.push_back( Vector3i(-1, 1,-1) ); |
64 |
+ |
cellOffsets_.push_back( Vector3i( 0, 1,-1) ); |
65 |
+ |
cellOffsets_.push_back( Vector3i( 1, 1,-1) ); |
66 |
+ |
cellOffsets_.push_back( Vector3i(-1,-1, 0) ); |
67 |
+ |
cellOffsets_.push_back( Vector3i( 0,-1, 0) ); |
68 |
+ |
cellOffsets_.push_back( Vector3i( 1,-1, 0) ); |
69 |
+ |
cellOffsets_.push_back( Vector3i(-1, 0, 0) ); |
70 |
+ |
cellOffsets_.push_back( Vector3i( 0, 0, 0) ); |
71 |
+ |
cellOffsets_.push_back( Vector3i( 1, 0, 0) ); |
72 |
+ |
cellOffsets_.push_back( Vector3i(-1, 1, 0) ); |
73 |
+ |
cellOffsets_.push_back( Vector3i( 0, 1, 0) ); |
74 |
+ |
cellOffsets_.push_back( Vector3i( 1, 1, 0) ); |
75 |
+ |
cellOffsets_.push_back( Vector3i(-1,-1, 1) ); |
76 |
+ |
cellOffsets_.push_back( Vector3i( 0,-1, 1) ); |
77 |
+ |
cellOffsets_.push_back( Vector3i( 1,-1, 1) ); |
78 |
+ |
cellOffsets_.push_back( Vector3i(-1, 0, 1) ); |
79 |
+ |
cellOffsets_.push_back( Vector3i( 0, 0, 1) ); |
80 |
+ |
cellOffsets_.push_back( Vector3i( 1, 0, 1) ); |
81 |
+ |
cellOffsets_.push_back( Vector3i(-1, 1, 1) ); |
82 |
+ |
cellOffsets_.push_back( Vector3i( 0, 1, 1) ); |
83 |
+ |
cellOffsets_.push_back( Vector3i( 1, 1, 1) ); |
84 |
|
#endif |
85 |
|
} |
86 |
|
|
168 |
|
AtomColToGlobal.resize(nAtomsInCol_); |
169 |
|
AtomPlanIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal); |
170 |
|
AtomPlanIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal); |
156 |
– |
|
157 |
– |
cerr << "Atoms in Local:\n"; |
158 |
– |
for (int i = 0; i < AtomLocalToGlobal.size(); i++) { |
159 |
– |
cerr << "i =\t" << i << "\t localAt =\t" << AtomLocalToGlobal[i] << "\n"; |
160 |
– |
} |
161 |
– |
cerr << "Atoms in Row:\n"; |
162 |
– |
for (int i = 0; i < AtomRowToGlobal.size(); i++) { |
163 |
– |
cerr << "i =\t" << i << "\t rowAt =\t" << AtomRowToGlobal[i] << "\n"; |
164 |
– |
} |
165 |
– |
cerr << "Atoms in Col:\n"; |
166 |
– |
for (int i = 0; i < AtomColToGlobal.size(); i++) { |
167 |
– |
cerr << "i =\t" << i << "\t colAt =\t" << AtomColToGlobal[i] << "\n"; |
168 |
– |
} |
171 |
|
|
172 |
|
cgRowToGlobal.resize(nGroupsInRow_); |
173 |
|
cgColToGlobal.resize(nGroupsInCol_); |
174 |
|
cgPlanIntRow->gather(cgLocalToGlobal, cgRowToGlobal); |
175 |
|
cgPlanIntColumn->gather(cgLocalToGlobal, cgColToGlobal); |
176 |
|
|
175 |
– |
cerr << "Gruops in Local:\n"; |
176 |
– |
for (int i = 0; i < cgLocalToGlobal.size(); i++) { |
177 |
– |
cerr << "i =\t" << i << "\t localCG =\t" << cgLocalToGlobal[i] << "\n"; |
178 |
– |
} |
179 |
– |
cerr << "Groups in Row:\n"; |
180 |
– |
for (int i = 0; i < cgRowToGlobal.size(); i++) { |
181 |
– |
cerr << "i =\t" << i << "\t rowCG =\t" << cgRowToGlobal[i] << "\n"; |
182 |
– |
} |
183 |
– |
cerr << "Groups in Col:\n"; |
184 |
– |
for (int i = 0; i < cgColToGlobal.size(); i++) { |
185 |
– |
cerr << "i =\t" << i << "\t colCG =\t" << cgColToGlobal[i] << "\n"; |
186 |
– |
} |
187 |
– |
|
188 |
– |
|
177 |
|
massFactorsRow.resize(nAtomsInRow_); |
178 |
|
massFactorsCol.resize(nAtomsInCol_); |
179 |
|
AtomPlanRealRow->gather(massFactors, massFactorsRow); |
233 |
|
} |
234 |
|
} |
235 |
|
|
236 |
< |
#endif |
249 |
< |
|
250 |
< |
// allocate memory for the parallel objects |
251 |
< |
atypesLocal.resize(nLocal_); |
252 |
< |
|
253 |
< |
for (int i = 0; i < nLocal_; i++) |
254 |
< |
atypesLocal[i] = ff_->getAtomType(idents[i]); |
255 |
< |
|
256 |
< |
groupList_.clear(); |
257 |
< |
groupList_.resize(nGroups_); |
258 |
< |
for (int i = 0; i < nGroups_; i++) { |
259 |
< |
int gid = cgLocalToGlobal[i]; |
260 |
< |
for (int j = 0; j < nLocal_; j++) { |
261 |
< |
int aid = AtomLocalToGlobal[j]; |
262 |
< |
if (globalGroupMembership[aid] == gid) { |
263 |
< |
groupList_[i].push_back(j); |
264 |
< |
} |
265 |
< |
} |
266 |
< |
} |
267 |
< |
|
236 |
> |
#else |
237 |
|
excludesForAtom.clear(); |
238 |
|
excludesForAtom.resize(nLocal_); |
239 |
|
toposForAtom.clear(); |
247 |
|
for (int j = 0; j < nLocal_; j++) { |
248 |
|
int jglob = AtomLocalToGlobal[j]; |
249 |
|
|
250 |
< |
if (excludes->hasPair(iglob, jglob)) |
250 |
> |
if (excludes->hasPair(iglob, jglob)) |
251 |
|
excludesForAtom[i].push_back(j); |
252 |
|
|
253 |
+ |
|
254 |
|
if (oneTwo->hasPair(iglob, jglob)) { |
255 |
|
toposForAtom[i].push_back(j); |
256 |
|
topoDist[i].push_back(1); |
267 |
|
} |
268 |
|
} |
269 |
|
} |
270 |
< |
|
270 |
> |
#endif |
271 |
> |
|
272 |
> |
// allocate memory for the parallel objects |
273 |
> |
atypesLocal.resize(nLocal_); |
274 |
> |
|
275 |
> |
for (int i = 0; i < nLocal_; i++) |
276 |
> |
atypesLocal[i] = ff_->getAtomType(idents[i]); |
277 |
> |
|
278 |
> |
groupList_.clear(); |
279 |
> |
groupList_.resize(nGroups_); |
280 |
> |
for (int i = 0; i < nGroups_; i++) { |
281 |
> |
int gid = cgLocalToGlobal[i]; |
282 |
> |
for (int j = 0; j < nLocal_; j++) { |
283 |
> |
int aid = AtomLocalToGlobal[j]; |
284 |
> |
if (globalGroupMembership[aid] == gid) { |
285 |
> |
groupList_[i].push_back(j); |
286 |
> |
} |
287 |
> |
} |
288 |
> |
} |
289 |
> |
|
290 |
> |
|
291 |
|
createGtypeCutoffMap(); |
292 |
|
|
293 |
|
} |
566 |
|
|
567 |
|
// gather up the cutoff group positions |
568 |
|
|
579 |
– |
cerr << "before gather\n"; |
580 |
– |
for (int i = 0; i < snap_->cgData.position.size(); i++) { |
581 |
– |
cerr << "cgpos = " << snap_->cgData.position[i] << "\n"; |
582 |
– |
} |
583 |
– |
|
569 |
|
cgPlanVectorRow->gather(snap_->cgData.position, |
570 |
|
cgRowData.position); |
571 |
|
|
587 |
– |
cerr << "after gather\n"; |
588 |
– |
for (int i = 0; i < cgRowData.position.size(); i++) { |
589 |
– |
cerr << "cgRpos = " << cgRowData.position[i] << "\n"; |
590 |
– |
} |
591 |
– |
|
572 |
|
cgPlanVectorColumn->gather(snap_->cgData.position, |
573 |
|
cgColData.position); |
594 |
– |
for (int i = 0; i < cgColData.position.size(); i++) { |
595 |
– |
cerr << "cgCpos = " << cgColData.position[i] << "\n"; |
596 |
– |
} |
574 |
|
|
575 |
|
|
576 |
|
// if needed, gather the atomic rotation matrices |
685 |
|
} |
686 |
|
|
687 |
|
AtomPlanRealColumn->scatter(atomColData.skippedCharge, skch_tmp); |
688 |
< |
for (int i = 0; i < ns; i++) |
688 |
> |
for (int i = 0; i < ns; i++) |
689 |
|
snap_->atomData.skippedCharge[i] += skch_tmp[i]; |
690 |
+ |
|
691 |
|
} |
692 |
|
|
693 |
|
nLocal_ = snap_->getNumberOfAtoms(); |
709 |
|
|
710 |
|
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
711 |
|
pairwisePot += pot_temp[ii]; |
712 |
+ |
|
713 |
+ |
for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
714 |
+ |
RealType ploc1 = pairwisePot[ii]; |
715 |
+ |
RealType ploc2 = 0.0; |
716 |
+ |
MPI::COMM_WORLD.Allreduce(&ploc1, &ploc2, 1, MPI::REALTYPE, MPI::SUM); |
717 |
+ |
pairwisePot[ii] = ploc2; |
718 |
+ |
} |
719 |
+ |
|
720 |
+ |
for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
721 |
+ |
RealType ploc1 = embeddingPot[ii]; |
722 |
+ |
RealType ploc2 = 0.0; |
723 |
+ |
MPI::COMM_WORLD.Allreduce(&ploc1, &ploc2, 1, MPI::REALTYPE, MPI::SUM); |
724 |
+ |
embeddingPot[ii] = ploc2; |
725 |
+ |
} |
726 |
+ |
|
727 |
|
#endif |
728 |
|
|
736 |
– |
cerr << "pairwisePot = " << pairwisePot << "\n"; |
729 |
|
} |
730 |
|
|
731 |
|
int ForceMatrixDecomposition::getNAtomsInRow() { |
760 |
|
|
761 |
|
#ifdef IS_MPI |
762 |
|
d = cgColData.position[cg2] - cgRowData.position[cg1]; |
771 |
– |
cerr << "cg1 = " << cg1 << "\tcg1p = " << cgRowData.position[cg1] << "\n"; |
772 |
– |
cerr << "cg2 = " << cg2 << "\tcg2p = " << cgColData.position[cg2] << "\n"; |
763 |
|
#else |
764 |
|
d = snap_->cgData.position[cg2] - snap_->cgData.position[cg1]; |
775 |
– |
cerr << "cg1 = " << cg1 << "\tcg1p = " << snap_->cgData.position[cg1] << "\n"; |
776 |
– |
cerr << "cg2 = " << cg2 << "\tcg2p = " << snap_->cgData.position[cg2] << "\n"; |
765 |
|
#endif |
766 |
|
|
767 |
|
snap_->wrapVector(d); |
837 |
|
bool ForceMatrixDecomposition::skipAtomPair(int atom1, int atom2) { |
838 |
|
int unique_id_1, unique_id_2; |
839 |
|
|
852 |
– |
|
853 |
– |
cerr << "sap with atom1, atom2 =\t" << atom1 << "\t" << atom2 << "\n"; |
840 |
|
#ifdef IS_MPI |
841 |
|
// in MPI, we have to look up the unique IDs for each atom |
842 |
|
unique_id_1 = AtomRowToGlobal[atom1]; |
843 |
|
unique_id_2 = AtomColToGlobal[atom2]; |
844 |
|
|
859 |
– |
cerr << "sap with uid1, uid2 =\t" << unique_id_1 << "\t" << unique_id_2 << "\n"; |
845 |
|
// this situation should only arise in MPI simulations |
846 |
|
if (unique_id_1 == unique_id_2) return true; |
847 |
|
|
865 |
|
* field) must still be handled for these pairs. |
866 |
|
*/ |
867 |
|
bool ForceMatrixDecomposition::excludeAtomPair(int atom1, int atom2) { |
868 |
< |
int unique_id_2; |
869 |
< |
#ifdef IS_MPI |
870 |
< |
// in MPI, we have to look up the unique IDs for the row atom. |
886 |
< |
unique_id_2 = AtomColToGlobal[atom2]; |
887 |
< |
#else |
888 |
< |
// in the normal loop, the atom numbers are unique |
889 |
< |
unique_id_2 = atom2; |
890 |
< |
#endif |
868 |
> |
|
869 |
> |
// excludesForAtom was constructed to use row/column indices in the MPI |
870 |
> |
// version, and to use local IDs in the non-MPI version: |
871 |
|
|
872 |
|
for (vector<int>::iterator i = excludesForAtom[atom1].begin(); |
873 |
|
i != excludesForAtom[atom1].end(); ++i) { |
874 |
< |
if ( (*i) == unique_id_2 ) return true; |
874 |
> |
if ( (*i) == atom2 ) return true; |
875 |
|
} |
876 |
|
|
877 |
|
return false; |
1111 |
|
// add this cutoff group to the list of groups in this cell; |
1112 |
|
cellListCol_[cellIndex].push_back(i); |
1113 |
|
} |
1114 |
+ |
|
1115 |
|
#else |
1116 |
|
for (int i = 0; i < nGroups_; i++) { |
1117 |
|
rs = snap_->cgData.position[i]; |
1137 |
|
// add this cutoff group to the list of groups in this cell; |
1138 |
|
cellList_[cellIndex].push_back(i); |
1139 |
|
} |
1140 |
+ |
|
1141 |
|
#endif |
1142 |
|
|
1143 |
|
for (int m1z = 0; m1z < nCells_.z(); m1z++) { |
1150 |
|
os != cellOffsets_.end(); ++os) { |
1151 |
|
|
1152 |
|
Vector3i m2v = m1v + (*os); |
1153 |
< |
|
1153 |
> |
|
1154 |
> |
|
1155 |
|
if (m2v.x() >= nCells_.x()) { |
1156 |
|
m2v.x() = 0; |
1157 |
|
} else if (m2v.x() < 0) { |
1169 |
|
} else if (m2v.z() < 0) { |
1170 |
|
m2v.z() = nCells_.z() - 1; |
1171 |
|
} |
1172 |
< |
|
1172 |
> |
|
1173 |
|
int m2 = Vlinear (m2v, nCells_); |
1174 |
|
|
1175 |
|
#ifdef IS_MPI |
1178 |
|
for (vector<int>::iterator j2 = cellListCol_[m2].begin(); |
1179 |
|
j2 != cellListCol_[m2].end(); ++j2) { |
1180 |
|
|
1181 |
< |
// In parallel, we need to visit *all* pairs of row & |
1182 |
< |
// column indicies and will truncate later on. |
1181 |
> |
// In parallel, we need to visit *all* pairs of row |
1182 |
> |
// & column indicies and will divide labor in the |
1183 |
> |
// force evaluation later. |
1184 |
|
dr = cgColData.position[(*j2)] - cgRowData.position[(*j1)]; |
1185 |
|
snap_->wrapVector(dr); |
1186 |
|
cuts = getGroupCutoffs( (*j1), (*j2) ); |