ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/parallel/ForceMatrixDecomposition.cpp
Revision: 1547
Committed: Mon Apr 11 18:44:16 2011 UTC (14 years ago) by gezelter
Original Path: branches/development/src/parallel/ForceDecomposition.cpp
File size: 8836 byte(s)
Log Message:
more minor changes to help build

File Contents

# Content
1 /*
2 * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3 *
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
7 * that the following conditions are met:
8 *
9 * 1. Redistributions of source code must retain the above copyright
10 * notice, this list of conditions and the following disclaimer.
11 *
12 * 2. Redistributions in binary form must reproduce the above copyright
13 * notice, this list of conditions and the following disclaimer in the
14 * documentation and/or other materials provided with the
15 * distribution.
16 *
17 * This software is provided "AS IS," without a warranty of any
18 * kind. All express or implied conditions, representations and
19 * warranties, including any implied warranty of merchantability,
20 * fitness for a particular purpose or non-infringement, are hereby
21 * excluded. The University of Notre Dame and its licensors shall not
22 * be liable for any damages suffered by licensee as a result of
23 * using, modifying or distributing the software or its
24 * derivatives. In no event will the University of Notre Dame or its
25 * licensors be liable for any lost revenue, profit or data, or for
26 * direct, indirect, special, consequential, incidental or punitive
27 * damages, however caused and regardless of the theory of liability,
28 * arising out of the use of or inability to use software, even if the
29 * University of Notre Dame has been advised of the possibility of
30 * such damages.
31 *
32 * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your
33 * research, please cite the appropriate papers when you publish your
34 * work. Good starting points are:
35 *
36 * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
37 * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
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 "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
56 Snapshot* snap = sman_->getCurrentSnapshot();
57 int nLocal = snap->getNumberOfAtoms();
58 int nGroups = snap->getNumberOfCutoffGroups();
59
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 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 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 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
111
112
113 void ForceDecomposition::distributeData() {
114 #ifdef IS_MPI
115 Snapshot* snap = sman_->getCurrentSnapshot();
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);
122
123 // gather up the cutoff group positions
124 cgCommVectorI->gather(snap->cgData.position,
125 snap->cgIData.position);
126 cgCommVectorJ->gather(snap->cgData.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);
133 AtomCommMatrixJ->gather(snap->atomData.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);
141 AtomCommMatrixJ->gather(snap->atomData.electroFrame,
142 snap->atomJData.electroFrame);
143 }
144 #endif
145 }
146
147 void ForceDecomposition::collectIntermediateData() {
148 #ifdef IS_MPI
149 Snapshot* snap = sman_->getCurrentSnapshot();
150
151 if (snap->atomData.getStorageLayout() & DataStorage::dslDensity) {
152
153 AtomCommRealI->scatter(snap->atomIData.density,
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];
161 }
162 #endif
163 }
164
165 void ForceDecomposition::distributeIntermediateData() {
166 #ifdef IS_MPI
167 Snapshot* snap = sman_->getCurrentSnapshot();
168 if (snap->atomData.getStorageLayout() & DataStorage::dslFunctional) {
169 AtomCommRealI->gather(snap->atomData.functional,
170 snap->atomIData.functional);
171 AtomCommRealJ->gather(snap->atomData.functional,
172 snap->atomJData.functional);
173 }
174
175 if (snap->atomData.getStorageLayout() & DataStorage::dslFunctionalDerivative) {
176 AtomCommRealI->gather(snap->atomData.functionalDerivative,
177 snap->atomIData.functionalDerivative);
178 AtomCommRealJ->gather(snap->atomData.functionalDerivative,
179 snap->atomJData.functionalDerivative);
180 }
181 #endif
182 }
183
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
233 } //end namespace OpenMD