ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-3.0/src/io/DumpReader.cpp
Revision: 1886
Committed: Wed Dec 15 16:13:35 2004 UTC (19 years, 9 months ago) by tim
File size: 16201 byte(s)
Log Message:
MPI in NVE is working

File Contents

# User Rev Content
1 gezelter 1490 #define _LARGEFILE_SOURCE64
2     #define _FILE_OFFSET_BITS 64
3    
4     #include <sys/types.h>
5     #include <sys/stat.h>
6    
7     #include <iostream>
8     #include <math.h>
9    
10     #include <stdio.h>
11     #include <stdlib.h>
12     #include <string.h>
13    
14 tim 1739 #include "io/DumpReader.hpp"
15 tim 1808 #include "primitives/Molecule.hpp"
16 tim 1492 #include "utils/simError.h"
17 tim 1784 #include "utils/MemoryUtils.hpp"
18 tim 1808 #include "utils/StringTokenizer.hpp"
19    
20 tim 1739 #ifdef IS_MPI
21 gezelter 1490
22     #include <mpi.h>
23     #define TAKE_THIS_TAG_CHAR 0
24     #define TAKE_THIS_TAG_INT 1
25 tim 1739
26 gezelter 1490 #endif // is_mpi
27    
28    
29 tim 1739 namespace oopse {
30 gezelter 1490
31 tim 1739 DumpReader::DumpReader(SimInfo* info, const std::string& filename)
32 tim 1885 : info_(info), filename_(filename), isScanned_(false), nframes_(0) {
33 gezelter 1490
34     #ifdef IS_MPI
35 tim 1739
36     if (worldRank == 0) {
37 gezelter 1490 #endif
38    
39 tim 1739 inFile_ = fopen(filename_.c_str(), "r");
40    
41     if (inFile_ == NULL) {
42 tim 1784 sprintf(painCave.errMsg, "Cannot open file: %s\n", filename_.c_str());
43 tim 1739 painCave.isFatal = 1;
44     simError();
45     }
46    
47 gezelter 1490 #ifdef IS_MPI
48 tim 1739
49     }
50    
51     strcpy(checkPointMsg, "Dump file opened for reading successfully.");
52     MPIcheckPoint();
53    
54 gezelter 1490 #endif
55 tim 1739
56     return;
57 gezelter 1490 }
58    
59 tim 1739 DumpReader::~DumpReader() {
60    
61 gezelter 1490 #ifdef IS_MPI
62 tim 1739
63     if (worldRank == 0) {
64 gezelter 1490 #endif
65    
66 tim 1739 int error;
67     error = fclose(inFile_);
68 gezelter 1490
69 tim 1739 if (error) {
70     sprintf(painCave.errMsg, "Error closing %s\n", filename_.c_str());
71     painCave.isFatal = 1;
72     simError();
73     }
74    
75     MemoryUtils::deleteVectorOfPointer(framePos_);
76    
77 gezelter 1490 #ifdef IS_MPI
78 tim 1739
79     }
80    
81     strcpy(checkPointMsg, "Dump file closed successfully.");
82     MPIcheckPoint();
83    
84 gezelter 1490 #endif
85    
86 tim 1739 return;
87 gezelter 1490 }
88    
89 tim 1739 int DumpReader::getNFrames(void) {
90 gezelter 1490
91 tim 1886 if (!isScanned_)
92     scanFile();
93 tim 1885
94 tim 1886 return nframes_;
95 gezelter 1490 }
96    
97 tim 1739 void DumpReader::scanFile(void) {
98 gezelter 1490 int i, j;
99     int lineNum = 0;
100 tim 1739 char readBuffer[maxBufferSize];
101     fpos_t * currPos;
102 gezelter 1490
103     #ifdef IS_MPI
104    
105 tim 1739 if (worldRank == 0) {
106 gezelter 1490 #endif // is_mpi
107    
108 tim 1739 rewind(inFile_);
109 gezelter 1490
110 tim 1739 currPos = new fpos_t;
111     fgetpos(inFile_, currPos);
112     fgets(readBuffer, sizeof(readBuffer), inFile_);
113     lineNum++;
114 gezelter 1490
115 tim 1739 if (feof(inFile_)) {
116     sprintf(painCave.errMsg,
117     "File \"%s\" ended unexpectedly at line %d\n",
118     filename_.c_str(),
119     lineNum);
120     painCave.isFatal = 1;
121     simError();
122     }
123 gezelter 1490
124 tim 1739 while (!feof(inFile_)) {
125     framePos_.push_back(currPos);
126 gezelter 1490
127 tim 1739 i = atoi(readBuffer);
128 gezelter 1490
129 tim 1739 fgets(readBuffer, sizeof(readBuffer), inFile_);
130     lineNum++;
131 gezelter 1490
132 tim 1739 if (feof(inFile_)) {
133     sprintf(painCave.errMsg,
134     "File \"%s\" ended unexpectedly at line %d\n",
135     filename_.c_str(),
136     lineNum);
137     painCave.isFatal = 1;
138     simError();
139     }
140 gezelter 1490
141 tim 1739 for(j = 0; j < i; j++) {
142     fgets(readBuffer, sizeof(readBuffer), inFile_);
143     lineNum++;
144 gezelter 1490
145 tim 1739 if (feof(inFile_)) {
146     sprintf(painCave.errMsg,
147     "File \"%s\" ended unexpectedly at line %d,"
148     " with atom %d\n", filename_.c_str(),
149     lineNum,
150     j);
151 gezelter 1490
152 tim 1739 painCave.isFatal = 1;
153     simError();
154     }
155     }
156 gezelter 1490
157 tim 1739 currPos = new fpos_t;
158     fgetpos(inFile_, currPos);
159     fgets(readBuffer, sizeof(readBuffer), inFile_);
160     lineNum++;
161     }
162 gezelter 1490
163 tim 1739 delete currPos;
164     rewind(inFile_);
165 tim 1886
166 tim 1885 nframes_ = framePos_.size();
167 tim 1739 #ifdef IS_MPI
168     }
169 gezelter 1490
170 tim 1886 MPI_Bcast(&nframes_, 1, MPI_INT, 0, MPI_COMM_WORLD);
171    
172 tim 1739 strcpy(checkPointMsg, "Successfully scanned DumpFile\n");
173     MPIcheckPoint();
174 gezelter 1490
175 tim 1739 #endif // is_mpi
176 gezelter 1490
177 tim 1886 isScanned_ = true;
178 tim 1739 }
179 gezelter 1490
180 tim 1739 void DumpReader::readFrame(int whichFrame) {
181     readSet(whichFrame);
182     }
183 gezelter 1490
184 tim 1739 void DumpReader::readSet(int whichFrame) {
185     int i;
186     int nTotObjs; // the number of atoms
187     char read_buffer[maxBufferSize]; //the line buffer for reading
188     char * eof_test; // ptr to see when we reach the end of the file
189 gezelter 1490
190 tim 1739 Molecule* mol;
191     StuntDouble* integrableObject;
192 tim 1808 SimInfo::MoleculeIterator mi;
193     Molecule::IntegrableObjectIterator ii;
194 gezelter 1490
195 tim 1739 #ifndef IS_MPI
196 gezelter 1490
197 tim 1739 fsetpos(inFile_, framePos_[whichFrame]);
198     eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
199 gezelter 1490
200 tim 1739 if (eof_test == NULL) {
201 gezelter 1490 sprintf(painCave.errMsg,
202 tim 1739 "DumpReader error: error reading 1st line of \"%s\"\n",
203     filename_.c_str());
204 gezelter 1490 painCave.isFatal = 1;
205     simError();
206     }
207    
208 tim 1739 nTotObjs = atoi(read_buffer);
209 gezelter 1490
210 tim 1739 if (nTotObjs != info_->getNGlobalIntegrableObjects()) {
211     sprintf(painCave.errMsg,
212     "DumpReader error. %s nIntegrable, %d, "
213     "does not match the meta-data file's nIntegrable, %d.\n",
214     filename_.c_str(),
215     nTotObjs,
216     info_->getNIntegrableObjects());
217 gezelter 1490
218 tim 1739 painCave.isFatal = 1;
219     simError();
220     }
221 gezelter 1490
222 tim 1739 //read the box mat from the comment line
223 gezelter 1490
224 tim 1739 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
225 gezelter 1490
226 tim 1739 if (eof_test == NULL) {
227     sprintf(painCave.errMsg, "error in reading commment in %s\n",
228     filename_.c_str());
229     painCave.isFatal = 1;
230     simError();
231 gezelter 1490 }
232    
233 tim 1808 parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
234 gezelter 1490
235 tim 1739 //parse dump lines
236 gezelter 1490
237 tim 1739 for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
238 gezelter 1490
239 tim 1808 for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
240     integrableObject = mol->nextIntegrableObject(ii)) {
241 gezelter 1490
242 tim 1739 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
243 gezelter 1490
244 tim 1739 if (eof_test == NULL) {
245     sprintf(painCave.errMsg,
246     "error in reading file %s\n"
247     "natoms = %d; index = %d\n"
248     "error reading the line from the file.\n",
249     filename_.c_str(),
250     nTotObjs,
251     i);
252 gezelter 1490
253 tim 1739 painCave.isFatal = 1;
254     simError();
255     }
256 gezelter 1490
257 tim 1808 parseDumpLine(read_buffer, integrableObject);
258 tim 1739
259 tim 1808 }
260 gezelter 1490 }
261    
262 tim 1739 // MPI Section of code..........
263 gezelter 1490
264 tim 1739 #else //IS_MPI
265 gezelter 1490
266 tim 1739 // first thing first, suspend fatalities.
267     int masterNode = 0;
268     int nCurObj;
269     painCave.isEventLoop = 1;
270 gezelter 1490
271 tim 1739 int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
272     int haveError;
273 gezelter 1490
274 tim 1739 MPI_Status istatus;
275     int nitems;
276 gezelter 1490
277 tim 1739 nTotObjs = info_->getNGlobalIntegrableObjects();
278     haveError = 0;
279 gezelter 1490
280 tim 1739 if (worldRank == masterNode) {
281     fsetpos(inFile_, framePos_[whichFrame]);
282 gezelter 1490
283 tim 1739 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
284 gezelter 1490
285 tim 1739 if (eof_test == NULL) {
286     sprintf(painCave.errMsg, "Error reading 1st line of %s \n ",
287     filename_.c_str());
288     painCave.isFatal = 1;
289     simError();
290     }
291 gezelter 1490
292 tim 1739 nitems = atoi(read_buffer);
293 gezelter 1490
294 tim 1739 // Check to see that the number of integrable objects in the
295     // intial configuration file is the same as derived from the
296     // meta-data file.
297 gezelter 1490
298 tim 1739 if (nTotObjs != nitems) {
299     sprintf(painCave.errMsg,
300     "DumpReader Error. %s nIntegrable, %d, "
301     "does not match the meta-data file's nIntegrable, %d.\n",
302     filename_.c_str(),
303     nTotObjs,
304 tim 1883 info_->getNGlobalIntegrableObjects());
305 gezelter 1490
306 tim 1739 painCave.isFatal = 1;
307     simError();
308 gezelter 1490 }
309    
310 tim 1739 //read the boxMat from the comment line
311 gezelter 1490
312 tim 1739 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
313 gezelter 1490
314 tim 1739 if (eof_test == NULL) {
315     sprintf(painCave.errMsg, "error in reading commment in %s\n",
316     filename_.c_str());
317     painCave.isFatal = 1;
318     simError();
319 gezelter 1490 }
320    
321 tim 1739 //Every single processor will parse the comment line by itself
322     //By using this way, we might lose some efficiency, but if we want to add
323     //more parameters into comment line, we only need to modify function
324     //parseCommentLine
325 gezelter 1490
326 tim 1739 MPI_Bcast(read_buffer, maxBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
327 tim 1808 parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
328 gezelter 1490
329 tim 1739 for(i = 0; i < info_->getNGlobalMolecules(); i++) {
330 tim 1883 int which_node = info_->getMolToProc(i);
331 gezelter 1490
332 tim 1739 if (which_node == masterNode) {
333     //molecules belong to master node
334 gezelter 1490
335 tim 1739 mol = info_->getMoleculeByGlobalIndex(i);
336 gezelter 1490
337 tim 1739 if (mol == NULL) {
338 tim 1884 sprintf(painCave.errMsg, "Molecule not found on node %d!", worldRank);
339 tim 1739 painCave.isFatal = 1;
340     simError();
341     }
342 gezelter 1490
343 tim 1739 for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
344     integrableObject = mol->nextIntegrableObject(ii)){
345    
346     eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
347 gezelter 1490
348 tim 1739 if (eof_test == NULL) {
349     sprintf(painCave.errMsg,
350     "error in reading file %s\n"
351     "natoms = %d; index = %d\n"
352     "error reading the line from the file.\n",
353     filename_.c_str(),
354     nTotObjs,
355     i);
356 gezelter 1490
357 tim 1739 painCave.isFatal = 1;
358     simError();
359     }
360 gezelter 1490
361 tim 1739 parseDumpLine(read_buffer, integrableObject);
362     }
363     } else {
364     //molecule belongs to slave nodes
365 gezelter 1490
366 tim 1739 MPI_Recv(&nCurObj, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT,
367     MPI_COMM_WORLD, &istatus);
368 gezelter 1490
369 tim 1739 for(int j = 0; j < nCurObj; j++) {
370     eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
371 gezelter 1490
372 tim 1739 if (eof_test == NULL) {
373     sprintf(painCave.errMsg,
374     "error in reading file %s\n"
375     "natoms = %d; index = %d\n"
376     "error reading the line from the file.\n",
377     filename_.c_str(),
378     nTotObjs,
379     i);
380 gezelter 1490
381 tim 1739 painCave.isFatal = 1;
382     simError();
383     }
384    
385     MPI_Send(read_buffer, maxBufferSize, MPI_CHAR, which_node,
386     TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
387     }
388     }
389     }
390     } else {
391     //actions taken at slave nodes
392     MPI_Bcast(read_buffer, maxBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
393 gezelter 1490
394 tim 1739 /**@todo*/
395 tim 1808 parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
396 gezelter 1490
397 tim 1739 for(i = 0; i < info_->getNGlobalMolecules(); i++) {
398 tim 1883 int which_node = info_->getMolToProc(i);
399 gezelter 1490
400 tim 1739 if (which_node == worldRank) {
401     //molecule with global index i belongs to this processor
402    
403     mol = info_->getMoleculeByGlobalIndex(i);
404     if (mol == NULL) {
405 tim 1884 sprintf(painCave.errMsg, "Molecule not found on node %d!", worldRank);
406 tim 1739 painCave.isFatal = 1;
407     simError();
408     }
409    
410     nCurObj = mol->getNIntegrableObjects();
411 gezelter 1490
412 tim 1739 MPI_Send(&nCurObj, 1, MPI_INT, masterNode, TAKE_THIS_TAG_INT,
413     MPI_COMM_WORLD);
414 gezelter 1490
415 tim 1739 for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
416     integrableObject = mol->nextIntegrableObject(ii)){
417    
418     MPI_Recv(read_buffer, maxBufferSize, MPI_CHAR, masterNode,
419     TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
420 gezelter 1490
421 tim 1739 parseDumpLine(read_buffer, integrableObject);
422     }
423    
424     }
425    
426     }
427    
428     }
429 gezelter 1490
430 tim 1739 #endif
431 gezelter 1490
432 tim 1739 }
433 gezelter 1490
434 tim 1739 void DumpReader::parseDumpLine(char *line, StuntDouble *integrableObject) {
435 gezelter 1490
436 tim 1739 Vector3d pos; // position place holders
437     Vector3d vel; // velocity placeholders
438     Quat4d q; // the quaternions
439     Vector3d ji; // angular velocity placeholders;
440     StringTokenizer tokenizer(line);
441     int nTokens;
442    
443     nTokens = tokenizer.countTokens();
444 gezelter 1490
445 tim 1808 if (nTokens < 14) {
446     sprintf(painCave.errMsg,
447     "Not enough Tokens.\n");
448     painCave.isFatal = 1;
449     simError();
450 gezelter 1490 }
451    
452 tim 1739 std::string name = tokenizer.nextToken();
453 gezelter 1490
454 tim 1739 if (name != integrableObject->getType()) {
455    
456 gezelter 1490 }
457    
458 tim 1739 pos[0] = tokenizer.nextTokenAsFloat();
459     pos[1] = tokenizer.nextTokenAsFloat();
460     pos[2] = tokenizer.nextTokenAsFloat();
461     integrableObject->setPos(pos);
462    
463     vel[0] = tokenizer.nextTokenAsFloat();
464     vel[1] = tokenizer.nextTokenAsFloat();
465     vel[2] = tokenizer.nextTokenAsFloat();
466     integrableObject->setVel(vel);
467 gezelter 1490
468 tim 1739 if (integrableObject->isDirectional()) {
469 tim 1808
470 tim 1739 q[0] = tokenizer.nextTokenAsFloat();
471     q[1] = tokenizer.nextTokenAsFloat();
472     q[2] = tokenizer.nextTokenAsFloat();
473     q[3] = tokenizer.nextTokenAsFloat();
474 gezelter 1490
475 tim 1739 double qlen = q.length();
476     if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
477    
478     sprintf(painCave.errMsg,
479     "initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2 ~ 0).\n");
480     painCave.isFatal = 1;
481     simError();
482    
483 tim 1864 }
484    
485     q.normalize();
486    
487 tim 1739 integrableObject->setQ(q);
488    
489     ji[0] = tokenizer.nextTokenAsFloat();
490     ji[1] = tokenizer.nextTokenAsFloat();
491     ji[2] = tokenizer.nextTokenAsFloat();
492     integrableObject->setJ(ji);
493 gezelter 1490 }
494    
495     }
496    
497    
498 tim 1808 void DumpReader::parseCommentLine(char* line, Snapshot* s) {
499 tim 1739 double currTime;
500     Mat3x3d hmat;
501     double chi;
502     double integralOfChiDt;
503     Mat3x3d eta;
504 gezelter 1490
505 tim 1739 StringTokenizer tokenizer(line);
506     int nTokens;
507 gezelter 1490
508 tim 1739 nTokens = tokenizer.countTokens();
509 gezelter 1490
510 tim 1739 //comment line should at least contain 10 tokens: current time(1 token) and h-matrix(9 tokens)
511     if (nTokens < 10) {
512     sprintf(painCave.errMsg,
513     "Not enough tokens in comment line: %s", line);
514     painCave.isFatal = 1;
515     simError();
516 gezelter 1490 }
517    
518 tim 1739 //read current time
519     currTime = tokenizer.nextTokenAsFloat();
520 tim 1808 s->setTime(currTime);
521 gezelter 1490
522 tim 1739 //read h-matrix
523     hmat(0, 0) = tokenizer.nextTokenAsFloat();
524     hmat(0, 1) = tokenizer.nextTokenAsFloat();
525     hmat(0, 2) = tokenizer.nextTokenAsFloat();
526     hmat(1, 0) = tokenizer.nextTokenAsFloat();
527     hmat(1, 1) = tokenizer.nextTokenAsFloat();
528     hmat(1, 2) = tokenizer.nextTokenAsFloat();
529     hmat(2, 0) = tokenizer.nextTokenAsFloat();
530     hmat(2, 1) = tokenizer.nextTokenAsFloat();
531     hmat(2, 2) = tokenizer.nextTokenAsFloat();
532     s->setHmat(hmat);
533 gezelter 1490
534 tim 1808 //read chi and integrablOfChidt, they should apprear in pair
535     if (tokenizer.countTokens() >= 2) {
536     chi = tokenizer.nextTokenAsFloat();
537     integralOfChiDt = tokenizer.nextTokenAsFloat();
538 gezelter 1490
539 tim 1808 s->setChi(chi);
540     s->setIntegralOfChiDt(integralOfChiDt);
541     }
542    
543     //read eta (eta is 3x3 matrix)
544     if (tokenizer.countTokens() >= 9) {
545     eta(0, 0) = tokenizer.nextTokenAsFloat();
546     eta(0, 1) = tokenizer.nextTokenAsFloat();
547     eta(0, 2) = tokenizer.nextTokenAsFloat();
548     eta(1, 0) = tokenizer.nextTokenAsFloat();
549     eta(1, 1) = tokenizer.nextTokenAsFloat();
550     eta(1, 2) = tokenizer.nextTokenAsFloat();
551     eta(2, 0) = tokenizer.nextTokenAsFloat();
552     eta(2, 1) = tokenizer.nextTokenAsFloat();
553     eta(2, 2) = tokenizer.nextTokenAsFloat();
554 gezelter 1490
555 tim 1808 s->setEta(eta);
556 tim 1739 }
557 gezelter 1490
558 tim 1739
559 gezelter 1490 }
560    
561 tim 1739 }//end namespace oopse