ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/io/DumpWriter.cpp
Revision: 1313
Committed: Wed Oct 22 20:01:49 2008 UTC (16 years, 6 months ago) by gezelter
Original Path: trunk/src/io/DumpWriter.cpp
File size: 11813 byte(s)
Log Message:
General bug-fixes and other changes to make particle pots work with
the Helfand Energy correlation function

File Contents

# User Rev Content
1 gezelter 507 /*
2 gezelter 246 * 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. Acknowledgement of the program authors must be made in any
10     * publication of scientific results based in part on use of the
11     * program. An acceptable form of acknowledgement is citation of
12     * the article in which the program was described (Matthew
13     * A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14     * J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15     * Parallel Simulation Engine for Molecular Dynamics,"
16     * J. Comput. Chem. 26, pp. 252-271 (2005))
17     *
18     * 2. Redistributions of source code must retain the above copyright
19     * notice, this list of conditions and the following disclaimer.
20     *
21     * 3. Redistributions in binary form must reproduce the above copyright
22     * notice, this list of conditions and the following disclaimer in the
23     * documentation and/or other materials provided with the
24     * distribution.
25     *
26     * This software is provided "AS IS," without a warranty of any
27     * kind. All express or implied conditions, representations and
28     * warranties, including any implied warranty of merchantability,
29     * fitness for a particular purpose or non-infringement, are hereby
30     * excluded. The University of Notre Dame and its licensors shall not
31     * be liable for any damages suffered by licensee as a result of
32     * using, modifying or distributing the software or its
33     * derivatives. In no event will the University of Notre Dame or its
34     * licensors be liable for any lost revenue, profit or data, or for
35     * direct, indirect, special, consequential, incidental or punitive
36     * damages, however caused and regardless of the theory of liability,
37     * arising out of the use of or inability to use software, even if the
38     * University of Notre Dame has been advised of the possibility of
39     * such damages.
40     */
41    
42     #include "io/DumpWriter.hpp"
43     #include "primitives/Molecule.hpp"
44     #include "utils/simError.h"
45 tim 376 #include "io/basic_teebuf.hpp"
46 tim 615 #include "io/gzstream.hpp"
47     #include "io/Globals.hpp"
48    
49 gezelter 2 #ifdef IS_MPI
50     #include <mpi.h>
51 gezelter 246 #endif //is_mpi
52 gezelter 2
53 gezelter 246 namespace oopse {
54 gezelter 2
55 gezelter 507 DumpWriter::DumpWriter(SimInfo* info)
56     : info_(info), filename_(info->getDumpFileName()), eorFilename_(info->getFinalConfigFileName()){
57 tim 615
58     Globals* simParams = info->getSimParams();
59     needCompression_ = simParams->getCompressDumpFile();
60 chrisfen 726 needForceVector_ = simParams->getOutputForceVector();
61 chuckv 791 createDumpFile_ = true;
62 tim 619 #ifdef HAVE_LIBZ
63 tim 615 if (needCompression_) {
64 tim 1024 filename_ += ".gz";
65     eorFilename_ += ".gz";
66 tim 615 }
67 tim 619 #endif
68 tim 615
69 tim 376 #ifdef IS_MPI
70    
71 tim 1024 if (worldRank == 0) {
72 tim 376 #endif // is_mpi
73 chuckv 791
74 tim 1024 dumpFile_ = createOStream(filename_);
75 tim 615
76 tim 1024 if (!dumpFile_) {
77     sprintf(painCave.errMsg, "Could not open \"%s\" for dump output.\n",
78     filename_.c_str());
79     painCave.isFatal = 1;
80     simError();
81     }
82 tim 376
83     #ifdef IS_MPI
84    
85 tim 1024 }
86 tim 376
87     #endif // is_mpi
88    
89 tim 1024 }
90 tim 376
91    
92 gezelter 507 DumpWriter::DumpWriter(SimInfo* info, const std::string& filename)
93     : info_(info), filename_(filename){
94 tim 615
95     Globals* simParams = info->getSimParams();
96     eorFilename_ = filename_.substr(0, filename_.rfind(".")) + ".eor";
97    
98     needCompression_ = simParams->getCompressDumpFile();
99 chrisfen 726 needForceVector_ = simParams->getOutputForceVector();
100 chuckv 791 createDumpFile_ = true;
101 tim 619 #ifdef HAVE_LIBZ
102 tim 615 if (needCompression_) {
103 tim 1024 filename_ += ".gz";
104     eorFilename_ += ".gz";
105 tim 615 }
106 tim 619 #endif
107 tim 615
108 gezelter 246 #ifdef IS_MPI
109 gezelter 2
110 tim 1024 if (worldRank == 0) {
111 gezelter 246 #endif // is_mpi
112 gezelter 2
113 chuckv 791
114 tim 1024 dumpFile_ = createOStream(filename_);
115 tim 615
116 tim 1024 if (!dumpFile_) {
117     sprintf(painCave.errMsg, "Could not open \"%s\" for dump output.\n",
118     filename_.c_str());
119     painCave.isFatal = 1;
120     simError();
121     }
122 gezelter 2
123     #ifdef IS_MPI
124    
125 tim 1024 }
126 gezelter 2
127     #endif // is_mpi
128 gezelter 246
129 tim 1024 }
130 chuckv 791
131     DumpWriter::DumpWriter(SimInfo* info, const std::string& filename, bool writeDumpFile)
132 tim 1024 : info_(info), filename_(filename){
133 chuckv 791
134     Globals* simParams = info->getSimParams();
135     eorFilename_ = filename_.substr(0, filename_.rfind(".")) + ".eor";
136    
137     needCompression_ = simParams->getCompressDumpFile();
138     needForceVector_ = simParams->getOutputForceVector();
139    
140     #ifdef HAVE_LIBZ
141     if (needCompression_) {
142     filename_ += ".gz";
143     eorFilename_ += ".gz";
144     }
145     #endif
146    
147     #ifdef IS_MPI
148    
149     if (worldRank == 0) {
150     #endif // is_mpi
151    
152     createDumpFile_ = writeDumpFile;
153     if (createDumpFile_) {
154     dumpFile_ = createOStream(filename_);
155    
156     if (!dumpFile_) {
157     sprintf(painCave.errMsg, "Could not open \"%s\" for dump output.\n",
158     filename_.c_str());
159     painCave.isFatal = 1;
160     simError();
161     }
162     }
163     #ifdef IS_MPI
164    
165     }
166 tim 1024
167 chuckv 791
168     #endif // is_mpi
169    
170     }
171 gezelter 2
172 gezelter 507 DumpWriter::~DumpWriter() {
173 gezelter 2
174     #ifdef IS_MPI
175 gezelter 246
176     if (worldRank == 0) {
177 gezelter 2 #endif // is_mpi
178 chuckv 791 if (createDumpFile_){
179 tim 1024 writeClosing(*dumpFile_);
180 chuckv 791 delete dumpFile_;
181     }
182 gezelter 2 #ifdef IS_MPI
183 gezelter 246
184     }
185    
186 gezelter 2 #endif // is_mpi
187 gezelter 246
188 gezelter 507 }
189 gezelter 2
190 tim 1024 void DumpWriter::writeFrameProperties(std::ostream& os, Snapshot* s) {
191 gezelter 2
192 tim 1024 char buffer[1024];
193    
194     os << " <FrameData>\n";
195    
196     RealType currentTime = s->getTime();
197 gezelter 1025 sprintf(buffer, " Time: %.10g\n", currentTime);
198 tim 1024 os << buffer;
199    
200 gezelter 246 Mat3x3d hmat;
201     hmat = s->getHmat();
202 tim 1024 sprintf(buffer, " Hmat: {{ %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }}\n",
203     hmat(0, 0), hmat(1, 0), hmat(2, 0),
204     hmat(0, 1), hmat(1, 1), hmat(2, 1),
205     hmat(0, 2), hmat(1, 2), hmat(2, 2));
206     os << buffer;
207 gezelter 2
208 tim 1024 RealType chi = s->getChi();
209     RealType integralOfChiDt = s->getIntegralOfChiDt();
210     sprintf(buffer, " Thermostat: %.10g , %.10g\n", chi, integralOfChiDt);
211     os << buffer;
212 gezelter 246
213 tim 1024 Mat3x3d eta;
214     eta = s->getEta();
215     sprintf(buffer, " Barostat: {{ %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }, { %.10g, %.10g, %.10g }}\n",
216     eta(0, 0), eta(1, 0), eta(2, 0),
217     eta(0, 1), eta(1, 1), eta(2, 1),
218     eta(0, 2), eta(1, 2), eta(2, 2));
219     os << buffer;
220 gezelter 246
221 tim 1024 os << " </FrameData>\n";
222 gezelter 507 }
223 gezelter 2
224 gezelter 507 void DumpWriter::writeFrame(std::ostream& os) {
225 gezelter 246
226 tim 1024 #ifdef IS_MPI
227     MPI_Status istatus;
228     #endif
229 gezelter 246
230     Molecule* mol;
231     StuntDouble* integrableObject;
232     SimInfo::MoleculeIterator mi;
233     Molecule::IntegrableObjectIterator ii;
234 gezelter 2
235 gezelter 246 #ifndef IS_MPI
236 tim 1024 os << " <Snapshot>\n";
237 gezelter 1025
238 tim 1024 writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot());
239 gezelter 2
240 tim 1024 os << " <StuntDoubles>\n";
241 gezelter 246 for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
242 gezelter 2
243 xsun 1217
244     for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
245 tim 1024 integrableObject = mol->nextIntegrableObject(ii)) {
246 xsun 1217 os << prepareDumpLine(integrableObject);
247    
248 tim 1024 }
249     }
250     os << " </StuntDoubles>\n";
251 xsun 1217
252 tim 1024 os << " </Snapshot>\n";
253 gezelter 2
254 tim 1024 os.flush();
255     #else
256     //every node prepares the dump lines for integrable objects belong to itself
257     std::string buffer;
258     for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
259 xsun 1217
260    
261 tim 1024 for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
262     integrableObject = mol->nextIntegrableObject(ii)) {
263 xsun 1217 buffer += prepareDumpLine(integrableObject);
264 gezelter 507 }
265 gezelter 2 }
266 tim 1024
267 gezelter 246 const int masterNode = 0;
268 gezelter 2
269 tim 1024 if (worldRank == masterNode) {
270     os << " <Snapshot>\n";
271     writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot());
272     os << " <StuntDoubles>\n";
273    
274 gezelter 1025 os << buffer;
275    
276 tim 1024 int nProc;
277 gezelter 507 MPI_Comm_size(MPI_COMM_WORLD, &nProc);
278 tim 1024 for (int i = 1; i < nProc; ++i) {
279 gezelter 2
280 tim 1024 // receive the length of the string buffer that was
281     // prepared by processor i
282 gezelter 2
283 tim 1024 int recvLength;
284     MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, &istatus);
285     char* recvBuffer = new char[recvLength];
286     if (recvBuffer == NULL) {
287     } else {
288 gezelter 1025 MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, &istatus);
289 tim 1024 os << recvBuffer;
290 gezelter 1313 delete [] recvBuffer;
291 tim 1024 }
292     }
293     os << " </StuntDoubles>\n";
294    
295     os << " </Snapshot>\n";
296 gezelter 507 os.flush();
297 gezelter 246 } else {
298 gezelter 1025 int sendBufferLength = buffer.size() + 1;
299 tim 1024 MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
300 gezelter 1025 MPI_Send((void *)buffer.c_str(), sendBufferLength, MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
301 tim 1024 }
302 gezelter 2
303 tim 1024 #endif // is_mpi
304 gezelter 2
305 tim 1024 }
306 gezelter 2
307 tim 1024 std::string DumpWriter::prepareDumpLine(StuntDouble* integrableObject) {
308    
309     int index = integrableObject->getGlobalIntegrableObjectIndex();
310     std::string type("pv");
311     std::string line;
312     char tempBuffer[4096];
313 gezelter 2
314 tim 1024 Vector3d pos;
315     Vector3d vel;
316     pos = integrableObject->getPos();
317     vel = integrableObject->getVel();
318 gezelter 1025 sprintf(tempBuffer, "%18.10g %18.10g %18.10g %13e %13e %13e",
319 tim 1024 pos[0], pos[1], pos[2],
320     vel[0], vel[1], vel[2]);
321     line += tempBuffer;
322 gezelter 2
323 tim 1024 if (integrableObject->isDirectional()) {
324     type += "qj";
325     Quat4d q;
326     Vector3d ji;
327     q = integrableObject->getQ();
328     ji = integrableObject->getJ();
329 gezelter 1025 sprintf(tempBuffer, " %13e %13e %13e %13e %13e %13e %13e",
330 tim 1024 q[0], q[1], q[2], q[3],
331     ji[0], ji[1], ji[2]);
332     line += tempBuffer;
333     }
334 gezelter 2
335 tim 1024 if (needForceVector_) {
336     type += "ft";
337     Vector3d frc;
338     Vector3d trq;
339     frc = integrableObject->getFrc();
340     trq = integrableObject->getTrq();
341    
342 gezelter 1025 sprintf(tempBuffer, " %13e %13e %13e %13e %13e %13e",
343 tim 1024 frc[0], frc[1], frc[2],
344     trq[0], trq[1], trq[2]);
345     line += tempBuffer;
346 gezelter 246 }
347 tim 1024
348 gezelter 1025 sprintf(tempBuffer, "%10d %7s %s\n", index, type.c_str(), line.c_str());
349 tim 1024 return std::string(tempBuffer);
350 gezelter 507 }
351 gezelter 2
352 gezelter 507 void DumpWriter::writeDump() {
353 tim 615 writeFrame(*dumpFile_);
354 gezelter 507 }
355 tim 376
356 gezelter 507 void DumpWriter::writeEor() {
357 tim 615 std::ostream* eorStream;
358 tim 376
359     #ifdef IS_MPI
360     if (worldRank == 0) {
361     #endif // is_mpi
362    
363 tim 615 eorStream = createOStream(eorFilename_);
364 tim 376
365     #ifdef IS_MPI
366     }
367     #endif // is_mpi
368    
369 tim 615 writeFrame(*eorStream);
370    
371     #ifdef IS_MPI
372     if (worldRank == 0) {
373     #endif // is_mpi
374 tim 1024 writeClosing(*eorStream);
375     delete eorStream;
376 tim 615 #ifdef IS_MPI
377     }
378     #endif // is_mpi
379    
380 gezelter 507 }
381 tim 376
382    
383 gezelter 507 void DumpWriter::writeDumpAndEor() {
384 tim 376 std::vector<std::streambuf*> buffers;
385 tim 615 std::ostream* eorStream;
386 tim 376 #ifdef IS_MPI
387     if (worldRank == 0) {
388     #endif // is_mpi
389    
390 tim 615 buffers.push_back(dumpFile_->rdbuf());
391 tim 376
392 tim 615 eorStream = createOStream(eorFilename_);
393 tim 376
394 tim 615 buffers.push_back(eorStream->rdbuf());
395 tim 376
396     #ifdef IS_MPI
397     }
398     #endif // is_mpi
399    
400     TeeBuf tbuf(buffers.begin(), buffers.end());
401     std::ostream os(&tbuf);
402    
403     writeFrame(os);
404 tim 615
405     #ifdef IS_MPI
406     if (worldRank == 0) {
407     #endif // is_mpi
408 tim 1024 writeClosing(*eorStream);
409     delete eorStream;
410 tim 615 #ifdef IS_MPI
411     }
412     #endif // is_mpi
413 tim 376
414 gezelter 507 }
415 tim 376
416 tim 1024 std::ostream* DumpWriter::createOStream(const std::string& filename) {
417 tim 619
418 tim 615 std::ostream* newOStream;
419 tim 619 #ifdef HAVE_LIBZ
420 tim 615 if (needCompression_) {
421 tim 1024 newOStream = new ogzstream(filename.c_str());
422 tim 615 } else {
423 tim 1024 newOStream = new std::ofstream(filename.c_str());
424 tim 615 }
425 tim 619 #else
426     newOStream = new std::ofstream(filename.c_str());
427     #endif
428 tim 1024 //write out MetaData first
429     (*newOStream) << "<OOPSE version=4>" << std::endl;
430     (*newOStream) << " <MetaData>" << std::endl;
431     (*newOStream) << info_->getRawMetaData();
432     (*newOStream) << " </MetaData>" << std::endl;
433 tim 615 return newOStream;
434 tim 1024 }
435 tim 376
436 tim 1024 void DumpWriter::writeClosing(std::ostream& os) {
437    
438     os << "</OOPSE>\n";
439     os.flush();
440     }
441    
442 gezelter 246 }//end namespace oopse