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 246 by gezelter, Wed Jan 12 22:41:40 2005 UTC vs.
Revision 1390 by gezelter, Wed Nov 25 20:02:06 2009 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]  Vardeman & Gezelter, in progress (2009).                        
40   */
41 +  
42 + #define _LARGEFILE_SOURCE64
43 + #define _FILE_OFFSET_BITS 64
44  
45 < #define _LARGEFILE_SOURCE64
46 < #define _FILE_OFFSET_BITS 64
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 OpenMD {
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 >        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 < #include <sys/types.h>
245 < #include <sys/stat.h>
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 < #include <iostream>
261 < #include <math.h>
260 >  }
261 >  
262 >  void DumpReader::readSet(int whichFrame) {    
263 >    std::string line;
264  
265 < #include <stdio.h>
266 < #include <stdlib.h>
267 < #include <string.h>
265 > #ifndef IS_MPI
266 >    inFile_->clear();  
267 >    inFile_->seekg(framePos_[whichFrame]);
268  
269 < #include "io/DumpReader.hpp"
56 < #include "primitives/Molecule.hpp"
57 < #include "utils/simError.h"
58 < #include "utils/MemoryUtils.hpp"
59 < #include "utils/StringTokenizer.hpp"
269 >    std::istream& inputStream = *inFile_;    
270  
271 < #ifdef IS_MPI
271 > #else
272 >    int masterNode = 0;
273 >    std::stringstream sstream;
274 >    if (worldRank == masterNode) {
275 >      std::string sendBuffer;
276  
277 < #include <mpi.h>
278 < #define TAKE_THIS_TAG_CHAR 0
279 < #define TAKE_THIS_TAG_INT 1
277 >      inFile_->clear();  
278 >      inFile_->seekg(framePos_[whichFrame]);
279 >      
280 >      while (inFile_->getline(buffer, bufferSize)) {
281  
282 < #endif // is_mpi
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 +      assert(recvBuffer);
300 +      recvBuffer[sendBufferSize] = '\0';
301 +      MPI_Bcast(recvBuffer, sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);    
302 +      sstream.str(recvBuffer);
303 +      delete [] recvBuffer;
304 +    }      
305  
306 < namespace oopse {
306 >    std::istream& inputStream = sstream;  
307 > #endif
308  
309 < DumpReader::DumpReader(SimInfo* info, const std::string& filename)
73 <                     : info_(info), filename_(filename), isScanned_(false), nframes_(0) {
309 >    inputStream.getline(buffer, bufferSize);
310  
311 < #ifdef IS_MPI
311 >    line = buffer;
312 >    if (line.find("<Snapshot>") == std::string::npos) {
313 >      sprintf(painCave.errMsg,
314 >              "DumpReader Error: can not find <Snapshot>\n");
315 >      painCave.isFatal = 1;
316 >      simError();
317 >    }
318 >    
319 >    //read frameData
320 >    readFrameProperties(inputStream);
321  
322 <    if (worldRank == 0) {
323 < #endif
322 >    //read StuntDoubles
323 >    readStuntDoubles(inputStream);    
324  
325 <        inFile_ = fopen(filename_.c_str(), "r");
325 >    inputStream.getline(buffer, bufferSize);
326 >    line = buffer;
327 >    if (line.find("</Snapshot>") == std::string::npos) {
328 >      sprintf(painCave.errMsg,
329 >              "DumpReader Error: can not find </Snapshot>\n");
330 >      painCave.isFatal = 1;
331 >      simError();
332 >    }        
333 >  
334 >  }
335 >  
336 >  void DumpReader::parseDumpLine(const std::string& line) {
337  
338 <        if (inFile_ == NULL) {
339 <            sprintf(painCave.errMsg, "Cannot open file: %s\n", filename_.c_str());
340 <            painCave.isFatal = 1;
341 <            simError();
342 <        }
338 >      
339 >    StringTokenizer tokenizer(line);
340 >    int nTokens;
341 >    
342 >    nTokens = tokenizer.countTokens();
343 >    
344 >    if (nTokens < 2) {
345 >      sprintf(painCave.errMsg,
346 >              "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
347 >      painCave.isFatal = 1;
348 >      simError();
349 >    }
350  
351 < #ifdef IS_MPI
351 >    int index = tokenizer.nextTokenAsInt();
352 >
353 >    StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index);
354  
355 +    if (integrableObject == NULL) {
356 +      return;
357      }
358 +    std::string type = tokenizer.nextToken();
359 +    int size = type.size();
360  
361 <    strcpy(checkPointMsg, "Dump file opened for reading successfully.");
362 <    MPIcheckPoint();
361 >    for(int i = 0; i < size; ++i) {
362 >      switch(type[i]) {
363 >        
364 >        case 'p': {
365 >            Vector3d pos;
366 >            pos[0] = tokenizer.nextTokenAsDouble();
367 >            pos[1] = tokenizer.nextTokenAsDouble();
368 >            pos[2] = tokenizer.nextTokenAsDouble();
369 >            if (needPos_) {
370 >              integrableObject->setPos(pos);
371 >            }            
372 >            break;
373 >        }
374 >        case 'v' : {
375 >            Vector3d vel;
376 >            vel[0] = tokenizer.nextTokenAsDouble();
377 >            vel[1] = tokenizer.nextTokenAsDouble();
378 >            vel[2] = tokenizer.nextTokenAsDouble();
379 >            if (needVel_) {
380 >              integrableObject->setVel(vel);
381 >            }
382 >            break;
383 >        }
384  
385 < #endif
385 >        case 'q' : {
386 >           Quat4d q;
387 >           if (integrableObject->isDirectional()) {
388 >              
389 >             q[0] = tokenizer.nextTokenAsDouble();
390 >             q[1] = tokenizer.nextTokenAsDouble();
391 >             q[2] = tokenizer.nextTokenAsDouble();
392 >             q[3] = tokenizer.nextTokenAsDouble();
393 >              
394 >             RealType qlen = q.length();
395 >             if (qlen < OpenMD::epsilon) { //check quaternion is not equal to 0
396 >                
397 >               sprintf(painCave.errMsg,
398 >                       "DumpReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n");
399 >               painCave.isFatal = 1;
400 >               simError();
401 >                
402 >             }  
403 >              
404 >             q.normalize();
405 >             if (needQuaternion_) {            
406 >               integrableObject->setQ(q);
407 >             }              
408 >           }            
409 >           break;
410 >        }  
411 >        case 'j' : {
412 >          Vector3d ji;
413 >          if (integrableObject->isDirectional()) {
414 >             ji[0] = tokenizer.nextTokenAsDouble();
415 >             ji[1] = tokenizer.nextTokenAsDouble();
416 >             ji[2] = tokenizer.nextTokenAsDouble();
417 >             if (needAngMom_) {
418 >               integrableObject->setJ(ji);
419 >             }
420 >          }
421 >          break;
422 >        }  
423 >        case 'f': {
424  
425 <    return;
426 < }
425 >          Vector3d force;
426 >          force[0] = tokenizer.nextTokenAsDouble();
427 >          force[1] = tokenizer.nextTokenAsDouble();
428 >          force[2] = tokenizer.nextTokenAsDouble();          
429 >          integrableObject->setFrc(force);
430 >          break;
431 >        }
432 >        case 't' : {
433  
434 < DumpReader::~DumpReader() {
435 <
436 < #ifdef IS_MPI
437 <
438 <    if (worldRank == 0) {
439 < #endif
106 <
107 <        int error;
108 <        error = fclose(inFile_);
109 <
110 <        if (error) {
111 <            sprintf(painCave.errMsg, "Error closing %s\n", filename_.c_str());
112 <            painCave.isFatal = 1;            
113 <            simError();
434 >           Vector3d torque;
435 >           torque[0] = tokenizer.nextTokenAsDouble();
436 >           torque[1] = tokenizer.nextTokenAsDouble();
437 >           torque[2] = tokenizer.nextTokenAsDouble();          
438 >           integrableObject->setTrq(torque);          
439 >           break;
440          }
441 +        default: {
442 +               sprintf(painCave.errMsg,
443 +                       "DumpReader Error: %s is an unrecognized type\n", type.c_str());
444 +               painCave.isFatal = 1;
445 +               simError();
446 +          break;  
447 +        }
448  
449 <        MemoryUtils::deleteVectorOfPointer(framePos_);
449 >      }
450 >    }
451 >    
452 >  }
453 >  
454  
455 < #ifdef IS_MPI
455 >  void  DumpReader::readStuntDoubles(std::istream& inputStream) {
456  
457 +    inputStream.getline(buffer, bufferSize);
458 +    std::string line(buffer);
459 +    
460 +    if (line.find("<StuntDoubles>") == std::string::npos) {
461 +      sprintf(painCave.errMsg,
462 +              "DumpReader Error: Missing <StuntDoubles>\n");
463 +      painCave.isFatal = 1;
464 +      simError();
465      }
466  
467 <    strcpy(checkPointMsg, "Dump file closed successfully.");
468 <    MPIcheckPoint();
467 >    while(inputStream.getline(buffer, bufferSize)) {
468 >      line = buffer;
469 >      
470 >      if(line.find("</StuntDoubles>") != std::string::npos) {
471 >        break;
472 >      }
473  
474 < #endif
474 >      parseDumpLine(line);
475 >    }
476 >  
477 >  }
478  
479 <    return;
128 < }
479 >  void DumpReader::readFrameProperties(std::istream& inputStream) {
480  
481 < int DumpReader::getNFrames(void) {
481 >    Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
482 >    inputStream.getline(buffer, bufferSize);
483 >    std::string line(buffer);
484  
485 <    if (!isScanned_)
486 <        scanFile();
485 >    if (line.find("<FrameData>") == std::string::npos) {
486 >      sprintf(painCave.errMsg,
487 >              "DumpReader Error: Missing <FrameData>\n");
488 >      painCave.isFatal = 1;
489 >      simError();
490 >    }
491  
492 <    return nframes_;
493 < }
494 <
495 < void DumpReader::scanFile(void) {
496 <  int i, j;
497 <  int lineNum = 0;
498 <  char readBuffer[maxBufferSize];
499 <  fpos_t * currPos;
500 <
501 < #ifdef IS_MPI
492 >    while(inputStream.getline(buffer, bufferSize)) {
493 >      line = buffer;
494 >      
495 >      if(line.find("</FrameData>") != std::string::npos) {
496 >        break;
497 >      }
498 >      
499 >      StringTokenizer tokenizer(line, " ;\t\n\r{}:,");
500 >      if (!tokenizer.hasMoreTokens()) {
501 >        sprintf(painCave.errMsg,
502 >                "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
503 >        painCave.isFatal = 1;
504 >        simError();      
505 >      }
506  
507 <    if (worldRank == 0) {
508 < #endif // is_mpi
509 <
510 <        rewind(inFile_);
511 <
512 <        currPos = new fpos_t;
513 <        fgetpos(inFile_, currPos);
514 <        fgets(readBuffer, sizeof(readBuffer), inFile_);
515 <        lineNum++;
516 <
517 <        if (feof(inFile_)) {
518 <            sprintf(painCave.errMsg,
519 <                    "File \"%s\" ended unexpectedly at line %d\n",
520 <                    filename_.c_str(),
521 <                    lineNum);
522 <            painCave.isFatal = 1;
523 <            simError();
524 <        }
525 <
526 <        while (!feof(inFile_)) {
527 <            framePos_.push_back(currPos);
528 <
529 <            i = atoi(readBuffer);
530 <
531 <            fgets(readBuffer, sizeof(readBuffer), inFile_);
532 <            lineNum++;
533 <
534 <            if (feof(inFile_)) {
535 <                sprintf(painCave.errMsg,
536 <                        "File \"%s\" ended unexpectedly at line %d\n",
537 <                        filename_.c_str(),
538 <                        lineNum);
539 <                painCave.isFatal = 1;
540 <                simError();
541 <            }
542 <
543 <            for(j = 0; j < i; j++) {
544 <                fgets(readBuffer, sizeof(readBuffer), inFile_);
545 <                lineNum++;
546 <
186 <                if (feof(inFile_)) {
187 <                    sprintf(painCave.errMsg,
188 <                            "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
507 >      std::string propertyName = tokenizer.nextToken();
508 >      if (propertyName == "Time") {
509 >        RealType currTime = tokenizer.nextTokenAsDouble();
510 >        s->setTime(currTime);
511 >      } else if (propertyName == "Hmat"){
512 >        Mat3x3d hmat;
513 >        hmat(0, 0) = tokenizer.nextTokenAsDouble();
514 >        hmat(0, 1) = tokenizer.nextTokenAsDouble();
515 >        hmat(0, 2) = tokenizer.nextTokenAsDouble();
516 >        hmat(1, 0) = tokenizer.nextTokenAsDouble();
517 >        hmat(1, 1) = tokenizer.nextTokenAsDouble();
518 >        hmat(1, 2) = tokenizer.nextTokenAsDouble();
519 >        hmat(2, 0) = tokenizer.nextTokenAsDouble();
520 >        hmat(2, 1) = tokenizer.nextTokenAsDouble();
521 >        hmat(2, 2) = tokenizer.nextTokenAsDouble();
522 >        s->setHmat(hmat);      
523 >      } else if (propertyName == "Thermostat") {
524 >        RealType chi = tokenizer.nextTokenAsDouble();
525 >        RealType integralOfChiDt = tokenizer.nextTokenAsDouble();
526 >        s->setChi(chi);
527 >        s->setIntegralOfChiDt(integralOfChiDt);        
528 >     } else if (propertyName == "Barostat") {
529 >        Mat3x3d eta;
530 >        eta(0, 0) = tokenizer.nextTokenAsDouble();
531 >        eta(0, 1) = tokenizer.nextTokenAsDouble();
532 >        eta(0, 2) = tokenizer.nextTokenAsDouble();
533 >        eta(1, 0) = tokenizer.nextTokenAsDouble();
534 >        eta(1, 1) = tokenizer.nextTokenAsDouble();
535 >        eta(1, 2) = tokenizer.nextTokenAsDouble();
536 >        eta(2, 0) = tokenizer.nextTokenAsDouble();
537 >        eta(2, 1) = tokenizer.nextTokenAsDouble();
538 >        eta(2, 2) = tokenizer.nextTokenAsDouble();
539 >        s->setEta(eta);
540 >      } else {
541 >        sprintf(painCave.errMsg,
542 >                "DumpReader Error: %s is an invalid property in <FrameData>\n", propertyName.c_str());
543 >        painCave.isFatal = 0;
544 >        simError();        
545 >      }
546 >      
547      }
548  
549 <    MPI_Bcast(&nframes_, 1, MPI_INT, 0, MPI_COMM_WORLD);
549 >  }
550  
551 <    strcpy(checkPointMsg, "Successfully scanned DumpFile\n");
552 <    MPIcheckPoint();
215 <
216 < #endif // is_mpi
217 <
218 <    isScanned_ = true;
219 < }
220 <
221 < void DumpReader::readFrame(int whichFrame) {
222 <    readSet(whichFrame);
223 < }
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();
247 <    }
248 <
249 <    nTotObjs = atoi(read_buffer);
250 <
251 <    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_->getNIntegrableObjects());
258 <
259 <        painCave.isFatal = 1;
260 <        simError();
261 <    }
262 <
263 <    //read the box mat from the comment line
264 <
265 <    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
266 <
267 <    if (eof_test == NULL) {
268 <        sprintf(painCave.errMsg, "error in reading commment in %s\n",
269 <                filename_.c_str());
270 <        painCave.isFatal = 1;
271 <        simError();
272 <    }
273 <
274 <    parseCommentLine(read_buffer, info_->getSnapshotManager()->getCurrentSnapshot());
275 <
276 <    //parse dump lines
277 <
278 <    for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
279 <
280 <        for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;
281 <            integrableObject = mol->nextIntegrableObject(ii)) {          
282 <
283 <            eof_test = fgets(read_buffer, sizeof(read_buffer), inFile_);
284 <
285 <            if (eof_test == NULL) {
286 <                sprintf(painCave.errMsg,
287 <                        "error in reading file %s\n"
288 <                            "natoms  = %d; index = %d\n"
289 <                            "error reading the line from the file.\n",
290 <                        filename_.c_str(),
291 <                        nTotObjs,
292 <                        i);
293 <
294 <                painCave.isFatal = 1;
295 <                simError();
296 <            }
297 <
298 <            parseDumpLine(read_buffer, integrableObject);
299 <            
300 <            }
301 <    }
302 <
303 <    // 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_);
325 <
326 <        if (eof_test == NULL) {
327 <            sprintf(painCave.errMsg, "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, "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, "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 <                                "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 <                                "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, "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 <                    "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 <    }
498 <
499 <    pos[0] = tokenizer.nextTokenAsDouble();
500 <    pos[1] = tokenizer.nextTokenAsDouble();
501 <    pos[2] = tokenizer.nextTokenAsDouble();
502 <    integrableObject->setPos(pos);
503 <    
504 <    vel[0] = tokenizer.nextTokenAsDouble();
505 <    vel[1] = tokenizer.nextTokenAsDouble();
506 <    vel[2] = tokenizer.nextTokenAsDouble();
507 <    integrableObject->setVel(vel);
508 <
509 <    if (integrableObject->isDirectional()) {
510 <        
511 <        q[0] = tokenizer.nextTokenAsDouble();
512 <        q[1] = tokenizer.nextTokenAsDouble();
513 <        q[2] = tokenizer.nextTokenAsDouble();
514 <        q[3] = tokenizer.nextTokenAsDouble();
515 <
516 <        double qlen = q.length();
517 <        if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
518 <            
519 <            sprintf(painCave.errMsg,
520 <                    "initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2 ~ 0).\n");
521 <            painCave.isFatal = 1;
522 <            simError();
523 <            
524 <        }
525 <
526 <        q.normalize();
527 <              
528 <        integrableObject->setQ(q);
529 <        
530 <        ji[0] = tokenizer.nextTokenAsDouble();
531 <        ji[1] = tokenizer.nextTokenAsDouble();
532 <        ji[2] = tokenizer.nextTokenAsDouble();
533 <        integrableObject->setJ(ji);
534 <    }
535 <
536 < }
537 <
538 <
539 < void DumpReader::parseCommentLine(char* line, Snapshot* s) {
540 <    double currTime;
541 <    Mat3x3d hmat;
542 <    double chi;
543 <    double integralOfChiDt;
544 <    Mat3x3d eta;
545 <
546 <    StringTokenizer tokenizer(line);
547 <    int nTokens;
548 <
549 <    nTokens = tokenizer.countTokens();
550 <
551 <    //comment line should at least contain 10 tokens: current time(1 token) and  h-matrix(9 tokens)
552 <    if (nTokens < 10) {
553 <            sprintf(painCave.errMsg,
554 <                    "Not enough tokens in comment line: %s", line);
555 <            painCave.isFatal = 1;
556 <            simError();  
557 <    }
558 <
559 <    //read current time
560 <    currTime = tokenizer.nextTokenAsDouble();
561 <    s->setTime(currTime);
562 <    
563 <    //read h-matrix
564 <    hmat(0, 0) = tokenizer.nextTokenAsDouble();
565 <    hmat(0, 1) = tokenizer.nextTokenAsDouble();
566 <    hmat(0, 2) = tokenizer.nextTokenAsDouble();
567 <    hmat(1, 0) = tokenizer.nextTokenAsDouble();
568 <    hmat(1, 1) = tokenizer.nextTokenAsDouble();
569 <    hmat(1, 2) = tokenizer.nextTokenAsDouble();
570 <    hmat(2, 0) = tokenizer.nextTokenAsDouble();
571 <    hmat(2, 1) = tokenizer.nextTokenAsDouble();
572 <    hmat(2, 2) = tokenizer.nextTokenAsDouble();
573 <    s->setHmat(hmat);
574 <    
575 <    //read chi and integrablOfChidt, they should apprear in pair
576 <    if (tokenizer.countTokens() >= 2) {
577 <        chi = tokenizer.nextTokenAsDouble();
578 <        integralOfChiDt = tokenizer.nextTokenAsDouble();            
579 <
580 <        s->setChi(chi);
581 <        s->setIntegralOfChiDt(integralOfChiDt);
582 <    }
583 <    
584 <    //read eta (eta is 3x3 matrix)
585 <    if (tokenizer.countTokens() >= 9) {
586 <        eta(0, 0) = tokenizer.nextTokenAsDouble();
587 <        eta(0, 1) = tokenizer.nextTokenAsDouble();
588 <        eta(0, 2) = tokenizer.nextTokenAsDouble();
589 <        eta(1, 0) = tokenizer.nextTokenAsDouble();
590 <        eta(1, 1) = tokenizer.nextTokenAsDouble();
591 <        eta(1, 2) = tokenizer.nextTokenAsDouble();
592 <        eta(2, 0) = tokenizer.nextTokenAsDouble();
593 <        eta(2, 1) = tokenizer.nextTokenAsDouble();
594 <        eta(2, 2) = tokenizer.nextTokenAsDouble();      
595 <
596 <        s->setEta(eta);
597 <    }
598 <
599 <    
600 < }
601 <
602 < }//end namespace oopse
551 >  
552 > }//end namespace OpenMD

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines