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

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

Comparing trunk/src/io/DumpReader.cpp (property svn:keywords):
Revision 721 by chrisfen, Tue Nov 8 13:32:27 2005 UTC vs.
Revision 2000 by gezelter, Sat May 31 22:35:05 2014 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines