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 1538 by chuckv, Tue Jan 11 18:58:12 2011 UTC vs.
Revision 1541 by gezelter, Fri Feb 4 20:04:56 2011 UTC

# Line 1 | Line 1
1 < /**
2 < * @file ForceDecomposition.cpp
3 < * @author Charles Vardeman <cvardema.at.nd.edu>
4 < * @date 08/18/2010
5 < * @time 11:56am
6 < * @version 1.0
1 > /*
2 > * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   *
8 – * @section LICENSE
9 – * Copyright (c) 2010 The University of Notre Dame. All Rights Reserved.
10 – *
4   * The University of Notre Dame grants you ("Licensee") a
5   * non-exclusive, royalty free, license to use, modify and
6   * redistribute this software in source and binary code form, provided
# Line 45 | Line 38
38   * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
39   * [4]  Vardeman & Gezelter, in progress (2009).                        
40   */
41 + #include "parallel/ForceDecomposition.hpp"
42 + #include "parallel/Communicator.hpp"
43 + #include "math/SquareMatrix3.hpp"
44  
45 + using namespace std;
46 + namespace OpenMD {
47  
48 <
51 < /*  -*- c++ -*-  */
52 < #include "config.h"
53 < #include <stdlib.h>
48 >  void ForceDecomposition::distributeInitialData() {
49   #ifdef IS_MPI
50 < #include <mpi.h>
51 < #endif
50 >    Snapshot* snap = sman_->getCurrentSnapshot();
51 >    int nAtoms = snap->getNumberOfAtoms();
52 >    int nGroups = snap->getNumberOfCutoffGroups();
53  
54 < #include <iostream>
55 < #include <vector>
56 < #include <algorithm>
61 < #include <cmath>
62 < #include "parallel/ForceDecomposition.hpp"
54 >    AtomCommRealI = new Communicator<Row,RealType>(nAtoms);
55 >    AtomCommVectorI = new Communicator<Row,Vector3d>(nAtoms);
56 >    AtomCommMatrixI = new Communicator<Row,Mat3x3d>(nAtoms);
57  
58 +    AtomCommRealJ = new Communicator<Column,RealType>(nAtoms);
59 +    AtomCommVectorJ = new Communicator<Column,Vector3d>(nAtoms);
60 +    AtomCommMatrixJ = new Communicator<Column,Mat3x3d>(nAtoms);
61  
62 < using namespace std;
63 < using namespace OpenMD;
62 >    cgCommVectorI = new Communicator<Row,Vector3d>(nGroups);
63 >    cgCommVectorJ = new Communicator<Column,Vector3d>(nGroups);
64  
65 < //__static
66 < #ifdef IS_MPI
70 < static vector<MPI:Comm> communictors;
71 < #endif
65 >    int nInRow = AtomCommRealI.getSize();
66 >    int nInCol = AtomCommRealJ.getSize();
67  
68 < //____ MPITypeTraits
69 < template<typename T>
70 < struct MPITypeTraits;
68 >    vector<vector<RealType> > pot_row(LR_POT_TYPES,
69 >                                      vector<RealType> (nInRow, 0.0));
70 >    vector<vector<RealType> > pot_col(LR_POT_TYPES,
71 >                                      vector<RealType> (nInCol, 0.0));
72  
73 < #ifdef IS_MPI
74 < template<>
79 < struct MPITypeTraits<RealType> {
80 <  static const MPI::Datatype datatype;
81 < };
82 < const MPI_Datatype MPITypeTraits<RealType>::datatype = MY_MPI_REAL;
73 >    vector<vector<RealType> > pot_local(LR_POT_TYPES,
74 >                                        vector<RealType> (nAtoms, 0.0));
75  
84 – template<>
85 – struct MPITypeTraits<int> {
86 –  static const MPI::Datatype datatype;
87 – };
88 – const MPI::Datatype MPITypeTraits<int>::datatype = MPI_INT;
76   #endif
77 +  }
78 +    
79  
91 – /**
92 – * Constructor for ForceDecomposition Parallel Decomposition Method
93 – * Will try to construct a symmetric grid of processors. Ideally, the
94 – * number of processors will be a square ex: 4, 9, 16, 25.
95 – *
96 – */
80  
81 < ForceDecomposition::ForceDecomposition() {
99 <
81 >  void ForceDecomposition::distributeData()  {
82   #ifdef IS_MPI
83 <  int nProcs = MPI::COMM_WORLD.Get_size();
84 <  int worldRank = MPI::COMM_WORLD.Get_rank();
85 < #endif
86 <
87 <  // First time through, construct column stride.
88 <  if (communicators.size() == 0)
89 <  {
90 <    int nColumnsMax = (int) round(sqrt((float) nProcs));
91 <    for (int i = 0; i < nProcs; ++i)
92 <    {
93 <      if (nProcs%i==0) nColumns=i;
83 >    Snapshot* snap = sman_->getCurrentSnapshot();
84 >    
85 >    // gather up the atomic positions
86 >    AtomCommVectorI->gather(snap->atomData.position,
87 >                            snap->atomIData.position);
88 >    AtomCommVectorJ->gather(snap->atomData.position,
89 >                            snap->atomJData.position);
90 >    
91 >    // gather up the cutoff group positions
92 >    cgCommVectorI->gather(snap->cgData.position,
93 >                          snap->cgIData.position);
94 >    cgCommVectorJ->gather(snap->cgData.position,
95 >                          snap->cgJData.position);
96 >    
97 >    // if needed, gather the atomic rotation matrices
98 >    if (snap->atomData.getStorageLayout() & DataStorage::dslAmat) {
99 >      AtomCommMatrixI->gather(snap->atomData.aMat,
100 >                              snap->atomIData.aMat);
101 >      AtomCommMatrixJ->gather(snap->atomData.aMat,
102 >                              snap->atomJData.aMat);
103      }
104 +    
105 +    // if needed, gather the atomic eletrostatic frames
106 +    if (snap->atomData.getStorageLayout() & DataStorage::dslElectroFrame) {
107 +      AtomCommMatrixI->gather(snap->atomData.electroFrame,
108 +                              snap->atomIData.electroFrame);
109 +      AtomCommMatrixJ->gather(snap->atomData.electroFrame,
110 +                              snap->atomJData.electroFrame);
111 +    }
112 + #endif      
113 +  }
114 +  
115 +  void ForceDecomposition::collectIntermediateData() {
116 + #ifdef IS_MPI
117 +    Snapshot* snap = sman_->getCurrentSnapshot();
118 +    
119 +    if (snap->atomData.getStorageLayout() & DataStorage::dslDensity) {
120  
121 <    int nRows = nProcs/nColumns;    
122 <    myRank_ = (int) worldRank%nColumns;
121 >      AtomCommRealI->scatter(snap->atomIData.density,
122 >                             snap->atomData.density);
123 >
124 >      int n = snap->atomData.density.size();
125 >      std::vector<RealType> rho_tmp(n, 0.0);
126 >      AtomCommRealJ->scatter(snap->atomJData.density, rho_tmp);
127 >      for (int i = 0; i < n; i++)
128 >        snap->atomData.density[i] += rho_tmp[i];
129 >    }
130 > #endif
131    }
132 <  else
133 <  {
134 <    myRank_ = myRank/nColumns;
132 >  
133 >  void ForceDecomposition::distributeIntermediateData() {
134 > #ifdef IS_MPI
135 >    Snapshot* snap = sman_->getCurrentSnapshot();
136 >    if (snap->atomData.getStorageLayout() & DataStorage::dslFunctional) {
137 >      AtomCommRealI->gather(snap->atomData.functional,
138 >                            snap->atomIData.functional);
139 >      AtomCommRealJ->gather(snap->atomData.functional,
140 >                            snap->atomJData.functional);
141 >    }
142 >    
143 >    if (snap->atomData.getStorageLayout() & DataStorage::dslFunctionalDerivative) {
144 >      AtomCommRealI->gather(snap->atomData.functionalDerivative,
145 >                            snap->atomIData.functionalDerivative);
146 >      AtomCommRealJ->gather(snap->atomData.functionalDerivative,
147 >                            snap->atomJData.functionalDerivative);
148 >    }
149 > #endif
150    }
121 –  MPI::Comm newComm = MPI:COMM_WORLD.Split(myRank_,0);
151    
123 –  isColumn_ = false;
152    
153 < }
153 >  void ForceDecomposition::collectData() {
154 > #ifdef IS_MPI
155 >    Snapshot* snap = sman_->getCurrentSnapshot();
156 >    
157 >    int n = snap->atomData.force.size();
158 >    std::vector<Vector3d> frc_tmp(n, 0.0);
159 >    
160 >    AtomCommVectorI->scatter(snap->atomIData.force, frc_tmp);
161 >    for (int i = 0; i < n; i++) {
162 >      snap->atomData.force[i] += frc_tmp[i];
163 >      frc_tmp[i] = 0.0;
164 >    }
165 >    
166 >    AtomCommVectorJ->scatter(snap->atomJData.force, frc_tmp);
167 >    for (int i = 0; i < n; i++)
168 >      snap->atomData.force[i] += frc_tmp[i];
169 >    
170 >    
171 >    if (snap->atomData.getStorageLayout() & DataStorage::dslTorque) {
172  
173 < ForceDecomposition::gather(sendbuf, receivebuf){
174 <  communicators(myIndex_).Allgatherv();
129 < }
173 >      int nt = snap->atomData.force.size();
174 >      std::vector<Vector3d> trq_tmp(nt, 0.0);
175  
176 +      AtomCommVectorI->scatter(snap->atomIData.torque, trq_tmp);
177 +      for (int i = 0; i < n; i++) {
178 +        snap->atomData.torque[i] += trq_tmp[i];
179 +        trq_tmp[i] = 0.0;
180 +      }
181 +      
182 +      AtomCommVectorJ->scatter(snap->atomJData.torque, trq_tmp);
183 +      for (int i = 0; i < n; i++)
184 +        snap->atomData.torque[i] += trq_tmp[i];
185 +    }
186 +    
187 +    
188 +    vector<vector<RealType> > pot_temp(LR_POT_TYPES,
189 +                                       vector<RealType> (nAtoms, 0.0));
190 +    
191 +    for (int i = 0; i < LR_POT_TYPES; i++) {
192 +      AtomCommRealI->scatter(pot_row[i], pot_temp[i]);
193 +      for (int ii = 0;  ii < pot_temp[i].size(); ii++ ) {
194 +        pot_local[i] += pot_temp[i][ii];
195 +      }
196 +    }
197 +    
198  
199  
200 < ForceDecomposition::scatter(sbuffer, rbuffer){
201 <  communicators(myIndex_).Reduce_scatter(sbuffer, recevbuf. recvcounts, MPI::DOUBLE, MPI::SUM);
202 < }
203 <
137 <
200 > #endif
201 >  }
202 >  
203 > } //end namespace OpenMD

Diff Legend

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