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 274 by tim, Tue Jan 25 21:59:18 2005 UTC vs.
branches/development/src/io/DumpReader.cpp (file contents), Revision 1752 by gezelter, Sun Jun 10 14:05:02 2012 UTC

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

Comparing:
trunk/src/io/DumpReader.cpp (property svn:keywords), Revision 274 by tim, Tue Jan 25 21:59:18 2005 UTC vs.
branches/development/src/io/DumpReader.cpp (property svn:keywords), Revision 1752 by gezelter, Sun Jun 10 14:05:02 2012 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines