42 |
|
#include "math/SquareMatrix3.hpp" |
43 |
|
#include "nonbonded/NonBondedInteraction.hpp" |
44 |
|
#include "brains/SnapshotManager.hpp" |
45 |
+ |
#include "brains/PairList.hpp" |
46 |
|
|
47 |
|
using namespace std; |
48 |
|
namespace OpenMD { |
55 |
|
void ForceMatrixDecomposition::distributeInitialData() { |
56 |
|
snap_ = sman_->getCurrentSnapshot(); |
57 |
|
storageLayout_ = sman_->getStorageLayout(); |
58 |
+ |
ff_ = info_->getForceField(); |
59 |
|
nLocal_ = snap_->getNumberOfAtoms(); |
60 |
|
nGroups_ = snap_->getNumberOfCutoffGroups(); |
61 |
|
|
62 |
+ |
// gather the information for atomtype IDs (atids): |
63 |
+ |
identsLocal = info_->getIdentArray(); |
64 |
+ |
AtomLocalToGlobal = info_->getGlobalAtomIndices(); |
65 |
+ |
cgLocalToGlobal = info_->getGlobalGroupIndices(); |
66 |
+ |
vector<int> globalGroupMembership = info_->getGlobalGroupMembership(); |
67 |
+ |
vector<RealType> massFactorsLocal = info_->getMassFactors(); |
68 |
+ |
PairList excludes = info_->getExcludedInteractions(); |
69 |
+ |
PairList oneTwo = info_->getOneTwoInteractions(); |
70 |
+ |
PairList oneThree = info_->getOneThreeInteractions(); |
71 |
+ |
PairList oneFour = info_->getOneFourInteractions(); |
72 |
+ |
|
73 |
|
#ifdef IS_MPI |
74 |
|
|
75 |
|
AtomCommIntRow = new Communicator<Row,int>(nLocal_); |
76 |
|
AtomCommRealRow = new Communicator<Row,RealType>(nLocal_); |
77 |
|
AtomCommVectorRow = new Communicator<Row,Vector3d>(nLocal_); |
78 |
|
AtomCommMatrixRow = new Communicator<Row,Mat3x3d>(nLocal_); |
79 |
+ |
AtomCommPotRow = new Communicator<Row,potVec>(nLocal_); |
80 |
|
|
81 |
|
AtomCommIntColumn = new Communicator<Column,int>(nLocal_); |
82 |
|
AtomCommRealColumn = new Communicator<Column,RealType>(nLocal_); |
83 |
|
AtomCommVectorColumn = new Communicator<Column,Vector3d>(nLocal_); |
84 |
|
AtomCommMatrixColumn = new Communicator<Column,Mat3x3d>(nLocal_); |
85 |
+ |
AtomCommPotColumn = new Communicator<Column,potVec>(nLocal_); |
86 |
|
|
87 |
|
cgCommIntRow = new Communicator<Row,int>(nGroups_); |
88 |
|
cgCommVectorRow = new Communicator<Row,Vector3d>(nGroups_); |
103 |
|
cgRowData.setStorageLayout(DataStorage::dslPosition); |
104 |
|
cgColData.resize(nGroupsInCol_); |
105 |
|
cgColData.setStorageLayout(DataStorage::dslPosition); |
106 |
< |
|
92 |
< |
vector<vector<RealType> > pot_row(N_INTERACTION_FAMILIES, |
93 |
< |
vector<RealType> (nAtomsInRow_, 0.0)); |
94 |
< |
vector<vector<RealType> > pot_col(N_INTERACTION_FAMILIES, |
95 |
< |
vector<RealType> (nAtomsInCol_, 0.0)); |
96 |
< |
|
97 |
< |
|
98 |
< |
vector<RealType> pot_local(N_INTERACTION_FAMILIES, 0.0); |
99 |
< |
|
100 |
< |
// gather the information for atomtype IDs (atids): |
101 |
< |
vector<int> identsLocal = info_->getIdentArray(); |
106 |
> |
|
107 |
|
identsRow.reserve(nAtomsInRow_); |
108 |
|
identsCol.reserve(nAtomsInCol_); |
109 |
|
|
110 |
|
AtomCommIntRow->gather(identsLocal, identsRow); |
111 |
|
AtomCommIntColumn->gather(identsLocal, identsCol); |
112 |
|
|
108 |
– |
AtomLocalToGlobal = info_->getGlobalAtomIndices(); |
113 |
|
AtomCommIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal); |
114 |
|
AtomCommIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal); |
115 |
|
|
112 |
– |
cgLocalToGlobal = info_->getGlobalGroupIndices(); |
116 |
|
cgCommIntRow->gather(cgLocalToGlobal, cgRowToGlobal); |
117 |
|
cgCommIntColumn->gather(cgLocalToGlobal, cgColToGlobal); |
118 |
|
|
119 |
< |
// still need: |
120 |
< |
// topoDist |
121 |
< |
// exclude |
119 |
> |
AtomCommRealRow->gather(massFactorsLocal, massFactorsRow); |
120 |
> |
AtomCommRealColumn->gather(massFactorsLocal, massFactorsCol); |
121 |
> |
|
122 |
> |
groupListRow_.clear(); |
123 |
> |
groupListRow_.reserve(nGroupsInRow_); |
124 |
> |
for (int i = 0; i < nGroupsInRow_; i++) { |
125 |
> |
int gid = cgRowToGlobal[i]; |
126 |
> |
for (int j = 0; j < nAtomsInRow_; j++) { |
127 |
> |
int aid = AtomRowToGlobal[j]; |
128 |
> |
if (globalGroupMembership[aid] == gid) |
129 |
> |
groupListRow_[i].push_back(j); |
130 |
> |
} |
131 |
> |
} |
132 |
> |
|
133 |
> |
groupListCol_.clear(); |
134 |
> |
groupListCol_.reserve(nGroupsInCol_); |
135 |
> |
for (int i = 0; i < nGroupsInCol_; i++) { |
136 |
> |
int gid = cgColToGlobal[i]; |
137 |
> |
for (int j = 0; j < nAtomsInCol_; j++) { |
138 |
> |
int aid = AtomColToGlobal[j]; |
139 |
> |
if (globalGroupMembership[aid] == gid) |
140 |
> |
groupListCol_[i].push_back(j); |
141 |
> |
} |
142 |
> |
} |
143 |
> |
|
144 |
> |
skipsForRowAtom.clear(); |
145 |
> |
skipsForRowAtom.reserve(nAtomsInRow_); |
146 |
> |
for (int i = 0; i < nAtomsInRow_; i++) { |
147 |
> |
int iglob = AtomRowToGlobal[i]; |
148 |
> |
for (int j = 0; j < nAtomsInCol_; j++) { |
149 |
> |
int jglob = AtomColToGlobal[j]; |
150 |
> |
if (excludes.hasPair(iglob, jglob)) |
151 |
> |
skipsForRowAtom[i].push_back(j); |
152 |
> |
} |
153 |
> |
} |
154 |
> |
|
155 |
> |
toposForRowAtom.clear(); |
156 |
> |
toposForRowAtom.reserve(nAtomsInRow_); |
157 |
> |
for (int i = 0; i < nAtomsInRow_; i++) { |
158 |
> |
int iglob = AtomRowToGlobal[i]; |
159 |
> |
int nTopos = 0; |
160 |
> |
for (int j = 0; j < nAtomsInCol_; j++) { |
161 |
> |
int jglob = AtomColToGlobal[j]; |
162 |
> |
if (oneTwo.hasPair(iglob, jglob)) { |
163 |
> |
toposForRowAtom[i].push_back(j); |
164 |
> |
topoDistRow[i][nTopos] = 1; |
165 |
> |
nTopos++; |
166 |
> |
} |
167 |
> |
if (oneThree.hasPair(iglob, jglob)) { |
168 |
> |
toposForRowAtom[i].push_back(j); |
169 |
> |
topoDistRow[i][nTopos] = 2; |
170 |
> |
nTopos++; |
171 |
> |
} |
172 |
> |
if (oneFour.hasPair(iglob, jglob)) { |
173 |
> |
toposForRowAtom[i].push_back(j); |
174 |
> |
topoDistRow[i][nTopos] = 3; |
175 |
> |
nTopos++; |
176 |
> |
} |
177 |
> |
} |
178 |
> |
} |
179 |
> |
|
180 |
|
#endif |
181 |
+ |
|
182 |
+ |
groupList_.clear(); |
183 |
+ |
groupList_.reserve(nGroups_); |
184 |
+ |
for (int i = 0; i < nGroups_; i++) { |
185 |
+ |
int gid = cgLocalToGlobal[i]; |
186 |
+ |
for (int j = 0; j < nLocal_; j++) { |
187 |
+ |
int aid = AtomLocalToGlobal[j]; |
188 |
+ |
if (globalGroupMembership[aid] == gid) |
189 |
+ |
groupList_[i].push_back(j); |
190 |
+ |
} |
191 |
+ |
} |
192 |
+ |
|
193 |
+ |
skipsForLocalAtom.clear(); |
194 |
+ |
skipsForLocalAtom.reserve(nLocal_); |
195 |
+ |
|
196 |
+ |
for (int i = 0; i < nLocal_; i++) { |
197 |
+ |
int iglob = AtomLocalToGlobal[i]; |
198 |
+ |
for (int j = 0; j < nLocal_; j++) { |
199 |
+ |
int jglob = AtomLocalToGlobal[j]; |
200 |
+ |
if (excludes.hasPair(iglob, jglob)) |
201 |
+ |
skipsForLocalAtom[i].push_back(j); |
202 |
+ |
} |
203 |
+ |
} |
204 |
+ |
|
205 |
+ |
toposForLocalAtom.clear(); |
206 |
+ |
toposForLocalAtom.reserve(nLocal_); |
207 |
+ |
for (int i = 0; i < nLocal_; i++) { |
208 |
+ |
int iglob = AtomLocalToGlobal[i]; |
209 |
+ |
int nTopos = 0; |
210 |
+ |
for (int j = 0; j < nLocal_; j++) { |
211 |
+ |
int jglob = AtomLocalToGlobal[j]; |
212 |
+ |
if (oneTwo.hasPair(iglob, jglob)) { |
213 |
+ |
toposForLocalAtom[i].push_back(j); |
214 |
+ |
topoDistLocal[i][nTopos] = 1; |
215 |
+ |
nTopos++; |
216 |
+ |
} |
217 |
+ |
if (oneThree.hasPair(iglob, jglob)) { |
218 |
+ |
toposForLocalAtom[i].push_back(j); |
219 |
+ |
topoDistLocal[i][nTopos] = 2; |
220 |
+ |
nTopos++; |
221 |
+ |
} |
222 |
+ |
if (oneFour.hasPair(iglob, jglob)) { |
223 |
+ |
toposForLocalAtom[i].push_back(j); |
224 |
+ |
topoDistLocal[i][nTopos] = 3; |
225 |
+ |
nTopos++; |
226 |
+ |
} |
227 |
+ |
} |
228 |
+ |
} |
229 |
|
} |
230 |
+ |
|
231 |
+ |
void ForceMatrixDecomposition::zeroWorkArrays() { |
232 |
+ |
|
233 |
+ |
for (int j = 0; j < N_INTERACTION_FAMILIES; j++) { |
234 |
+ |
longRangePot_[j] = 0.0; |
235 |
+ |
} |
236 |
+ |
|
237 |
+ |
#ifdef IS_MPI |
238 |
+ |
if (storageLayout_ & DataStorage::dslForce) { |
239 |
+ |
fill(atomRowData.force.begin(), atomRowData.force.end(), V3Zero); |
240 |
+ |
fill(atomColData.force.begin(), atomColData.force.end(), V3Zero); |
241 |
+ |
} |
242 |
+ |
|
243 |
+ |
if (storageLayout_ & DataStorage::dslTorque) { |
244 |
+ |
fill(atomRowData.torque.begin(), atomRowData.torque.end(), V3Zero); |
245 |
+ |
fill(atomColData.torque.begin(), atomColData.torque.end(), V3Zero); |
246 |
+ |
} |
247 |
|
|
248 |
+ |
fill(pot_row.begin(), pot_row.end(), |
249 |
+ |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
250 |
|
|
251 |
+ |
fill(pot_col.begin(), pot_col.end(), |
252 |
+ |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
253 |
+ |
|
254 |
+ |
pot_local = Vector<RealType, N_INTERACTION_FAMILIES>(0.0); |
255 |
|
|
256 |
+ |
if (storageLayout_ & DataStorage::dslParticlePot) { |
257 |
+ |
fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), 0.0); |
258 |
+ |
fill(atomColData.particlePot.begin(), atomColData.particlePot.end(), 0.0); |
259 |
+ |
} |
260 |
+ |
|
261 |
+ |
if (storageLayout_ & DataStorage::dslDensity) { |
262 |
+ |
fill(atomRowData.density.begin(), atomRowData.density.end(), 0.0); |
263 |
+ |
fill(atomColData.density.begin(), atomColData.density.end(), 0.0); |
264 |
+ |
} |
265 |
+ |
|
266 |
+ |
if (storageLayout_ & DataStorage::dslFunctional) { |
267 |
+ |
fill(atomRowData.functional.begin(), atomRowData.functional.end(), 0.0); |
268 |
+ |
fill(atomColData.functional.begin(), atomColData.functional.end(), 0.0); |
269 |
+ |
} |
270 |
+ |
|
271 |
+ |
if (storageLayout_ & DataStorage::dslFunctionalDerivative) { |
272 |
+ |
fill(atomRowData.functionalDerivative.begin(), |
273 |
+ |
atomRowData.functionalDerivative.end(), 0.0); |
274 |
+ |
fill(atomColData.functionalDerivative.begin(), |
275 |
+ |
atomColData.functionalDerivative.end(), 0.0); |
276 |
+ |
} |
277 |
+ |
|
278 |
+ |
#else |
279 |
+ |
|
280 |
+ |
if (storageLayout_ & DataStorage::dslParticlePot) { |
281 |
+ |
fill(snap_->atomData.particlePot.begin(), |
282 |
+ |
snap_->atomData.particlePot.end(), 0.0); |
283 |
+ |
} |
284 |
+ |
|
285 |
+ |
if (storageLayout_ & DataStorage::dslDensity) { |
286 |
+ |
fill(snap_->atomData.density.begin(), |
287 |
+ |
snap_->atomData.density.end(), 0.0); |
288 |
+ |
} |
289 |
+ |
if (storageLayout_ & DataStorage::dslFunctional) { |
290 |
+ |
fill(snap_->atomData.functional.begin(), |
291 |
+ |
snap_->atomData.functional.end(), 0.0); |
292 |
+ |
} |
293 |
+ |
if (storageLayout_ & DataStorage::dslFunctionalDerivative) { |
294 |
+ |
fill(snap_->atomData.functionalDerivative.begin(), |
295 |
+ |
snap_->atomData.functionalDerivative.end(), 0.0); |
296 |
+ |
} |
297 |
+ |
#endif |
298 |
+ |
|
299 |
+ |
} |
300 |
+ |
|
301 |
+ |
|
302 |
|
void ForceMatrixDecomposition::distributeData() { |
303 |
|
snap_ = sman_->getCurrentSnapshot(); |
304 |
|
storageLayout_ = sman_->getStorageLayout(); |
334 |
|
#endif |
335 |
|
} |
336 |
|
|
337 |
+ |
/* collects information obtained during the pre-pair loop onto local |
338 |
+ |
* data structures. |
339 |
+ |
*/ |
340 |
|
void ForceMatrixDecomposition::collectIntermediateData() { |
341 |
|
snap_ = sman_->getCurrentSnapshot(); |
342 |
|
storageLayout_ = sman_->getStorageLayout(); |
348 |
|
snap_->atomData.density); |
349 |
|
|
350 |
|
int n = snap_->atomData.density.size(); |
351 |
< |
std::vector<RealType> rho_tmp(n, 0.0); |
351 |
> |
vector<RealType> rho_tmp(n, 0.0); |
352 |
|
AtomCommRealColumn->scatter(atomColData.density, rho_tmp); |
353 |
|
for (int i = 0; i < n; i++) |
354 |
|
snap_->atomData.density[i] += rho_tmp[i]; |
355 |
|
} |
356 |
|
#endif |
357 |
|
} |
358 |
< |
|
358 |
> |
|
359 |
> |
/* |
360 |
> |
* redistributes information obtained during the pre-pair loop out to |
361 |
> |
* row and column-indexed data structures |
362 |
> |
*/ |
363 |
|
void ForceMatrixDecomposition::distributeIntermediateData() { |
364 |
|
snap_ = sman_->getCurrentSnapshot(); |
365 |
|
storageLayout_ = sman_->getStorageLayout(); |
417 |
|
|
418 |
|
nLocal_ = snap_->getNumberOfAtoms(); |
419 |
|
|
420 |
< |
vector<vector<RealType> > pot_temp(N_INTERACTION_FAMILIES, |
421 |
< |
vector<RealType> (nLocal_, 0.0)); |
420 |
> |
vector<potVec> pot_temp(nLocal_, |
421 |
> |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
422 |
> |
|
423 |
> |
// scatter/gather pot_row into the members of my column |
424 |
> |
|
425 |
> |
AtomCommPotRow->scatter(pot_row, pot_temp); |
426 |
> |
|
427 |
> |
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
428 |
> |
pot_local += pot_temp[ii]; |
429 |
|
|
430 |
< |
for (int i = 0; i < N_INTERACTION_FAMILIES; i++) { |
431 |
< |
AtomCommRealRow->scatter(pot_row[i], pot_temp[i]); |
432 |
< |
for (int ii = 0; ii < pot_temp[i].size(); ii++ ) { |
433 |
< |
pot_local[i] += pot_temp[i][ii]; |
434 |
< |
} |
435 |
< |
} |
430 |
> |
fill(pot_temp.begin(), pot_temp.end(), |
431 |
> |
Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
432 |
> |
|
433 |
> |
AtomCommPotColumn->scatter(pot_col, pot_temp); |
434 |
> |
|
435 |
> |
for (int ii = 0; ii < pot_temp.size(); ii++ ) |
436 |
> |
pot_local += pot_temp[ii]; |
437 |
> |
|
438 |
|
#endif |
439 |
|
} |
440 |
|
|
441 |
+ |
int ForceMatrixDecomposition::getNAtomsInRow() { |
442 |
+ |
#ifdef IS_MPI |
443 |
+ |
return nAtomsInRow_; |
444 |
+ |
#else |
445 |
+ |
return nLocal_; |
446 |
+ |
#endif |
447 |
+ |
} |
448 |
+ |
|
449 |
+ |
/** |
450 |
+ |
* returns the list of atoms belonging to this group. |
451 |
+ |
*/ |
452 |
+ |
vector<int> ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){ |
453 |
+ |
#ifdef IS_MPI |
454 |
+ |
return groupListRow_[cg1]; |
455 |
+ |
#else |
456 |
+ |
return groupList_[cg1]; |
457 |
+ |
#endif |
458 |
+ |
} |
459 |
+ |
|
460 |
+ |
vector<int> ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){ |
461 |
+ |
#ifdef IS_MPI |
462 |
+ |
return groupListCol_[cg2]; |
463 |
+ |
#else |
464 |
+ |
return groupList_[cg2]; |
465 |
+ |
#endif |
466 |
+ |
} |
467 |
|
|
468 |
|
Vector3d ForceMatrixDecomposition::getIntergroupVector(int cg1, int cg2){ |
469 |
|
Vector3d d; |
505 |
|
snap_->wrapVector(d); |
506 |
|
return d; |
507 |
|
} |
508 |
+ |
|
509 |
+ |
RealType ForceMatrixDecomposition::getMassFactorRow(int atom1) { |
510 |
+ |
#ifdef IS_MPI |
511 |
+ |
return massFactorsRow[atom1]; |
512 |
+ |
#else |
513 |
+ |
return massFactorsLocal[atom1]; |
514 |
+ |
#endif |
515 |
+ |
} |
516 |
+ |
|
517 |
+ |
RealType ForceMatrixDecomposition::getMassFactorColumn(int atom2) { |
518 |
+ |
#ifdef IS_MPI |
519 |
+ |
return massFactorsCol[atom2]; |
520 |
+ |
#else |
521 |
+ |
return massFactorsLocal[atom2]; |
522 |
+ |
#endif |
523 |
+ |
|
524 |
+ |
} |
525 |
|
|
526 |
|
Vector3d ForceMatrixDecomposition::getInteratomicVector(int atom1, int atom2){ |
527 |
|
Vector3d d; |
536 |
|
return d; |
537 |
|
} |
538 |
|
|
539 |
+ |
vector<int> ForceMatrixDecomposition::getSkipsForRowAtom(int atom1) { |
540 |
+ |
#ifdef IS_MPI |
541 |
+ |
return skipsForRowAtom[atom1]; |
542 |
+ |
#else |
543 |
+ |
return skipsForLocalAtom[atom1]; |
544 |
+ |
#endif |
545 |
+ |
} |
546 |
+ |
|
547 |
+ |
/** |
548 |
+ |
* There are a number of reasons to skip a pair or a |
549 |
+ |
* particle. Mostly we do this to exclude atoms who are involved in |
550 |
+ |
* short range interactions (bonds, bends, torsions), but we also |
551 |
+ |
* need to exclude some overcounted interactions that result from |
552 |
+ |
* the parallel decomposition. |
553 |
+ |
*/ |
554 |
+ |
bool ForceMatrixDecomposition::skipAtomPair(int atom1, int atom2) { |
555 |
+ |
int unique_id_1, unique_id_2; |
556 |
+ |
|
557 |
+ |
#ifdef IS_MPI |
558 |
+ |
// in MPI, we have to look up the unique IDs for each atom |
559 |
+ |
unique_id_1 = AtomRowToGlobal[atom1]; |
560 |
+ |
unique_id_2 = AtomColToGlobal[atom2]; |
561 |
+ |
|
562 |
+ |
// this situation should only arise in MPI simulations |
563 |
+ |
if (unique_id_1 == unique_id_2) return true; |
564 |
+ |
|
565 |
+ |
// this prevents us from doing the pair on multiple processors |
566 |
+ |
if (unique_id_1 < unique_id_2) { |
567 |
+ |
if ((unique_id_1 + unique_id_2) % 2 == 0) return true; |
568 |
+ |
} else { |
569 |
+ |
if ((unique_id_1 + unique_id_2) % 2 == 1) return true; |
570 |
+ |
} |
571 |
+ |
#else |
572 |
+ |
// in the normal loop, the atom numbers are unique |
573 |
+ |
unique_id_1 = atom1; |
574 |
+ |
unique_id_2 = atom2; |
575 |
+ |
#endif |
576 |
+ |
|
577 |
+ |
#ifdef IS_MPI |
578 |
+ |
for (vector<int>::iterator i = skipsForRowAtom[atom1].begin(); |
579 |
+ |
i != skipsForRowAtom[atom1].end(); ++i) { |
580 |
+ |
if ( (*i) == unique_id_2 ) return true; |
581 |
+ |
} |
582 |
+ |
#else |
583 |
+ |
for (vector<int>::iterator i = skipsForLocalAtom[atom1].begin(); |
584 |
+ |
i != skipsForLocalAtom[atom1].end(); ++i) { |
585 |
+ |
if ( (*i) == unique_id_2 ) return true; |
586 |
+ |
} |
587 |
+ |
#endif |
588 |
+ |
} |
589 |
+ |
|
590 |
+ |
int ForceMatrixDecomposition::getTopoDistance(int atom1, int atom2) { |
591 |
+ |
|
592 |
+ |
#ifdef IS_MPI |
593 |
+ |
for (int i = 0; i < toposForRowAtom[atom1].size(); i++) { |
594 |
+ |
if ( toposForRowAtom[atom1][i] == atom2 ) return topoDistRow[atom1][i]; |
595 |
+ |
} |
596 |
+ |
#else |
597 |
+ |
for (int i = 0; i < toposForLocalAtom[atom1].size(); i++) { |
598 |
+ |
if ( toposForLocalAtom[atom1][i] == atom2 ) return topoDistLocal[atom1][i]; |
599 |
+ |
} |
600 |
+ |
#endif |
601 |
+ |
|
602 |
+ |
// zero is default for unconnected (i.e. normal) pair interactions |
603 |
+ |
return 0; |
604 |
+ |
} |
605 |
+ |
|
606 |
|
void ForceMatrixDecomposition::addForceToAtomRow(int atom1, Vector3d fg){ |
607 |
|
#ifdef IS_MPI |
608 |
|
atomRowData.force[atom1] += fg; |
624 |
|
InteractionData idat; |
625 |
|
|
626 |
|
#ifdef IS_MPI |
627 |
+ |
|
628 |
+ |
idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), |
629 |
+ |
ff_->getAtomType(identsCol[atom2]) ); |
630 |
+ |
|
631 |
+ |
|
632 |
|
if (storageLayout_ & DataStorage::dslAmat) { |
633 |
|
idat.A1 = &(atomRowData.aMat[atom1]); |
634 |
|
idat.A2 = &(atomColData.aMat[atom2]); |
649 |
|
idat.rho2 = &(atomColData.density[atom2]); |
650 |
|
} |
651 |
|
|
652 |
+ |
if (storageLayout_ & DataStorage::dslFunctional) { |
653 |
+ |
idat.frho1 = &(atomRowData.functional[atom1]); |
654 |
+ |
idat.frho2 = &(atomColData.functional[atom2]); |
655 |
+ |
} |
656 |
+ |
|
657 |
|
if (storageLayout_ & DataStorage::dslFunctionalDerivative) { |
658 |
|
idat.dfrho1 = &(atomRowData.functionalDerivative[atom1]); |
659 |
|
idat.dfrho2 = &(atomColData.functionalDerivative[atom2]); |
660 |
|
} |
661 |
+ |
|
662 |
+ |
if (storageLayout_ & DataStorage::dslParticlePot) { |
663 |
+ |
idat.particlePot1 = &(atomRowData.particlePot[atom1]); |
664 |
+ |
idat.particlePot2 = &(atomColData.particlePot[atom2]); |
665 |
+ |
} |
666 |
+ |
|
667 |
|
#else |
668 |
+ |
|
669 |
+ |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
670 |
+ |
ff_->getAtomType(identsLocal[atom2]) ); |
671 |
+ |
|
672 |
|
if (storageLayout_ & DataStorage::dslAmat) { |
673 |
|
idat.A1 = &(snap_->atomData.aMat[atom1]); |
674 |
|
idat.A2 = &(snap_->atomData.aMat[atom2]); |
689 |
|
idat.rho2 = &(snap_->atomData.density[atom2]); |
690 |
|
} |
691 |
|
|
692 |
+ |
if (storageLayout_ & DataStorage::dslFunctional) { |
693 |
+ |
idat.frho1 = &(snap_->atomData.functional[atom1]); |
694 |
+ |
idat.frho2 = &(snap_->atomData.functional[atom2]); |
695 |
+ |
} |
696 |
+ |
|
697 |
|
if (storageLayout_ & DataStorage::dslFunctionalDerivative) { |
698 |
|
idat.dfrho1 = &(snap_->atomData.functionalDerivative[atom1]); |
699 |
|
idat.dfrho2 = &(snap_->atomData.functionalDerivative[atom2]); |
700 |
|
} |
701 |
+ |
|
702 |
+ |
if (storageLayout_ & DataStorage::dslParticlePot) { |
703 |
+ |
idat.particlePot1 = &(snap_->atomData.particlePot[atom1]); |
704 |
+ |
idat.particlePot2 = &(snap_->atomData.particlePot[atom2]); |
705 |
+ |
} |
706 |
+ |
|
707 |
|
#endif |
708 |
|
return idat; |
709 |
+ |
} |
710 |
+ |
|
711 |
+ |
|
712 |
+ |
void ForceMatrixDecomposition::unpackInteractionData(InteractionData idat, int atom1, int atom2) { |
713 |
+ |
#ifdef IS_MPI |
714 |
+ |
pot_row[atom1] += 0.5 * *(idat.pot); |
715 |
+ |
pot_col[atom2] += 0.5 * *(idat.pot); |
716 |
+ |
|
717 |
+ |
atomRowData.force[atom1] += *(idat.f1); |
718 |
+ |
atomColData.force[atom2] -= *(idat.f1); |
719 |
+ |
#else |
720 |
+ |
longRangePot_ += *(idat.pot); |
721 |
+ |
|
722 |
+ |
snap_->atomData.force[atom1] += *(idat.f1); |
723 |
+ |
snap_->atomData.force[atom2] -= *(idat.f1); |
724 |
+ |
#endif |
725 |
+ |
|
726 |
|
} |
727 |
|
|
728 |
+ |
|
729 |
|
InteractionData ForceMatrixDecomposition::fillSkipData(int atom1, int atom2){ |
730 |
|
|
731 |
|
InteractionData idat; |
732 |
|
#ifdef IS_MPI |
733 |
+ |
idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), |
734 |
+ |
ff_->getAtomType(identsCol[atom2]) ); |
735 |
+ |
|
736 |
|
if (storageLayout_ & DataStorage::dslElectroFrame) { |
737 |
|
idat.eFrame1 = &(atomRowData.electroFrame[atom1]); |
738 |
|
idat.eFrame2 = &(atomColData.electroFrame[atom2]); |
741 |
|
idat.t1 = &(atomRowData.torque[atom1]); |
742 |
|
idat.t2 = &(atomColData.torque[atom2]); |
743 |
|
} |
388 |
– |
if (storageLayout_ & DataStorage::dslForce) { |
389 |
– |
idat.t1 = &(atomRowData.force[atom1]); |
390 |
– |
idat.t2 = &(atomColData.force[atom2]); |
391 |
– |
} |
744 |
|
#else |
745 |
+ |
idat.atypes = make_pair( ff_->getAtomType(identsLocal[atom1]), |
746 |
+ |
ff_->getAtomType(identsLocal[atom2]) ); |
747 |
+ |
|
748 |
|
if (storageLayout_ & DataStorage::dslElectroFrame) { |
749 |
|
idat.eFrame1 = &(snap_->atomData.electroFrame[atom1]); |
750 |
|
idat.eFrame2 = &(snap_->atomData.electroFrame[atom2]); |
753 |
|
idat.t1 = &(snap_->atomData.torque[atom1]); |
754 |
|
idat.t2 = &(snap_->atomData.torque[atom2]); |
755 |
|
} |
756 |
< |
if (storageLayout_ & DataStorage::dslForce) { |
402 |
< |
idat.t1 = &(snap_->atomData.force[atom1]); |
403 |
< |
idat.t2 = &(snap_->atomData.force[atom2]); |
404 |
< |
} |
405 |
< |
#endif |
406 |
< |
|
756 |
> |
#endif |
757 |
|
} |
758 |
|
|
409 |
– |
|
410 |
– |
|
411 |
– |
|
759 |
|
/* |
760 |
|
* buildNeighborList |
761 |
|
* |
854 |
|
} |
855 |
|
#endif |
856 |
|
|
510 |
– |
|
511 |
– |
|
857 |
|
for (int m1z = 0; m1z < nCells_.z(); m1z++) { |
858 |
|
for (int m1y = 0; m1y < nCells_.y(); m1y++) { |
859 |
|
for (int m1x = 0; m1x < nCells_.x(); m1x++) { |