ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/io/DumpReader.cpp
(Generate patch)

Comparing trunk/src/io/DumpReader.cpp (file contents):
Revision 381 by tim, Tue Mar 1 14:45:45 2005 UTC vs.
Revision 1249 by gezelter, Mon Jun 2 02:42:50 2008 UTC

# Line 1 | Line 1
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. 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 < */
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. 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 > #define _LARGEFILE_SOURCE64
43 > #define _FILE_OFFSET_BITS 64
44  
45 < #define _LARGEFILE_SOURCE64
46 < #define _FILE_OFFSET_BITS 64
47 <
48 < #include <sys/types.h>
49 < #include <sys/stat.h>
50 <
51 < #include <iostream>
52 < #include <math.h>
53 <
54 < #include <stdio.h>
55 < #include <stdlib.h>
56 < #include <string.h>
57 <
58 < #include "io/DumpReader.hpp"
59 < #include "primitives/Molecule.hpp"
60 < #include "utils/simError.h"
61 < #include "utils/MemoryUtils.hpp"
62 < #include "utils/StringTokenizer.hpp"
63 <
64 < #ifdef IS_MPI
65 <
66 < #include <mpi.h>
67 < #define TAKE_THIS_TAG_CHAR 0
68 < #define TAKE_THIS_TAG_INT 1
69 <
70 < #endif // is_mpi
71 <
72 <
73 < namespace oopse {
74 <
75 < DumpReader::DumpReader(SimInfo* info, const std::string& filename)
76 <                     : info_(info), filename_(filename), isScanned_(false), nframes_(0) {
77 <
78 < #ifdef IS_MPI
79 <
80 <    if (worldRank == 0) {
81 < #endif
82 <
83 <        inFile_ = fopen(filename_.c_str(), "r");
84 <
85 <        if (inFile_ == NULL) {
86 <            sprintf(painCave.errMsg, "DumpReader: Cannot open file: %s\n", filename_.c_str());
87 <            painCave.isFatal = 1;
88 <            simError();
45 > #include <sys/types.h>
46 > #include <sys/stat.h>
47 >
48 > #include <iostream>
49 > #include <math.h>
50 >
51 > #include <stdio.h>
52 > #include <stdlib.h>
53 > #include <string.h>
54 >
55 > #include "io/DumpReader.hpp"
56 > #include "primitives/Molecule.hpp"
57 > #include "utils/simError.h"
58 > #include "utils/MemoryUtils.hpp"
59 > #include "utils/StringTokenizer.hpp"
60 >
61 > #ifdef IS_MPI
62 >
63 > #include <mpi.h>
64 > #define TAKE_THIS_TAG_CHAR 0
65 > #define TAKE_THIS_TAG_INT 1
66 >
67 > #endif // is_mpi
68 >
69 >
70 > namespace oopse {
71 >  
72 >  DumpReader::DumpReader(SimInfo* info, const std::string& filename)
73 >    : info_(info), filename_(filename), isScanned_(false), nframes_(0), needCOMprops_(false) {
74 >    
75 > #ifdef IS_MPI
76 >    
77 >    if (worldRank == 0) {
78 > #endif
79 >      
80 >      inFile_ = new std::ifstream(filename_.c_str());
81 >      
82 >      if (inFile_->fail()) {
83 >        sprintf(painCave.errMsg,
84 >                "DumpReader: Cannot open file: %s\n",
85 >                filename_.c_str());
86 >        painCave.isFatal = 1;
87 >        simError();
88 >      }
89 >      
90 > #ifdef IS_MPI
91 >      
92 >    }
93 >    
94 >    strcpy(checkPointMsg, "Dump file opened for reading successfully.");
95 >    errorCheckPoint();
96 >    
97 > #endif
98 >    
99 >    return;
100 >  }
101 >  
102 >  DumpReader::~DumpReader() {
103 >    
104 > #ifdef IS_MPI
105 >    
106 >    if (worldRank == 0) {
107 > #endif
108 >      
109 >      delete inFile_;
110 >      
111 > #ifdef IS_MPI
112 >      
113 >    }
114 >    
115 >    strcpy(checkPointMsg, "Dump file closed successfully.");
116 >    errorCheckPoint();
117 >    
118 > #endif
119 >    
120 >    return;
121 >  }
122 >  
123 >  int DumpReader::getNFrames(void) {
124 >    
125 >    if (!isScanned_)
126 >      scanFile();
127 >    
128 >    return nframes_;
129 >  }
130 >  
131 >  void DumpReader::scanFile(void) {
132 >    int lineNo = 0;
133 >    std::streampos prevPos;
134 >    std::streampos  currPos;
135 >    
136 > #ifdef IS_MPI
137 >    
138 >    if (worldRank == 0) {
139 > #endif // is_mpi
140 >      
141 >      currPos = inFile_->tellg();
142 >      prevPos = currPos;
143 >      bool foundOpenSnapshotTag = false;
144 >      bool foundClosedSnapshotTag = false;
145 >      while(inFile_->getline(buffer, bufferSize)) {
146 >        ++lineNo;
147 >        
148 >        std::string line = buffer;
149 >        currPos = inFile_->tellg();
150 >        if (line.find("<Snapshot>")!= std::string::npos) {
151 >          if (foundOpenSnapshotTag) {
152 >            sprintf(painCave.errMsg,
153 >                    "DumpReader:<Snapshot> is multiply nested at line %d in %s \n", lineNo,
154 >                    filename_.c_str());
155 >            painCave.isFatal = 1;
156 >            simError();          
157 >          }
158 >          foundOpenSnapshotTag = true;
159 >          foundClosedSnapshotTag = false;
160 >          framePos_.push_back(prevPos);
161 >          
162 >        } else if (line.find("</Snapshot>") != std::string::npos){
163 >          if (!foundOpenSnapshotTag) {
164 >            sprintf(painCave.errMsg,
165 >                    "DumpReader:</Snapshot> appears before <Snapshot> at line %d in %s \n", lineNo,
166 >                    filename_.c_str());
167 >            painCave.isFatal = 1;
168 >            simError();
169 >          }
170 >          
171 >          if (foundClosedSnapshotTag) {
172 >            sprintf(painCave.errMsg,
173 >                    "DumpReader:</Snapshot> appears multiply nested at line %d in %s \n", lineNo,
174 >                    filename_.c_str());
175 >            painCave.isFatal = 1;
176 >            simError();
177 >          }
178 >          foundClosedSnapshotTag = true;
179 >          foundOpenSnapshotTag = false;
180          }
181 <
182 < #ifdef IS_MPI
183 <
184 <    }
185 <
186 <    strcpy(checkPointMsg, "Dump file opened for reading successfully.");
187 <    MPIcheckPoint();
188 <
189 < #endif
190 <
191 <    return;
192 < }
193 <
194 < DumpReader::~DumpReader() {
195 <
196 < #ifdef IS_MPI
197 <
198 <    if (worldRank == 0) {
199 < #endif
200 <
201 <        int error;
202 <        error = fclose(inFile_);
181 >        prevPos = currPos;
182 >      }
183 >      
184 >      // only found <Snapshot> for the last frame means the file is corrupted, we should discard
185 >      // it and give a warning message
186 >      if (foundOpenSnapshotTag) {
187 >        sprintf(painCave.errMsg,
188 >                "DumpReader: last frame in %s is invalid\n", filename_.c_str());
189 >        painCave.isFatal = 0;
190 >        simError();      
191 >        framePos_.pop_back();
192 >      }
193 >      
194 >      nframes_ = framePos_.size();
195 >      
196 >      if (nframes_ == 0) {
197 >        sprintf(painCave.errMsg,
198 >                "DumpReader: %s does not contain a valid frame\n", filename_.c_str());
199 >        painCave.isFatal = 1;
200 >        simError();      
201 >      }
202 > #ifdef IS_MPI
203 >    }
204 >    
205 >    MPI_Bcast(&nframes_, 1, MPI_INT, 0, MPI_COMM_WORLD);
206 >    
207 > #endif // is_mpi
208 >    
209 >    isScanned_ = true;
210 >  }
211 >  
212 >  void DumpReader::readFrame(int whichFrame) {
213 >    if (!isScanned_)
214 >      scanFile();
215 >        
216 >    int storageLayout = info_->getSnapshotManager()->getStorageLayout();
217 >    
218 >    if (storageLayout & DataStorage::dslPosition) {
219 >      needPos_ = true;
220 >    } else {
221 >      needPos_ = false;
222 >    }
223 >    
224 >    if (storageLayout & DataStorage::dslVelocity) {
225 >      needVel_ = true;
226 >    } else {
227 >      needVel_ = false;
228 >    }
229 >    
230 >    if (storageLayout & DataStorage::dslAmat || storageLayout & DataStorage::dslElectroFrame) {
231 >      needQuaternion_ = true;
232 >    } else {
233 >      needQuaternion_ = false;
234 >    }
235 >    
236 >    if (storageLayout & DataStorage::dslAngularMomentum) {
237 >      needAngMom_ = true;
238 >    } else {
239 >      needAngMom_ = false;    
240 >    }
241 >    
242 >    readSet(whichFrame);
243  
244 <        if (error) {
245 <            sprintf(painCave.errMsg, "DumpReader Error: Error closing %s\n", filename_.c_str());
246 <            painCave.isFatal = 1;            
247 <            simError();
244 >    if (needCOMprops_) {
245 >      Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
246 >      Vector3d com;
247 >      Vector3d comvel;
248 >      Vector3d comw;
249 >      if (needPos_ && needVel_){
250 >        info_->getComAll(com, comvel);
251 >        comw = info_->getAngularMomentum();
252 >      }else{
253 >        com = info_->getCom();
254 >        comvel = 0.0;
255 >        comw   = 0.0;
256 >      }
257 >      s->setCOMprops(com, comvel, comw);      
258 >    }
259 >
260 >  }
261 >  
262 >  void DumpReader::readSet(int whichFrame) {    
263 >    std::string line;
264 >
265 > #ifndef IS_MPI
266 >    inFile_->clear();  
267 >    inFile_->seekg(framePos_[whichFrame]);
268 >
269 >    std::istream& inputStream = *inFile_;    
270 >
271 > #else
272 >    int masterNode = 0;
273 >    std::stringstream sstream;
274 >    if (worldRank == masterNode) {
275 >      std::string sendBuffer;
276 >
277 >      inFile_->clear();  
278 >      inFile_->seekg(framePos_[whichFrame]);
279 >      
280 >      while (inFile_->getline(buffer, bufferSize)) {
281 >
282 >        line = buffer;
283 >        sendBuffer += line;
284 >        sendBuffer += '\n';
285 >        if (line.find("</Snapshot>") != std::string::npos) {
286 >          break;
287 >        }        
288 >      }
289 >
290 >      int sendBufferSize = sendBuffer.size();
291 >      MPI_Bcast(&sendBufferSize, 1, MPI_INT, masterNode, MPI_COMM_WORLD);    
292 >      MPI_Bcast((void *)sendBuffer.c_str(), sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);    
293 >      
294 >      sstream.str(sendBuffer);
295 >    } else {
296 >      int sendBufferSize;
297 >      MPI_Bcast(&sendBufferSize, 1, MPI_INT, masterNode, MPI_COMM_WORLD);    
298 >      char * recvBuffer = new char[sendBufferSize+1];
299 >      MPI_Bcast(recvBuffer, sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);    
300 >      sstream.str(recvBuffer);
301 >    }      
302 >
303 >    std::istream& inputStream = sstream;  
304 > #endif
305 >
306 >    inputStream.getline(buffer, bufferSize);
307 >
308 >    line = buffer;
309 >    if (line.find("<Snapshot>") == std::string::npos) {
310 >      sprintf(painCave.errMsg,
311 >              "DumpReader Error: can not find <Snapshot>\n");
312 >      painCave.isFatal = 1;
313 >      simError();
314 >    }
315 >    
316 >    //read frameData
317 >    readFrameProperties(inputStream);
318 >
319 >    //read StuntDoubles
320 >    readStuntDoubles(inputStream);    
321 >
322 >    inputStream.getline(buffer, bufferSize);
323 >    line = buffer;
324 >    if (line.find("</Snapshot>") == std::string::npos) {
325 >      sprintf(painCave.errMsg,
326 >              "DumpReader Error: can not find </Snapshot>\n");
327 >      painCave.isFatal = 1;
328 >      simError();
329 >    }        
330 >  
331 >  }
332 >  
333 >  void DumpReader::parseDumpLine(const std::string& line) {
334 >
335 >      
336 >    StringTokenizer tokenizer(line);
337 >    int nTokens;
338 >    
339 >    nTokens = tokenizer.countTokens();
340 >    
341 >    if (nTokens < 2) {
342 >      sprintf(painCave.errMsg,
343 >              "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
344 >      painCave.isFatal = 1;
345 >      simError();
346 >    }
347 >
348 >    int index = tokenizer.nextTokenAsInt();
349 >
350 >    StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index);
351 >
352 >    if (integrableObject == NULL) {
353 >      return;
354 >    }
355 >    std::string type = tokenizer.nextToken();
356 >    int size = type.size();
357 >
358 >    for(int i = 0; i < size; ++i) {
359 >      switch(type[i]) {
360 >        
361 >        case 'p': {
362 >            Vector3d pos;
363 >            pos[0] = tokenizer.nextTokenAsDouble();
364 >            pos[1] = tokenizer.nextTokenAsDouble();
365 >            pos[2] = tokenizer.nextTokenAsDouble();
366 >            if (needPos_) {
367 >              integrableObject->setPos(pos);
368 >            }            
369 >            break;
370          }
371 +        case 'v' : {
372 +            Vector3d vel;
373 +            vel[0] = tokenizer.nextTokenAsDouble();
374 +            vel[1] = tokenizer.nextTokenAsDouble();
375 +            vel[2] = tokenizer.nextTokenAsDouble();
376 +            if (needVel_) {
377 +              integrableObject->setVel(vel);
378 +            }
379 +            break;
380 +        }
381  
382 <        MemoryUtils::deleteVectorOfPointer(framePos_);
382 >        case 'q' : {
383 >           Quat4d q;
384 >           if (integrableObject->isDirectional()) {
385 >              
386 >             q[0] = tokenizer.nextTokenAsDouble();
387 >             q[1] = tokenizer.nextTokenAsDouble();
388 >             q[2] = tokenizer.nextTokenAsDouble();
389 >             q[3] = tokenizer.nextTokenAsDouble();
390 >              
391 >             RealType qlen = q.length();
392 >             if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
393 >                
394 >               sprintf(painCave.errMsg,
395 >                       "DumpReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n");
396 >               painCave.isFatal = 1;
397 >               simError();
398 >                
399 >             }  
400 >              
401 >             q.normalize();
402 >             if (needQuaternion_) {            
403 >               integrableObject->setQ(q);
404 >             }              
405 >           }            
406 >           break;
407 >        }  
408 >        case 'j' : {
409 >          Vector3d ji;
410 >          if (integrableObject->isDirectional()) {
411 >             ji[0] = tokenizer.nextTokenAsDouble();
412 >             ji[1] = tokenizer.nextTokenAsDouble();
413 >             ji[2] = tokenizer.nextTokenAsDouble();
414 >             if (needAngMom_) {
415 >               integrableObject->setJ(ji);
416 >             }
417 >          }
418 >          break;
419 >        }  
420 >        case 'f': {
421  
422 < #ifdef IS_MPI
422 >          Vector3d force;
423 >          force[0] = tokenizer.nextTokenAsDouble();
424 >          force[1] = tokenizer.nextTokenAsDouble();
425 >          force[2] = tokenizer.nextTokenAsDouble();          
426 >          integrableObject->setFrc(force);
427 >          break;
428 >        }
429 >        case 't' : {
430  
431 +           Vector3d torque;
432 +           torque[0] = tokenizer.nextTokenAsDouble();
433 +           torque[1] = tokenizer.nextTokenAsDouble();
434 +           torque[2] = tokenizer.nextTokenAsDouble();          
435 +           integrableObject->setTrq(torque);          
436 +           break;
437 +        }
438 +        default: {
439 +               sprintf(painCave.errMsg,
440 +                       "DumpReader Error: %s is an unrecognized type\n", type.c_str());
441 +               painCave.isFatal = 1;
442 +               simError();
443 +          break;  
444 +        }
445 +
446 +      }
447      }
448 +    
449 +  }
450 +  
451  
452 <    strcpy(checkPointMsg, "Dump file closed successfully.");
123 <    MPIcheckPoint();
452 >  void  DumpReader::readStuntDoubles(std::istream& inputStream) {
453  
454 < #endif
454 >    inputStream.getline(buffer, bufferSize);
455 >    std::string line(buffer);
456 >    
457 >    if (line.find("<StuntDoubles>") == std::string::npos) {
458 >      sprintf(painCave.errMsg,
459 >              "DumpReader Error: Missing <StuntDoubles>\n");
460 >      painCave.isFatal = 1;
461 >      simError();
462 >    }
463  
464 <    return;
465 < }
464 >    while(inputStream.getline(buffer, bufferSize)) {
465 >      line = buffer;
466 >      
467 >      if(line.find("</StuntDoubles>") != std::string::npos) {
468 >        break;
469 >      }
470  
471 < int DumpReader::getNFrames(void) {
471 >      parseDumpLine(line);
472 >    }
473 >  
474 >  }
475  
476 <    if (!isScanned_)
133 <        scanFile();
476 >  void DumpReader::readFrameProperties(std::istream& inputStream) {
477  
478 <    return nframes_;
479 < }
478 >    Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
479 >    inputStream.getline(buffer, bufferSize);
480 >    std::string line(buffer);
481 >
482 >    if (line.find("<FrameData>") == std::string::npos) {
483 >      sprintf(painCave.errMsg,
484 >              "DumpReader Error: Missing <FrameData>\n");
485 >      painCave.isFatal = 1;
486 >      simError();
487 >    }
488 >
489 >    while(inputStream.getline(buffer, bufferSize)) {
490 >      line = buffer;
491 >      
492 >      if(line.find("</FrameData>") != std::string::npos) {
493 >        break;
494 >      }
495 >      
496 >      StringTokenizer tokenizer(line, " ;\t\n\r{}:,");
497 >      if (!tokenizer.hasMoreTokens()) {
498 >        sprintf(painCave.errMsg,
499 >                "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
500 >        painCave.isFatal = 1;
501 >        simError();      
502 >      }
503  
504 < void DumpReader::scanFile(void) {
505 <  int i, j;
506 <  int lineNum = 0;
507 <  char readBuffer[maxBufferSize];
508 <  fpos_t * currPos;
509 <
510 < #ifdef IS_MPI
511 <
512 <    if (worldRank == 0) {
513 < #endif // is_mpi
514 <
515 <        rewind(inFile_);
516 <
517 <        currPos = new fpos_t;
518 <        fgetpos(inFile_, currPos);
519 <        fgets(readBuffer, sizeof(readBuffer), inFile_);
520 <        lineNum++;
521 <
522 <        if (feof(inFile_)) {
523 <            sprintf(painCave.errMsg,
524 <                    "DumpReader Error: File \"%s\" ended unexpectedly at line %d\n",
525 <                    filename_.c_str(),
526 <                    lineNum);
527 <            painCave.isFatal = 1;
528 <            simError();
529 <        }
530 <
531 <        while (!feof(inFile_)) {
532 <            framePos_.push_back(currPos);
533 <
534 <            i = atoi(readBuffer);
535 <
536 <            fgets(readBuffer, sizeof(readBuffer), inFile_);
537 <            lineNum++;
538 <
539 <            if (feof(inFile_)) {
540 <                sprintf(painCave.errMsg,
541 <                        "DumpReader Error: File \"%s\" ended unexpectedly at line %d\n",
542 <                        filename_.c_str(),
543 <                        lineNum);
178 <                painCave.isFatal = 1;
179 <                simError();
180 <            }
181 <
182 <            for(j = 0; j < i; j++) {
183 <                fgets(readBuffer, sizeof(readBuffer), inFile_);
184 <                lineNum++;
185 <
186 <                if (feof(inFile_)) {
187 <                    sprintf(painCave.errMsg,
188 <                            "DumpReader Error: File \"%s\" ended unexpectedly at line %d,"
189 <                                " with atom %d\n", filename_.c_str(),
190 <                            lineNum,
191 <                            j);
192 <
193 <                    painCave.isFatal = 1;
194 <                    simError();
195 <                }
196 <            }
197 <
198 <            currPos = new fpos_t;
199 <            fgetpos(inFile_, currPos);
200 <            fgets(readBuffer, sizeof(readBuffer), inFile_);
201 <            lineNum++;
202 <        }
203 <
204 <        delete currPos;
205 <        rewind(inFile_);
206 <        
207 <        nframes_ = framePos_.size();
208 < #ifdef IS_MPI
504 >      std::string propertyName = tokenizer.nextToken();
505 >      if (propertyName == "Time") {
506 >        RealType currTime = tokenizer.nextTokenAsDouble();
507 >        s->setTime(currTime);
508 >      } else if (propertyName == "Hmat"){
509 >        Mat3x3d hmat;
510 >        hmat(0, 0) = tokenizer.nextTokenAsDouble();
511 >        hmat(0, 1) = tokenizer.nextTokenAsDouble();
512 >        hmat(0, 2) = tokenizer.nextTokenAsDouble();
513 >        hmat(1, 0) = tokenizer.nextTokenAsDouble();
514 >        hmat(1, 1) = tokenizer.nextTokenAsDouble();
515 >        hmat(1, 2) = tokenizer.nextTokenAsDouble();
516 >        hmat(2, 0) = tokenizer.nextTokenAsDouble();
517 >        hmat(2, 1) = tokenizer.nextTokenAsDouble();
518 >        hmat(2, 2) = tokenizer.nextTokenAsDouble();
519 >        s->setHmat(hmat);      
520 >      } else if (propertyName == "Thermostat") {
521 >        RealType chi = tokenizer.nextTokenAsDouble();
522 >        RealType integralOfChiDt = tokenizer.nextTokenAsDouble();
523 >        s->setChi(chi);
524 >        s->setIntegralOfChiDt(integralOfChiDt);        
525 >     } else if (propertyName == "Barostat") {
526 >        Mat3x3d eta;
527 >        eta(0, 0) = tokenizer.nextTokenAsDouble();
528 >        eta(0, 1) = tokenizer.nextTokenAsDouble();
529 >        eta(0, 2) = tokenizer.nextTokenAsDouble();
530 >        eta(1, 0) = tokenizer.nextTokenAsDouble();
531 >        eta(1, 1) = tokenizer.nextTokenAsDouble();
532 >        eta(1, 2) = tokenizer.nextTokenAsDouble();
533 >        eta(2, 0) = tokenizer.nextTokenAsDouble();
534 >        eta(2, 1) = tokenizer.nextTokenAsDouble();
535 >        eta(2, 2) = tokenizer.nextTokenAsDouble();
536 >        s->setEta(eta);
537 >      } else {
538 >        sprintf(painCave.errMsg,
539 >                "DumpReader Error: %s is an invalid property in <FrameData>\n", propertyName.c_str());
540 >        painCave.isFatal = 0;
541 >        simError();        
542 >      }
543 >      
544      }
545  
546 <    MPI_Bcast(&nframes_, 1, MPI_INT, 0, MPI_COMM_WORLD);
212 <
213 <    strcpy(checkPointMsg, "Successfully scanned DumpFile\n");
214 <    MPIcheckPoint();
215 <
216 < #endif // is_mpi
217 <
218 <    isScanned_ = true;
219 < }
220 <
221 < void DumpReader::readFrame(int whichFrame) {
222 <    int storageLayout = info_->getSnapshotManager()->getStorageLayout();
223 <
224 <    if (storageLayout & DataStorage::dslPosition) {
225 <        needPos_ = true;
226 <    } else {
227 <        needPos_ = false;
228 <    }
229 <
230 <    if (storageLayout & DataStorage::dslVelocity) {
231 <        needVel_ = true;
232 <    } else {
233 <        needVel_ = false;
234 <    }
235 <
236 <    if (storageLayout & DataStorage::dslAmat || storageLayout & DataStorage::dslElectroFrame) {
237 <        needQuaternion_ = true;
238 <    } else {
239 <        needQuaternion_ = false;
240 <    }
241 <
242 <    if (storageLayout & DataStorage::dslAngularMomentum) {
243 <        needAngMom_ = true;
244 <    } else {
245 <        needAngMom_ = false;    
246 <    }
247 <
248 <    readSet(whichFrame);
249 < }
250 <
251 < void DumpReader::readSet(int whichFrame) {
252 <  int i;
253 <  int nTotObjs;                  // the number of atoms
254 <  char read_buffer[maxBufferSize];  //the line buffer for reading
255 <  char * eof_test;               // ptr to see when we reach the end of the file
256 <
257 <  Molecule* mol;
258 <  StuntDouble* integrableObject;
259 <  SimInfo::MoleculeIterator mi;
260 <  Molecule::IntegrableObjectIterator ii;
261 <
262 < #ifndef IS_MPI
263 <
264 <    fsetpos(inFile_, framePos_[whichFrame]);
265 <    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
266 <
267 <    if (eof_test == NULL) {
268 <        sprintf(painCave.errMsg,
269 <                "DumpReader error: error reading 1st line of \"%s\"\n",
270 <                filename_.c_str());
271 <        painCave.isFatal = 1;
272 <        simError();
273 <    }
274 <
275 <    nTotObjs = atoi(read_buffer);
276 <
277 <    if (nTotObjs != info_->getNGlobalIntegrableObjects()) {
278 <        sprintf(painCave.errMsg,
279 <                "DumpReader error. %s nIntegrable, %d, "
280 <                    "does not match the meta-data file's nIntegrable, %d.\n",
281 <                filename_.c_str(),
282 <                nTotObjs,
283 <                info_->getNGlobalIntegrableObjects());
284 <
285 <        painCave.isFatal = 1;
286 <        simError();
287 <    }
288 <
289 <    //read the box mat from the comment line
290 <
291 <    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
292 <
293 <    if (eof_test == NULL) {
294 <        sprintf(painCave.errMsg, "DumpReader Error: error in reading commment in %s\n",
295 <                filename_.c_str());
296 <        painCave.isFatal = 1;
297 <        simError();
298 <    }
299 <
300 <    parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
301 <
302 <    //parse dump lines
303 <
304 <    for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
305 <
306 <        for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
307 <            integrableObject = mol->nextIntegrableObject(ii)) {          
308 <
309 <            eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
310 <
311 <            if (eof_test == NULL) {
312 <                sprintf(painCave.errMsg,
313 <                        "DumpReader Error: error in reading file %s\n"
314 <                            "natoms  = %d; index = %d\n"
315 <                            "error reading the line from the file.\n",
316 <                        filename_.c_str(),
317 <                        nTotObjs,
318 <                        i);
319 <
320 <                painCave.isFatal = 1;
321 <                simError();
322 <            }
323 <
324 <            parseDumpLine(read_buffer, integrableObject);
325 <            
326 <            }
327 <    }
328 <
329 <    // MPI Section of code..........
330 <
331 < #else //IS_MPI
332 <
333 <    // first thing first, suspend fatalities.
334 <    int masterNode = 0;
335 <    int nCurObj;
336 <    painCave.isEventLoop = 1;
337 <
338 <    int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
339 <    int haveError;
340 <
341 <    MPI_Status istatus;
342 <    int nitems;
343 <
344 <    nTotObjs = info_->getNGlobalIntegrableObjects();
345 <    haveError = 0;
346 <
347 <    if (worldRank == masterNode) {
348 <        fsetpos(inFile_, framePos_[whichFrame]);
349 <
350 <        eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
351 <
352 <        if (eof_test == NULL) {
353 <            sprintf(painCave.errMsg, "DumpReader Error: Error reading 1st line of %s \n ",
354 <                    filename_.c_str());
355 <            painCave.isFatal = 1;
356 <            simError();
357 <        }
358 <
359 <        nitems = atoi(read_buffer);
360 <
361 <        // Check to see that the number of integrable objects in the
362 <        // intial configuration file is the same as derived from the
363 <        // meta-data file.
364 <
365 <        if (nTotObjs != nitems) {
366 <            sprintf(painCave.errMsg,
367 <                    "DumpReader Error. %s nIntegrable, %d, "
368 <                        "does not match the meta-data file's nIntegrable, %d.\n",
369 <                    filename_.c_str(),
370 <                    nTotObjs,
371 <                    info_->getNGlobalIntegrableObjects());
546 >  }
547  
548 <            painCave.isFatal = 1;
549 <            simError();
375 <        }
376 <
377 <        //read the boxMat from the comment line
378 <
379 <        eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
380 <
381 <        if (eof_test == NULL) {
382 <            sprintf(painCave.errMsg, "DumpReader Error: error in reading commment in %s\n",
383 <                    filename_.c_str());
384 <            painCave.isFatal = 1;
385 <            simError();
386 <        }
387 <
388 <        //Every single processor will parse the comment line by itself
389 <        //By using this way, we might lose some efficiency, but if we want to add
390 <        //more parameters into comment line, we only need to modify function
391 <        //parseCommentLine
392 <
393 <        MPI_Bcast(read_buffer, maxBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
394 <        parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
395 <
396 <        for(i = 0; i < info_->getNGlobalMolecules(); i++) {
397 <            int which_node = info_->getMolToProc(i);
398 <
399 <            if (which_node == masterNode) {
400 <                //molecules belong to master node
401 <
402 <                mol = info_->getMoleculeByGlobalIndex(i);
403 <
404 <                if (mol == NULL) {
405 <                    sprintf(painCave.errMsg, "DumpReader Error: Molecule not found on node %d!", worldRank);
406 <                        painCave.isFatal = 1;
407 <                    simError();
408 <                }
409 <
410 <                for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
411 <                       integrableObject = mol->nextIntegrableObject(ii)){
412 <                        
413 <                    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
414 <
415 <                    if (eof_test == NULL) {
416 <                        sprintf(painCave.errMsg,
417 <                                "DumpReader Error: error in reading file %s\n"
418 <                                    "natoms  = %d; index = %d\n"
419 <                                    "error reading the line from the file.\n",
420 <                                filename_.c_str(),
421 <                                nTotObjs,
422 <                                i);
423 <
424 <                        painCave.isFatal = 1;
425 <                        simError();
426 <                    }
427 <
428 <                    parseDumpLine(read_buffer, integrableObject);
429 <                }
430 <            } else {
431 <                //molecule belongs to slave nodes
432 <
433 <                MPI_Recv(&nCurObj, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT,
434 <                         MPI_COMM_WORLD, &istatus);
435 <
436 <                for(int j = 0; j < nCurObj; j++) {
437 <                    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
438 <
439 <                    if (eof_test == NULL) {
440 <                        sprintf(painCave.errMsg,
441 <                                "DumpReader Error: error in reading file %s\n"
442 <                                    "natoms  = %d; index = %d\n"
443 <                                    "error reading the line from the file.\n",
444 <                                filename_.c_str(),
445 <                                nTotObjs,
446 <                                i);
447 <
448 <                        painCave.isFatal = 1;
449 <                        simError();
450 <                    }
451 <                    
452 <                    MPI_Send(read_buffer, maxBufferSize, MPI_CHAR, which_node,
453 <                             TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
454 <                }
455 <            }
456 <        }
457 <    } else {
458 <        //actions taken at slave nodes
459 <        MPI_Bcast(read_buffer, maxBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
460 <
461 <        /**@todo*/
462 <        parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
463 <
464 <        for(i = 0; i < info_->getNGlobalMolecules(); i++) {
465 <            int which_node = info_->getMolToProc(i);
466 <
467 <            if (which_node == worldRank) {
468 <                //molecule with global index i belongs to this processor
469 <                
470 <                mol = info_->getMoleculeByGlobalIndex(i);
471 <                if (mol == NULL) {
472 <                    sprintf(painCave.errMsg, "DumpReader Error: Molecule not found on node %d!", worldRank);
473 <                    painCave.isFatal = 1;
474 <                    simError();
475 <                }
476 <                
477 <                nCurObj = mol->getNIntegrableObjects();
478 <
479 <                MPI_Send(&nCurObj, 1, MPI_INT, masterNode, TAKE_THIS_TAG_INT,
480 <                         MPI_COMM_WORLD);
481 <
482 <                for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
483 <                       integrableObject = mol->nextIntegrableObject(ii)){
484 <                        
485 <                    MPI_Recv(read_buffer, maxBufferSize, MPI_CHAR, masterNode,
486 <                             TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
487 <
488 <                    parseDumpLine(read_buffer, integrableObject);
489 <                }
490 <                      
491 <            }
492 <            
493 <        }
494 <        
495 <    }
496 <
497 < #endif
498 <
499 < }
500 <
501 < void DumpReader::parseDumpLine(char *line, StuntDouble *integrableObject) {
502 <
503 <    Vector3d pos;  // position place holders
504 <    Vector3d vel;  // velocity placeholders
505 <    Quat4d q;    // the quaternions
506 <    Vector3d ji;   // angular velocity placeholders;
507 <    StringTokenizer tokenizer(line);
508 <    int nTokens;
509 <    
510 <    nTokens = tokenizer.countTokens();
511 <
512 <    if (nTokens < 14) {
513 <            sprintf(painCave.errMsg,
514 <                    "DumpReader Error: Not enough Tokens.\n%s\n", line);
515 <            painCave.isFatal = 1;
516 <            simError();
517 <    }
518 <
519 <    std::string name = tokenizer.nextToken();
520 <
521 <    if (name != integrableObject->getType()) {
522 <
523 <            sprintf(painCave.errMsg,
524 <                    "DumpReader Error: Atom type [%s] in %s does not match Atom Type [%s] in .md file.\n",
525 <                    name.c_str(), filename_.c_str(), integrableObject->getType().c_str());
526 <            painCave.isFatal = 1;
527 <            simError();        
528 <    }
529 <
530 <    pos[0] = tokenizer.nextTokenAsDouble();
531 <    pos[1] = tokenizer.nextTokenAsDouble();
532 <    pos[2] = tokenizer.nextTokenAsDouble();
533 <    if (needPos_) {
534 <        integrableObject->setPos(pos);
535 <    }
536 <    
537 <    vel[0] = tokenizer.nextTokenAsDouble();
538 <    vel[1] = tokenizer.nextTokenAsDouble();
539 <    vel[2] = tokenizer.nextTokenAsDouble();
540 <    if (needVel_) {
541 <        integrableObject->setVel(vel);
542 <    }
543 <    
544 <    if (integrableObject->isDirectional()) {
545 <        
546 <        q[0] = tokenizer.nextTokenAsDouble();
547 <        q[1] = tokenizer.nextTokenAsDouble();
548 <        q[2] = tokenizer.nextTokenAsDouble();
549 <        q[3] = tokenizer.nextTokenAsDouble();
550 <
551 <        double qlen = q.length();
552 <        if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
553 <            
554 <            sprintf(painCave.errMsg,
555 <                    "DumpReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2 ~ 0).\n");
556 <            painCave.isFatal = 1;
557 <            simError();
558 <            
559 <        }
560 <
561 <        q.normalize();
562 <        if (needQuaternion_) {          
563 <            integrableObject->setQ(q);
564 <        }
565 <
566 <        ji[0] = tokenizer.nextTokenAsDouble();
567 <        ji[1] = tokenizer.nextTokenAsDouble();
568 <        ji[2] = tokenizer.nextTokenAsDouble();
569 <        if (needAngMom_) {
570 <            integrableObject->setJ(ji);
571 <        }
572 <    }
573 <
574 < }
575 <
576 <
577 < void DumpReader::parseCommentLine(char* line, Snapshot* s) {
578 <    double currTime;
579 <    Mat3x3d hmat;
580 <    double chi;
581 <    double integralOfChiDt;
582 <    Mat3x3d eta;
583 <
584 <    StringTokenizer tokenizer(line);
585 <    int nTokens;
586 <
587 <    nTokens = tokenizer.countTokens();
588 <
589 <    //comment line should at least contain 10 tokens: current time(1 token) and  h-matrix(9 tokens)
590 <    if (nTokens < 10) {
591 <            sprintf(painCave.errMsg,
592 <                    "DumpReader Error: Not enough tokens in comment line: %s", line);
593 <            painCave.isFatal = 1;
594 <            simError();  
595 <    }
596 <
597 <    //read current time
598 <    currTime = tokenizer.nextTokenAsDouble();
599 <    s->setTime(currTime);
600 <    
601 <    //read h-matrix
602 <    hmat(0, 0) = tokenizer.nextTokenAsDouble();
603 <    hmat(0, 1) = tokenizer.nextTokenAsDouble();
604 <    hmat(0, 2) = tokenizer.nextTokenAsDouble();
605 <    hmat(1, 0) = tokenizer.nextTokenAsDouble();
606 <    hmat(1, 1) = tokenizer.nextTokenAsDouble();
607 <    hmat(1, 2) = tokenizer.nextTokenAsDouble();
608 <    hmat(2, 0) = tokenizer.nextTokenAsDouble();
609 <    hmat(2, 1) = tokenizer.nextTokenAsDouble();
610 <    hmat(2, 2) = tokenizer.nextTokenAsDouble();
611 <    s->setHmat(hmat);
612 <    
613 <    //read chi and integralOfChidt, they should apprear in pair
614 <    if (tokenizer.countTokens() >= 2) {
615 <        chi = tokenizer.nextTokenAsDouble();
616 <        integralOfChiDt = tokenizer.nextTokenAsDouble();            
617 <
618 <        s->setChi(chi);
619 <        s->setIntegralOfChiDt(integralOfChiDt);
620 <    }
621 <    
622 <    //read eta (eta is 3x3 matrix)
623 <    if (tokenizer.countTokens() >= 9) {
624 <        eta(0, 0) = tokenizer.nextTokenAsDouble();
625 <        eta(0, 1) = tokenizer.nextTokenAsDouble();
626 <        eta(0, 2) = tokenizer.nextTokenAsDouble();
627 <        eta(1, 0) = tokenizer.nextTokenAsDouble();
628 <        eta(1, 1) = tokenizer.nextTokenAsDouble();
629 <        eta(1, 2) = tokenizer.nextTokenAsDouble();
630 <        eta(2, 0) = tokenizer.nextTokenAsDouble();
631 <        eta(2, 1) = tokenizer.nextTokenAsDouble();
632 <        eta(2, 2) = tokenizer.nextTokenAsDouble();      
633 <
634 <        s->setEta(eta);
635 <    }
636 <
637 <    
638 < }
639 <
640 < }//end namespace oopse
548 >  
549 > }//end namespace oopse

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines