ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/parallel/ForceMatrixDecomposition.cpp
(Generate patch)

Comparing branches/development/src/parallel/ForceMatrixDecomposition.cpp (file contents):
Revision 1569 by gezelter, Thu May 26 13:55:04 2011 UTC vs.
Revision 1570 by gezelter, Thu May 26 21:56:04 2011 UTC

# Line 42 | Line 42
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 {
# Line 63 | Line 64 | namespace OpenMD {
64      cgLocalToGlobal = info_->getGlobalGroupIndices();
65      vector<int> globalGroupMembership = info_->getGlobalGroupMembership();
66      vector<RealType> massFactorsLocal = info_->getMassFactors();
67 +    PairList excludes = info_->getExcludedInteractions();
68 +    PairList oneTwo = info_->getOneTwoInteractions();
69 +    PairList oneThree = info_->getOneThreeInteractions();
70 +    PairList oneFour = info_->getOneFourInteractions();
71      vector<RealType> pot_local(N_INTERACTION_FAMILIES, 0.0);
72  
73   #ifdef IS_MPI
# Line 139 | Line 144 | namespace OpenMD {
144        }      
145      }
146  
147 +    skipsForRowAtom.clear();
148 +    skipsForRowAtom.reserve(nAtomsInRow_);
149 +    for (int i = 0; i < nAtomsInRow_; i++) {
150 +      int iglob = AtomColToGlobal[i];
151 +      for (int j = 0; j < nAtomsInCol_; j++) {
152 +        int jglob = AtomRowToGlobal[j];        
153 +        if (excludes.hasPair(iglob, jglob))
154 +          skipsForRowAtom[i].push_back(j);      
155 +      }      
156 +    }
157 +
158 +    toposForRowAtom.clear();
159 +    toposForRowAtom.reserve(nAtomsInRow_);
160 +    for (int i = 0; i < nAtomsInRow_; i++) {
161 +      int iglob = AtomColToGlobal[i];
162 +      int nTopos = 0;
163 +      for (int j = 0; j < nAtomsInCol_; j++) {
164 +        int jglob = AtomRowToGlobal[j];        
165 +        if (oneTwo.hasPair(iglob, jglob)) {
166 +          toposForRowAtom[i].push_back(j);
167 +          topoDistRow[i][nTopos] = 1;
168 +          nTopos++;
169 +        }
170 +        if (oneThree.hasPair(iglob, jglob)) {
171 +          toposForRowAtom[i].push_back(j);
172 +          topoDistRow[i][nTopos] = 2;
173 +          nTopos++;
174 +        }
175 +        if (oneFour.hasPair(iglob, jglob)) {
176 +          toposForRowAtom[i].push_back(j);
177 +          topoDistRow[i][nTopos] = 3;
178 +          nTopos++;
179 +        }
180 +      }      
181 +    }
182 +
183   #endif
184  
185      groupList_.clear();
# Line 152 | Line 193 | namespace OpenMD {
193        }      
194      }
195  
196 <  
197 <    // still need:
157 <    // topoDist
158 <    // exclude
196 >    skipsForLocalAtom.clear();
197 >    skipsForLocalAtom.reserve(nLocal_);
198  
199 <  }
200 <    
199 >    for (int i = 0; i < nLocal_; i++) {
200 >      int iglob = AtomLocalToGlobal[i];
201 >      for (int j = 0; j < nLocal_; j++) {
202 >        int jglob = AtomLocalToGlobal[j];        
203 >        if (excludes.hasPair(iglob, jglob))
204 >          skipsForLocalAtom[i].push_back(j);      
205 >      }      
206 >    }
207  
208 <
208 >    toposForLocalAtom.clear();
209 >    toposForLocalAtom.reserve(nLocal_);
210 >    for (int i = 0; i < nLocal_; i++) {
211 >      int iglob = AtomLocalToGlobal[i];
212 >      int nTopos = 0;
213 >      for (int j = 0; j < nLocal_; j++) {
214 >        int jglob = AtomLocalToGlobal[j];        
215 >        if (oneTwo.hasPair(iglob, jglob)) {
216 >          toposForLocalAtom[i].push_back(j);
217 >          topoDistLocal[i][nTopos] = 1;
218 >          nTopos++;
219 >        }
220 >        if (oneThree.hasPair(iglob, jglob)) {
221 >          toposForLocalAtom[i].push_back(j);
222 >          topoDistLocal[i][nTopos] = 2;
223 >          nTopos++;
224 >        }
225 >        if (oneFour.hasPair(iglob, jglob)) {
226 >          toposForLocalAtom[i].push_back(j);
227 >          topoDistLocal[i][nTopos] = 3;
228 >          nTopos++;
229 >        }
230 >      }      
231 >    }
232 >  }
233 >  
234    void ForceMatrixDecomposition::distributeData()  {
235      snap_ = sman_->getCurrentSnapshot();
236      storageLayout_ = sman_->getStorageLayout();
# Line 284 | Line 354 | namespace OpenMD {
354   #endif
355    }
356  
357 +  int ForceMatrixDecomposition::getNAtomsInRow() {  
358 + #ifdef IS_MPI
359 +    return nAtomsInRow_;
360 + #else
361 +    return nLocal_;
362 + #endif
363 +  }
364 +
365    /**
366     * returns the list of atoms belonging to this group.  
367     */
# Line 374 | Line 452 | namespace OpenMD {
452      return d;    
453    }
454  
455 +  vector<int> ForceMatrixDecomposition::getSkipsForRowAtom(int atom1) {
456 + #ifdef IS_MPI
457 +    return skipsForRowAtom[atom1];
458 + #else
459 +    return skipsForLocalAtom[atom1];
460 + #endif
461 +  }
462 +
463 +  /**
464 +   * there are a number of reasons to skip a pair or a particle mostly
465 +   * we do this to exclude atoms who are involved in short range
466 +   * interactions (bonds, bends, torsions), but we also need to
467 +   * exclude some overcounted interactions that result from the
468 +   * parallel decomposition.
469 +   */
470 +  bool ForceMatrixDecomposition::skipAtomPair(int atom1, int atom2) {
471 +    int unique_id_1, unique_id_2;
472 +
473 + #ifdef IS_MPI
474 +    // in MPI, we have to look up the unique IDs for each atom
475 +    unique_id_1 = AtomRowToGlobal[atom1];
476 +    unique_id_2 = AtomColToGlobal[atom2];
477 +
478 +    // this situation should only arise in MPI simulations
479 +    if (unique_id_1 == unique_id_2) return true;
480 +    
481 +    // this prevents us from doing the pair on multiple processors
482 +    if (unique_id_1 < unique_id_2) {
483 +      if ((unique_id_1 + unique_id_2) % 2 == 0) return true;
484 +    } else {
485 +      if ((unique_id_1 + unique_id_2) % 2 == 1) return true;
486 +    }
487 + #else
488 +    // in the normal loop, the atom numbers are unique
489 +    unique_id_1 = atom1;
490 +    unique_id_2 = atom2;
491 + #endif
492 +    
493 + #ifdef IS_MPI
494 +    for (vector<int>::iterator i = skipsForRowAtom[atom1].begin();
495 +         i != skipsForRowAtom[atom1].end(); ++i) {
496 +      if ( (*i) == unique_id_2 ) return true;
497 +    }    
498 + #else
499 +    for (vector<int>::iterator i = skipsForLocalAtom[atom1].begin();
500 +         i != skipsForLocalAtom[atom1].end(); ++i) {
501 +      if ( (*i) == unique_id_2 ) return true;
502 +    }    
503 + #endif
504 +  }
505 +
506 +  int ForceMatrixDecomposition::getTopoDistance(int atom1, int atom2) {
507 +    
508 + #ifdef IS_MPI
509 +    for (int i = 0; i < toposForRowAtom[atom1].size(); i++) {
510 +      if ( toposForRowAtom[atom1][i] == atom2 ) return topoDistRow[atom1][i];
511 +    }
512 + #else
513 +    for (int i = 0; i < toposForLocalAtom[atom1].size(); i++) {
514 +      if ( toposForLocalAtom[atom1][i] == atom2 ) return topoDistLocal[atom1][i];
515 +    }
516 + #endif
517 +
518 +    // zero is default for unconnected (i.e. normal) pair interactions
519 +    return 0;
520 +  }
521 +
522    void ForceMatrixDecomposition::addForceToAtomRow(int atom1, Vector3d fg){
523   #ifdef IS_MPI
524      atomRowData.force[atom1] += fg;
# Line 419 | Line 564 | namespace OpenMD {
564        idat.dfrho1 = &(atomRowData.functionalDerivative[atom1]);
565        idat.dfrho2 = &(atomColData.functionalDerivative[atom2]);
566      }
567 +
568   #else
569      if (storageLayout_ & DataStorage::dslAmat) {
570        idat.A1 = &(snap_->atomData.aMat[atom1]);

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines