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/ForceDecomposition.cpp (file contents):
Revision 1539 by gezelter, Fri Jan 14 22:31:31 2011 UTC vs.
Revision 1547 by gezelter, Mon Apr 11 18:44:16 2011 UTC

# Line 39 | Line 39
39   * [4]  Vardeman & Gezelter, in progress (2009).                        
40   */
41   #include "parallel/ForceDecomposition.hpp"
42 #include "parallel/Communicator.hpp"
42   #include "math/SquareMatrix3.hpp"
43 + #include "nonbonded/NonBondedInteraction.hpp"
44 + #include "brains/SnapshotManager.hpp"
45  
46 + using namespace std;
47   namespace OpenMD {
48  
49 +  /**
50 +   * distributeInitialData is essentially a copy of the older fortran
51 +   * SimulationSetup
52 +   */
53 +  
54    void ForceDecomposition::distributeInitialData() {
55 < #ifdef IS_MPI
55 > #ifdef IS_MPI    
56 >    Snapshot* snap = sman_->getCurrentSnapshot();
57 >    int nLocal = snap->getNumberOfAtoms();
58 >    int nGroups = snap->getNumberOfCutoffGroups();
59  
60 <    int nAtoms;
61 <    int nGroups;
60 >    AtomCommIntI = new Communicator<Row,int>(nLocal);
61 >    AtomCommRealI = new Communicator<Row,RealType>(nLocal);
62 >    AtomCommVectorI = new Communicator<Row,Vector3d>(nLocal);
63 >    AtomCommMatrixI = new Communicator<Row,Mat3x3d>(nLocal);
64  
65 <    AtomCommRealI = new Comm<I,RealType>(nAtoms);
66 <    AtomCommVectorI = new Comm<I,Vector3d>(nAtoms);
67 <    AtomCommMatrixI = new Comm<I,Mat3x3d>(nAtoms);
65 >    AtomCommIntJ = new Communicator<Column,int>(nLocal);
66 >    AtomCommRealJ = new Communicator<Column,RealType>(nLocal);
67 >    AtomCommVectorJ = new Communicator<Column,Vector3d>(nLocal);
68 >    AtomCommMatrixJ = new Communicator<Column,Mat3x3d>(nLocal);
69  
70 <    AtomCommRealJ = new Comm<J,RealType>(nAtoms);
71 <    AtomCommVectorJ = new Comm<J,Vector3d>(nAtoms);
72 <    AtomCommMatrixJ = new Comm<J,Mat3x3d>(nAtoms);
70 >    cgCommIntI = new Communicator<Row,int>(nGroups);
71 >    cgCommVectorI = new Communicator<Row,Vector3d>(nGroups);
72 >    cgCommIntJ = new Communicator<Column,int>(nGroups);
73 >    cgCommVectorJ = new Communicator<Column,Vector3d>(nGroups);
74  
75 <    cgCommVectorI = new Comm<I,Vector3d>(nGroups);
76 <    cgCommVectorJ = new Comm<J,Vector3d>(nGroups);
77 <    // more to come
75 >    int nAtomsInRow = AtomCommIntI->getSize();
76 >    int nAtomsInCol = AtomCommIntJ->getSize();
77 >    int nGroupsInRow = cgCommIntI->getSize();
78 >    int nGroupsInCol = cgCommIntJ->getSize();
79 >
80 >    vector<vector<RealType> > pot_row(N_INTERACTION_FAMILIES,
81 >                                      vector<RealType> (nAtomsInRow, 0.0));
82 >    vector<vector<RealType> > pot_col(N_INTERACTION_FAMILIES,
83 >                                      vector<RealType> (nAtomsInCol, 0.0));
84 >    
85 >    vector<RealType> pot_local(N_INTERACTION_FAMILIES, 0.0);
86 >
87 >    // gather the information for atomtype IDs (atids):
88 >    vector<int> identsLocal = info_->getIdentArray();
89 >    identsRow.reserve(nAtomsInRow);
90 >    identsCol.reserve(nAtomsInCol);
91 >
92 >    AtomCommIntI->gather(identsLocal, identsRow);
93 >    AtomCommIntJ->gather(identsLocal, identsCol);
94 >
95 >    AtomLocalToGlobal = info_->getLocalToGlobalAtomIndex();
96 >    AtomCommIntI->gather(AtomLocalToGlobal, AtomRowToGlobal);
97 >    AtomCommIntJ->gather(AtomLocalToGlobal, AtomColToGlobal);
98 >
99 >    cgLocalToGlobal = info_->getLocalToGlobalCutoffGroupIndex();
100 >    cgCommIntI->gather(cgLocalToGlobal, cgRowToGlobal);
101 >    cgCommIntJ->gather(cgLocalToGlobal, cgColToGlobal);
102 >
103 >    
104 >
105 >    // still need:
106 >    // topoDist
107 >    // exclude
108   #endif
109    }
110      
# Line 69 | Line 113 | namespace OpenMD {
113    void ForceDecomposition::distributeData()  {
114   #ifdef IS_MPI
115      Snapshot* snap = sman_->getCurrentSnapshot();
116 <
116 >    
117      // gather up the atomic positions
118      AtomCommVectorI->gather(snap->atomData.position,
119                              snap->atomIData.position);
120      AtomCommVectorJ->gather(snap->atomData.position,
121 <                           snap->atomJData.position);
121 >                            snap->atomJData.position);
122      
123      // gather up the cutoff group positions
124      cgCommVectorI->gather(snap->cgData.position,
125 <                         snap->cgIData.position);
125 >                          snap->cgIData.position);
126      cgCommVectorJ->gather(snap->cgData.position,
127 <                         snap->cgJData.position);
127 >                          snap->cgJData.position);
128      
129      // if needed, gather the atomic rotation matrices
130      if (snap->atomData.getStorageLayout() & DataStorage::dslAmat) {
131        AtomCommMatrixI->gather(snap->atomData.aMat,
132 <                             snap->atomIData.aMat);
132 >                              snap->atomIData.aMat);
133        AtomCommMatrixJ->gather(snap->atomData.aMat,
134 <                             snap->atomJData.aMat);
134 >                              snap->atomJData.aMat);
135      }
136      
137      // if needed, gather the atomic eletrostatic frames
138      if (snap->atomData.getStorageLayout() & DataStorage::dslElectroFrame) {
139        AtomCommMatrixI->gather(snap->atomData.electroFrame,
140 <                             snap->atomIData.electroFrame);
140 >                              snap->atomIData.electroFrame);
141        AtomCommMatrixJ->gather(snap->atomData.electroFrame,
142 <                             snap->atomJData.electroFrame);
142 >                              snap->atomJData.electroFrame);
143      }
144   #endif      
145    }
# Line 103 | Line 147 | namespace OpenMD {
147    void ForceDecomposition::collectIntermediateData() {
148   #ifdef IS_MPI
149      Snapshot* snap = sman_->getCurrentSnapshot();
106    // gather up the atomic positions
150      
151      if (snap->atomData.getStorageLayout() & DataStorage::dslDensity) {
152 +
153        AtomCommRealI->scatter(snap->atomIData.density,
154 <                            snap->atomData.density);
155 <      std::vector<RealType> rho_tmp;
156 <      int n = snap->getNumberOfAtoms();
157 <      rho_tmp.reserve( n );
154 >                             snap->atomData.density);
155 >
156 >      int n = snap->atomData.density.size();
157 >      std::vector<RealType> rho_tmp(n, 0.0);
158        AtomCommRealJ->scatter(snap->atomJData.density, rho_tmp);
159        for (int i = 0; i < n; i++)
160          snap->atomData.density[i] += rho_tmp[i];
# Line 123 | Line 167 | namespace OpenMD {
167      Snapshot* snap = sman_->getCurrentSnapshot();
168      if (snap->atomData.getStorageLayout() & DataStorage::dslFunctional) {
169        AtomCommRealI->gather(snap->atomData.functional,
170 <                           snap->atomIData.functional);
170 >                            snap->atomIData.functional);
171        AtomCommRealJ->gather(snap->atomData.functional,
172 <                           snap->atomJData.functional);
172 >                            snap->atomJData.functional);
173      }
174      
175      if (snap->atomData.getStorageLayout() & DataStorage::dslFunctionalDerivative) {
176        AtomCommRealI->gather(snap->atomData.functionalDerivative,
177 <                           snap->atomIData.functionalDerivative);
177 >                            snap->atomIData.functionalDerivative);
178        AtomCommRealJ->gather(snap->atomData.functionalDerivative,
179 <                           snap->atomJData.functionalDerivative);
179 >                            snap->atomJData.functionalDerivative);
180      }
181   #endif
182    }
# Line 140 | Line 184 | namespace OpenMD {
184    
185    void ForceDecomposition::collectData() {
186   #ifdef IS_MPI
187 +    Snapshot* snap = sman_->getCurrentSnapshot();
188 +    
189 +    int n = snap->atomData.force.size();
190 +    vector<Vector3d> frc_tmp(n, V3Zero);
191 +    
192 +    AtomCommVectorI->scatter(snap->atomIData.force, frc_tmp);
193 +    for (int i = 0; i < n; i++) {
194 +      snap->atomData.force[i] += frc_tmp[i];
195 +      frc_tmp[i] = 0.0;
196 +    }
197 +    
198 +    AtomCommVectorJ->scatter(snap->atomJData.force, frc_tmp);
199 +    for (int i = 0; i < n; i++)
200 +      snap->atomData.force[i] += frc_tmp[i];
201 +    
202 +    
203 +    if (snap->atomData.getStorageLayout() & DataStorage::dslTorque) {
204 +
205 +      int nt = snap->atomData.force.size();
206 +      vector<Vector3d> trq_tmp(nt, V3Zero);
207 +
208 +      AtomCommVectorI->scatter(snap->atomIData.torque, trq_tmp);
209 +      for (int i = 0; i < n; i++) {
210 +        snap->atomData.torque[i] += trq_tmp[i];
211 +        trq_tmp[i] = 0.0;
212 +      }
213 +      
214 +      AtomCommVectorJ->scatter(snap->atomJData.torque, trq_tmp);
215 +      for (int i = 0; i < n; i++)
216 +        snap->atomData.torque[i] += trq_tmp[i];
217 +    }
218 +    
219 +    int nLocal = snap->getNumberOfAtoms();
220 +
221 +    vector<vector<RealType> > pot_temp(N_INTERACTION_FAMILIES,
222 +                                       vector<RealType> (nLocal, 0.0));
223 +    
224 +    for (int i = 0; i < N_INTERACTION_FAMILIES; i++) {
225 +      AtomCommRealI->scatter(pot_row[i], pot_temp[i]);
226 +      for (int ii = 0;  ii < pot_temp[i].size(); ii++ ) {
227 +        pot_local[i] += pot_temp[i][ii];
228 +      }
229 +    }
230   #endif
231    }
232    

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines