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 1567 by gezelter, Tue May 24 21:24:45 2011 UTC vs.
Revision 1569 by gezelter, Thu May 26 13:55:04 2011 UTC

# Line 57 | Line 57 | namespace OpenMD {
57      nLocal_ = snap_->getNumberOfAtoms();
58      nGroups_ = snap_->getNumberOfCutoffGroups();
59  
60 +    // gather the information for atomtype IDs (atids):
61 +    vector<int> identsLocal = info_->getIdentArray();
62 +    AtomLocalToGlobal = info_->getGlobalAtomIndices();
63 +    cgLocalToGlobal = info_->getGlobalGroupIndices();
64 +    vector<int> globalGroupMembership = info_->getGlobalGroupMembership();
65 +    vector<RealType> massFactorsLocal = info_->getMassFactors();
66 +    vector<RealType> pot_local(N_INTERACTION_FAMILIES, 0.0);
67 +
68   #ifdef IS_MPI
69  
70      AtomCommIntRow = new Communicator<Row,int>(nLocal_);
# Line 93 | Line 101 | namespace OpenMD {
101                                        vector<RealType> (nAtomsInRow_, 0.0));
102      vector<vector<RealType> > pot_col(N_INTERACTION_FAMILIES,
103                                        vector<RealType> (nAtomsInCol_, 0.0));
96 –
97 –
98 –    vector<RealType> pot_local(N_INTERACTION_FAMILIES, 0.0);
104      
100 –    // gather the information for atomtype IDs (atids):
101 –    vector<int> identsLocal = info_->getIdentArray();
105      identsRow.reserve(nAtomsInRow_);
106      identsCol.reserve(nAtomsInCol_);
107      
108      AtomCommIntRow->gather(identsLocal, identsRow);
109      AtomCommIntColumn->gather(identsLocal, identsCol);
110      
108 –    AtomLocalToGlobal = info_->getGlobalAtomIndices();
111      AtomCommIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal);
112      AtomCommIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal);
113      
112 –    cgLocalToGlobal = info_->getGlobalGroupIndices();
114      cgCommIntRow->gather(cgLocalToGlobal, cgRowToGlobal);
115      cgCommIntColumn->gather(cgLocalToGlobal, cgColToGlobal);
116 +
117 +    AtomCommRealRow->gather(massFactorsLocal, massFactorsRow);
118 +    AtomCommRealColumn->gather(massFactorsLocal, massFactorsCol);
119 +
120 +    groupListRow_.clear();
121 +    groupListRow_.reserve(nGroupsInRow_);
122 +    for (int i = 0; i < nGroupsInRow_; i++) {
123 +      int gid = cgRowToGlobal[i];
124 +      for (int j = 0; j < nAtomsInRow_; j++) {
125 +        int aid = AtomRowToGlobal[j];
126 +        if (globalGroupMembership[aid] == gid)
127 +          groupListRow_[i].push_back(j);
128 +      }      
129 +    }
130 +
131 +    groupListCol_.clear();
132 +    groupListCol_.reserve(nGroupsInCol_);
133 +    for (int i = 0; i < nGroupsInCol_; i++) {
134 +      int gid = cgColToGlobal[i];
135 +      for (int j = 0; j < nAtomsInCol_; j++) {
136 +        int aid = AtomColToGlobal[j];
137 +        if (globalGroupMembership[aid] == gid)
138 +          groupListCol_[i].push_back(j);
139 +      }      
140 +    }
141  
142 + #endif
143 +
144 +    groupList_.clear();
145 +    groupList_.reserve(nGroups_);
146 +    for (int i = 0; i < nGroups_; i++) {
147 +      int gid = cgLocalToGlobal[i];
148 +      for (int j = 0; j < nLocal_; j++) {
149 +        int aid = AtomLocalToGlobal[j];
150 +        if (globalGroupMembership[aid] == gid)
151 +          groupList_[i].push_back(j);
152 +      }      
153 +    }
154 +
155 +  
156      // still need:
157      // topoDist
158      // exclude
159 < #endif
159 >
160    }
161      
162  
# Line 241 | Line 281 | namespace OpenMD {
281          pot_local[i] += pot_temp[i][ii];
282        }
283      }
284 + #endif
285 +  }
286 +
287 +  /**
288 +   * returns the list of atoms belonging to this group.  
289 +   */
290 +  vector<int> ForceMatrixDecomposition::getAtomsInGroupRow(int cg1){
291 + #ifdef IS_MPI
292 +    return groupListRow_[cg1];
293 + #else
294 +    return groupList_[cg1];
295   #endif
296    }
297  
298 +  vector<int> ForceMatrixDecomposition::getAtomsInGroupColumn(int cg2){
299 + #ifdef IS_MPI
300 +    return groupListCol_[cg2];
301 + #else
302 +    return groupList_[cg2];
303 + #endif
304 +  }
305    
306    Vector3d ForceMatrixDecomposition::getIntergroupVector(int cg1, int cg2){
307      Vector3d d;
# Line 285 | Line 343 | namespace OpenMD {
343      snap_->wrapVector(d);
344      return d;    
345    }
346 +
347 +  RealType ForceMatrixDecomposition::getMassFactorRow(int atom1) {
348 + #ifdef IS_MPI
349 +    return massFactorsRow[atom1];
350 + #else
351 +    return massFactorsLocal[atom1];
352 + #endif
353 +  }
354 +
355 +  RealType ForceMatrixDecomposition::getMassFactorColumn(int atom2) {
356 + #ifdef IS_MPI
357 +    return massFactorsCol[atom2];
358 + #else
359 +    return massFactorsLocal[atom2];
360 + #endif
361 +
362 +  }
363      
364    Vector3d ForceMatrixDecomposition::getInteratomicVector(int atom1, int atom2){
365      Vector3d d;
# Line 406 | Line 481 | namespace OpenMD {
481      
482    }
483  
409 –  SelfData ForceMatrixDecomposition::fillSelfData(int atom1) {
410 –    SelfData sdat;
411 –    // Still Missing atype, skippedCharge, potVec pot,
412 –    if (storageLayout_ & DataStorage::dslElectroFrame) {
413 –      sdat.eFrame = &(snap_->atomData.electroFrame[atom1]);
414 –    }
415 –    
416 –    if (storageLayout_ & DataStorage::dslTorque) {
417 –      sdat.t = &(snap_->atomData.torque[atom1]);
418 –    }
419 –    
420 –    if (storageLayout_ & DataStorage::dslDensity) {
421 –      sdat.rho = &(snap_->atomData.density[atom1]);
422 –    }
423 –    
424 –    if (storageLayout_ & DataStorage::dslFunctional) {
425 –      sdat.frho = &(snap_->atomData.functional[atom1]);
426 –    }
427 –    
428 –    if (storageLayout_ & DataStorage::dslFunctionalDerivative) {
429 –      sdat.dfrhodrho = &(snap_->atomData.functionalDerivative[atom1]);
430 –    }
484  
432 –    return sdat;    
433 –  }
485  
486  
436 –
487    /*
488     * buildNeighborList
489     *
# Line 444 | Line 494 | namespace OpenMD {
494        
495      vector<pair<int, int> > neighborList;
496   #ifdef IS_MPI
497 <    CellListRow.clear();
498 <    CellListCol.clear();
497 >    cellListRow_.clear();
498 >    cellListCol_.clear();
499   #else
500 <    CellList.clear();
500 >    cellList_.clear();
501   #endif
502  
503      // dangerous to not do error checking.
454 –    RealType skinThickness_ = info_->getSimParams()->getSkinThickness();
504      RealType rCut_;
505  
506      RealType rList_ = (rCut_ + skinThickness_);
# Line 461 | Line 510 | namespace OpenMD {
510      Vector3d Hx = Hmat.getColumn(0);
511      Vector3d Hy = Hmat.getColumn(1);
512      Vector3d Hz = Hmat.getColumn(2);
464 –    Vector3i nCells;
513  
514 <    nCells.x() = (int) ( Hx.length() )/ rList_;
515 <    nCells.y() = (int) ( Hy.length() )/ rList_;
516 <    nCells.z() = (int) ( Hz.length() )/ rList_;
514 >    nCells_.x() = (int) ( Hx.length() )/ rList_;
515 >    nCells_.y() = (int) ( Hy.length() )/ rList_;
516 >    nCells_.z() = (int) ( Hz.length() )/ rList_;
517  
518      Mat3x3d invHmat = snap_->getInvHmat();
519      Vector3d rs, scaled, dr;
# Line 483 | Line 531 | namespace OpenMD {
531          scaled[j] -= roundMe(scaled[j]);
532      
533        // find xyz-indices of cell that cutoffGroup is in.
534 <      whichCell.x() = nCells.x() * scaled.x();
535 <      whichCell.y() = nCells.y() * scaled.y();
536 <      whichCell.z() = nCells.z() * scaled.z();
534 >      whichCell.x() = nCells_.x() * scaled.x();
535 >      whichCell.y() = nCells_.y() * scaled.y();
536 >      whichCell.z() = nCells_.z() * scaled.z();
537  
538        // find single index of this cell:
539 <      cellIndex = Vlinear(whichCell, nCells);
539 >      cellIndex = Vlinear(whichCell, nCells_);
540        // add this cutoff group to the list of groups in this cell;
541 <      CellListRow[cellIndex].push_back(i);
541 >      cellListRow_[cellIndex].push_back(i);
542      }
543  
544      for (int i = 0; i < nGroupsInCol_; i++) {
# Line 503 | Line 551 | namespace OpenMD {
551          scaled[j] -= roundMe(scaled[j]);
552  
553        // find xyz-indices of cell that cutoffGroup is in.
554 <      whichCell.x() = nCells.x() * scaled.x();
555 <      whichCell.y() = nCells.y() * scaled.y();
556 <      whichCell.z() = nCells.z() * scaled.z();
554 >      whichCell.x() = nCells_.x() * scaled.x();
555 >      whichCell.y() = nCells_.y() * scaled.y();
556 >      whichCell.z() = nCells_.z() * scaled.z();
557  
558        // find single index of this cell:
559 <      cellIndex = Vlinear(whichCell, nCells);
559 >      cellIndex = Vlinear(whichCell, nCells_);
560        // add this cutoff group to the list of groups in this cell;
561 <      CellListCol[cellIndex].push_back(i);
561 >      cellListCol_[cellIndex].push_back(i);
562      }
563   #else
564      for (int i = 0; i < nGroups_; i++) {
# Line 523 | Line 571 | namespace OpenMD {
571          scaled[j] -= roundMe(scaled[j]);
572  
573        // find xyz-indices of cell that cutoffGroup is in.
574 <      whichCell.x() = nCells.x() * scaled.x();
575 <      whichCell.y() = nCells.y() * scaled.y();
576 <      whichCell.z() = nCells.z() * scaled.z();
574 >      whichCell.x() = nCells_.x() * scaled.x();
575 >      whichCell.y() = nCells_.y() * scaled.y();
576 >      whichCell.z() = nCells_.z() * scaled.z();
577  
578        // find single index of this cell:
579 <      cellIndex = Vlinear(whichCell, nCells);
579 >      cellIndex = Vlinear(whichCell, nCells_);
580        // add this cutoff group to the list of groups in this cell;
581 <      CellList[cellIndex].push_back(i);
581 >      cellList_[cellIndex].push_back(i);
582      }
583   #endif
584  
585  
586  
587 <    for (int m1z = 0; m1z < nCells.z(); m1z++) {
588 <      for (int m1y = 0; m1y < nCells.y(); m1y++) {
589 <        for (int m1x = 0; m1x < nCells.x(); m1x++) {
587 >    for (int m1z = 0; m1z < nCells_.z(); m1z++) {
588 >      for (int m1y = 0; m1y < nCells_.y(); m1y++) {
589 >        for (int m1x = 0; m1x < nCells_.x(); m1x++) {
590            Vector3i m1v(m1x, m1y, m1z);
591 <          int m1 = Vlinear(m1v, nCells);
544 <          for (int offset = 0; offset < nOffset_; offset++) {
545 <            Vector3i m2v = m1v + cellOffsets_[offset];
591 >          int m1 = Vlinear(m1v, nCells_);
592  
593 <            if (m2v.x() >= nCells.x()) {
593 >          for (vector<Vector3i>::iterator os = cellOffsets_.begin();
594 >               os != cellOffsets_.end(); ++os) {
595 >            
596 >            Vector3i m2v = m1v + (*os);
597 >            
598 >            if (m2v.x() >= nCells_.x()) {
599                m2v.x() = 0;          
600              } else if (m2v.x() < 0) {
601 <              m2v.x() = nCells.x() - 1;
601 >              m2v.x() = nCells_.x() - 1;
602              }
603 <
604 <            if (m2v.y() >= nCells.y()) {
603 >            
604 >            if (m2v.y() >= nCells_.y()) {
605                m2v.y() = 0;          
606              } else if (m2v.y() < 0) {
607 <              m2v.y() = nCells.y() - 1;
607 >              m2v.y() = nCells_.y() - 1;
608              }
609 <
610 <            if (m2v.z() >= nCells.z()) {
609 >            
610 >            if (m2v.z() >= nCells_.z()) {
611                m2v.z() = 0;          
612              } else if (m2v.z() < 0) {
613 <              m2v.z() = nCells.z() - 1;
613 >              m2v.z() = nCells_.z() - 1;
614              }
615 +            
616 +            int m2 = Vlinear (m2v, nCells_);
617  
565 –            int m2 = Vlinear (m2v, nCells);
566 –
618   #ifdef IS_MPI
619 <            for (vector<int>::iterator j1 = CellListRow[m1].begin();
620 <                 j1 != CellListRow[m1].end(); ++j1) {
621 <              for (vector<int>::iterator j2 = CellListCol[m2].begin();
622 <                   j2 != CellListCol[m2].end(); ++j2) {
619 >            for (vector<int>::iterator j1 = cellListRow_[m1].begin();
620 >                 j1 != cellListRow_[m1].end(); ++j1) {
621 >              for (vector<int>::iterator j2 = cellListCol_[m2].begin();
622 >                   j2 != cellListCol_[m2].end(); ++j2) {
623                                
624                  // Always do this if we're in different cells or if
625                  // we're in the same cell and the global index of the
# Line 584 | Line 635 | namespace OpenMD {
635                }
636              }
637   #else
638 <            for (vector<int>::iterator j1 = CellList[m1].begin();
639 <                 j1 != CellList[m1].end(); ++j1) {
640 <              for (vector<int>::iterator j2 = CellList[m2].begin();
641 <                   j2 != CellList[m2].end(); ++j2) {
638 >            for (vector<int>::iterator j1 = cellList_[m1].begin();
639 >                 j1 != cellList_[m1].end(); ++j1) {
640 >              for (vector<int>::iterator j2 = cellList_[m2].begin();
641 >                   j2 != cellList_[m2].end(); ++j2) {
642                                
643                  // Always do this if we're in different cells or if
644                  // we're in the same cell and the global index of the
# Line 607 | Line 658 | namespace OpenMD {
658          }
659        }
660      }
661 +
662 +    // save the local cutoff group positions for the check that is
663 +    // done on each loop:
664 +    saved_CG_positions_.clear();
665 +    for (int i = 0; i < nGroups_; i++)
666 +      saved_CG_positions_.push_back(snap_->cgData.position[i]);
667 +
668      return neighborList;
669    }
670   } //end namespace OpenMD

Diff Legend

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