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 1025 by gezelter, Wed Aug 30 20:33:44 2006 UTC vs.
Revision 2069 by gezelter, Thu Mar 5 16:30:23 2015 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) {
71 >    : info_(info), filename_(filename), isScanned_(false), nframes_(0), needCOMprops_(false) {
72      
73   #ifdef IS_MPI
74      
75      if (worldRank == 0) {
76   #endif
77        
78 <      inFile_ = new std::ifstream(filename_.c_str());
78 >      inFile_ = new std::ifstream(filename_.c_str(),  
79 >                                  ifstream::in | ifstream::binary);
80        
81        if (inFile_->fail()) {
82          sprintf(painCave.errMsg,
# Line 92 | Line 91 | namespace oopse {
91      }
92      
93      strcpy(checkPointMsg, "Dump file opened for reading successfully.");
94 <    MPIcheckPoint();
94 >    errorCheckPoint();
95      
96   #endif
97      
# Line 104 | Line 103 | namespace oopse {
103   #ifdef IS_MPI
104      
105      if (worldRank == 0) {
106 < #endif
107 <      
106 > #endif
107 >
108        delete inFile_;
109        
110   #ifdef IS_MPI
# Line 113 | Line 112 | namespace oopse {
112      }
113      
114      strcpy(checkPointMsg, "Dump file closed successfully.");
115 <    MPIcheckPoint();
115 >    errorCheckPoint();
116      
117   #endif
118      
# Line 129 | Line 128 | namespace oopse {
128    }
129    
130    void DumpReader::scanFile(void) {
131 <    int lineNo = 0;
131 >
132      std::streampos prevPos;
133      std::streampos  currPos;
134      
# Line 142 | Line 141 | namespace oopse {
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          
# Line 227 | 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 240 | Line 243 | namespace oopse {
243      }
244      
245      readSet(whichFrame);
246 +
247 +    if (needCOMprops_) {
248 +      Thermo thermo(info_);
249 +      Vector3d com;
250 +
251 +      if (needPos_ && needVel_) {
252 +        Vector3d comvel;
253 +        Vector3d comw;
254 +        thermo.getComAll(com, comvel);
255 +        comw = thermo.getAngularMomentum();
256 +      } else {
257 +        com = thermo.getCom();
258 +      }                    
259 +    }
260    }
261    
262    void DumpReader::readSet(int whichFrame) {    
263      std::string line;
264  
265   #ifndef IS_MPI
249
266      inFile_->clear();  
267      inFile_->seekg(framePos_[whichFrame]);
268  
# Line 273 | Line 289 | namespace oopse {
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);    
292 >      MPI_Bcast((void *)sendBuffer.c_str(), sendBufferSize,
293 >                MPI_CHAR, masterNode, MPI_COMM_WORLD);    
294        
295        sstream.str(sendBuffer);
296      } else {
297        int sendBufferSize;
298 <      MPI_Bcast(&sendBufferSize, 1, MPI_INT, masterNode, MPI_COMM_WORLD);    
298 >      MPI_Bcast(&sendBufferSize, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
299        char * recvBuffer = new char[sendBufferSize+1];
300 <      MPI_Bcast(recvBuffer, sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);    
300 >      assert(recvBuffer);
301 >      recvBuffer[sendBufferSize] = '\0';
302 >      MPI_Bcast(recvBuffer, sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
303        sstream.str(recvBuffer);
304 +      delete [] recvBuffer;
305      }      
306  
307      std::istream& inputStream = sstream;  
# Line 305 | Line 325 | namespace oopse {
325  
326      inputStream.getline(buffer, bufferSize);
327      line = buffer;
328 <    if (line.find("</Snapshot>") == std::string::npos) {
329 <      sprintf(painCave.errMsg,
330 <              "DumpReader Error: can not find </Snapshot>\n");
331 <      painCave.isFatal = 1;
332 <      simError();
333 <    }        
334 <  
328 >
329 >    if (line.find("<SiteData>") != std::string::npos) {
330 >      //read SiteData
331 >      readSiteData(inputStream);        
332 >    } else {
333 >      if (line.find("</Snapshot>") == std::string::npos) {
334 >        sprintf(painCave.errMsg,
335 >                "DumpReader Error: can not find </Snapshot>\n");
336 >        painCave.isFatal = 1;
337 >        simError();
338 >      }        
339 >    }
340    }
341    
342    void DumpReader::parseDumpLine(const std::string& line) {
# Line 322 | Line 347 | namespace oopse {
347      
348      nTokens = tokenizer.countTokens();
349      
350 <    if (nTokens < 2) {
350 >    if (nTokens < 2) {  
351        sprintf(painCave.errMsg,
352                "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
353        painCave.isFatal = 1;
# Line 331 | Line 356 | namespace oopse {
356  
357      int index = tokenizer.nextTokenAsInt();
358  
359 <    StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index);
359 >    StuntDouble* sd = info_->getIOIndexToIntegrableObject(index);
360  
361 <    if (integrableObject == NULL) {
361 >    if (sd == NULL) {
362        return;
363      }
364      std::string type = tokenizer.nextToken();
365      int size = type.size();
366 +
367 +    size_t found;
368 +    
369 +    if (needPos_) {
370 +      found = type.find("p");      
371 +      if (found == std::string::npos) {
372 +        sprintf(painCave.errMsg,
373 +                "DumpReader Error: StuntDouble %d has no Position\n"
374 +                "\tField (\"p\") specified.\n%s\n", index,
375 +                line.c_str());  
376 +        painCave.isFatal = 1;
377 +        simError();
378 +      }
379 +    }
380 +    
381 +    if (sd->isDirectional()) {
382 +      if (needQuaternion_) {
383 +        found = type.find("q");      
384 +        if (found == std::string::npos) {
385 +          sprintf(painCave.errMsg,
386 +                  "DumpReader Error: Directional StuntDouble %d has no\n"
387 +                  "\tQuaternion Field (\"q\") specified.\n%s\n", index,
388 +                  line.c_str());  
389 +          painCave.isFatal = 1;
390 +          simError();
391 +        }
392 +      }      
393 +    }
394 +
395      for(int i = 0; i < size; ++i) {
396        switch(type[i]) {
397          
# Line 347 | Line 401 | namespace oopse {
401              pos[1] = tokenizer.nextTokenAsDouble();
402              pos[2] = tokenizer.nextTokenAsDouble();
403              if (needPos_) {
404 <              integrableObject->setPos(pos);
404 >              sd->setPos(pos);
405              }            
406              break;
407          }
# Line 357 | Line 411 | namespace oopse {
411              vel[1] = tokenizer.nextTokenAsDouble();
412              vel[2] = tokenizer.nextTokenAsDouble();
413              if (needVel_) {
414 <              integrableObject->setVel(vel);
414 >              sd->setVel(vel);
415              }
416              break;
417          }
418  
419          case 'q' : {
420             Quat4d q;
421 <           if (integrableObject->isDirectional()) {
421 >           if (sd->isDirectional()) {
422                
423               q[0] = tokenizer.nextTokenAsDouble();
424               q[1] = tokenizer.nextTokenAsDouble();
# Line 372 | Line 426 | namespace oopse {
426               q[3] = tokenizer.nextTokenAsDouble();
427                
428               RealType qlen = q.length();
429 <             if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
429 >             if (qlen < OpenMD::epsilon) { //check quaternion is not equal to 0
430                  
431                 sprintf(painCave.errMsg,
432                         "DumpReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n");
# Line 383 | Line 437 | namespace oopse {
437                
438               q.normalize();
439               if (needQuaternion_) {            
440 <               integrableObject->setQ(q);
440 >               sd->setQ(q);
441               }              
442             }            
443             break;
444          }  
445          case 'j' : {
446            Vector3d ji;
447 <          if (integrableObject->isDirectional()) {
447 >          if (sd->isDirectional()) {
448               ji[0] = tokenizer.nextTokenAsDouble();
449               ji[1] = tokenizer.nextTokenAsDouble();
450               ji[2] = tokenizer.nextTokenAsDouble();
451               if (needAngMom_) {
452 <               integrableObject->setJ(ji);
452 >               sd->setJ(ji);
453               }
454            }
455 +          break;
456          }  
457          case 'f': {
458  
# Line 405 | Line 460 | namespace oopse {
460            force[0] = tokenizer.nextTokenAsDouble();
461            force[1] = tokenizer.nextTokenAsDouble();
462            force[2] = tokenizer.nextTokenAsDouble();          
463 <          integrableObject->setFrc(force);
463 >          sd->setFrc(force);
464            break;
465          }
466          case 't' : {
# Line 414 | Line 469 | namespace oopse {
469             torque[0] = tokenizer.nextTokenAsDouble();
470             torque[1] = tokenizer.nextTokenAsDouble();
471             torque[2] = tokenizer.nextTokenAsDouble();          
472 <           integrableObject->setTrq(torque);          
472 >           sd->setTrq(torque);          
473             break;
474          }
475 +        case 'u' : {
476 +
477 +           RealType particlePot;
478 +           particlePot = tokenizer.nextTokenAsDouble();
479 +           sd->setParticlePot(particlePot);          
480 +           break;
481 +        }
482 +        case 'c' : {
483 +
484 +           RealType flucQPos;
485 +           flucQPos = tokenizer.nextTokenAsDouble();
486 +           sd->setFlucQPos(flucQPos);          
487 +           break;
488 +        }
489 +        case 'w' : {
490 +
491 +           RealType flucQVel;
492 +           flucQVel = tokenizer.nextTokenAsDouble();
493 +           sd->setFlucQVel(flucQVel);          
494 +           break;
495 +        }
496 +        case 'g' : {
497 +
498 +           RealType flucQFrc;
499 +           flucQFrc = tokenizer.nextTokenAsDouble();
500 +           sd->setFlucQFrc(flucQFrc);          
501 +           break;
502 +        }
503 +        case 'e' : {
504 +
505 +           Vector3d eField;
506 +           eField[0] = tokenizer.nextTokenAsDouble();
507 +           eField[1] = tokenizer.nextTokenAsDouble();
508 +           eField[2] = tokenizer.nextTokenAsDouble();          
509 +           sd->setElectricField(eField);          
510 +           break;
511 +        }
512 +        case 's' : {
513 +
514 +           RealType sPot;
515 +           sPot = tokenizer.nextTokenAsDouble();
516 +           sd->setSitePotential(sPot);          
517 +           break;
518 +        }
519          default: {
520                 sprintf(painCave.errMsg,
521                         "DumpReader Error: %s is an unrecognized type\n", type.c_str());
# Line 427 | Line 526 | namespace oopse {
526  
527        }
528      }
430    
529    }
530    
531  
532 <  void  DumpReader::readStuntDoubles(std::istream& inputStream) {
532 >  void DumpReader::parseSiteLine(const std::string& line) {
533  
534 +    StringTokenizer tokenizer(line);
535 +    int nTokens;
536 +        
537 +    nTokens = tokenizer.countTokens();
538 +    
539 +    if (nTokens < 2) {  
540 +      sprintf(painCave.errMsg,
541 +              "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
542 +      painCave.isFatal = 1;
543 +      simError();
544 +    }
545 +
546 +    /**
547 +     * The first token is the global integrable object index.
548 +     */
549 +
550 +    int index = tokenizer.nextTokenAsInt();
551 +    StuntDouble* sd = info_->getIOIndexToIntegrableObject(index);
552 +    if (sd == NULL) {
553 +      return;
554 +    }
555 +
556 +    /**
557 +     * Test to see if the next token is an integer or not.  If not,
558 +     * we've got data on the integrable object itself.  If there is an
559 +     * integer, we're parsing data for a site on a rigid body.
560 +     */
561 +
562 +    std::string indexTest = tokenizer.peekNextToken();
563 +    std::istringstream i(indexTest);
564 +    int siteIndex;
565 +    if (i >> siteIndex) {
566 +      // chew up this token and parse as an int:
567 +      siteIndex = tokenizer.nextTokenAsInt();
568 +
569 +      if (sd->isRigidBody()) {
570 +        RigidBody* rb = static_cast<RigidBody*>(sd);
571 +        sd = rb->getAtoms()[siteIndex];
572 +      }
573 +    }
574 +
575 +    /**
576 +     * The next token contains information on what follows.
577 +     */
578 +    std::string type = tokenizer.nextToken();
579 +    int size = type.size();
580 +    
581 +    for(int i = 0; i < size; ++i) {
582 +      switch(type[i]) {
583 +        
584 +      case 'u' : {
585 +        
586 +        RealType particlePot;
587 +        particlePot = tokenizer.nextTokenAsDouble();
588 +        sd->setParticlePot(particlePot);
589 +        break;
590 +      }
591 +      case 'c' : {
592 +        
593 +        RealType flucQPos;
594 +        flucQPos = tokenizer.nextTokenAsDouble();
595 +        sd->setFlucQPos(flucQPos);
596 +        break;
597 +      }
598 +      case 'w' : {
599 +        
600 +        RealType flucQVel;
601 +        flucQVel = tokenizer.nextTokenAsDouble();
602 +        sd->setFlucQVel(flucQVel);
603 +        break;
604 +      }
605 +      case 'g' : {
606 +        
607 +        RealType flucQFrc;
608 +        flucQFrc = tokenizer.nextTokenAsDouble();
609 +        sd->setFlucQFrc(flucQFrc);
610 +        break;
611 +      }
612 +      case 'e' : {
613 +        
614 +        Vector3d eField;
615 +        eField[0] = tokenizer.nextTokenAsDouble();
616 +        eField[1] = tokenizer.nextTokenAsDouble();
617 +        eField[2] = tokenizer.nextTokenAsDouble();  
618 +        sd->setElectricField(eField);          
619 +        break;
620 +      }
621 +      case 's' : {
622 +        
623 +        RealType sPot;
624 +        sPot = tokenizer.nextTokenAsDouble();
625 +        sd->setSitePotential(sPot);          
626 +        break;
627 +      }
628 +      default: {
629 +        sprintf(painCave.errMsg,
630 +                "DumpReader Error: %s is an unrecognized type\n", type.c_str());
631 +        painCave.isFatal = 1;
632 +        simError();
633 +        break;  
634 +      }
635 +      }
636 +    }    
637 +  }
638 +  
639 +  
640 +    void  DumpReader::readStuntDoubles(std::istream& inputStream) {
641 +    
642      inputStream.getline(buffer, bufferSize);
643      std::string line(buffer);
644      
# Line 455 | Line 661 | namespace oopse {
661    
662    }
663  
664 +  void  DumpReader::readSiteData(std::istream& inputStream) {
665 +
666 +    std::string line(buffer);
667 +
668 +    // We already found the starting <SiteData> tag or we wouldn't be
669 +    // here, so just start parsing until we get to the ending
670 +    // </SiteData> tag:
671 +    
672 +    while(inputStream.getline(buffer, bufferSize)) {
673 +      line = buffer;
674 +      
675 +      if(line.find("</SiteData>") != std::string::npos) {
676 +        break;
677 +      }
678 +
679 +      parseSiteLine(line);
680 +    }
681 +  
682 +  }
683 +
684    void DumpReader::readFrameProperties(std::istream& inputStream) {
685  
686      Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
# Line 500 | Line 726 | namespace oopse {
726          hmat(2, 2) = tokenizer.nextTokenAsDouble();
727          s->setHmat(hmat);      
728        } else if (propertyName == "Thermostat") {
729 <        RealType chi = tokenizer.nextTokenAsDouble();
730 <        RealType integralOfChiDt = tokenizer.nextTokenAsDouble();
731 <        s->setChi(chi);
732 <        s->setIntegralOfChiDt(integralOfChiDt);        
729 >        pair<RealType, RealType> thermostat;
730 >        thermostat.first = tokenizer.nextTokenAsDouble();
731 >        thermostat.second = tokenizer.nextTokenAsDouble();
732 >        s->setThermostat(thermostat);
733       } else if (propertyName == "Barostat") {
734          Mat3x3d eta;
735          eta(0, 0) = tokenizer.nextTokenAsDouble();
# Line 515 | Line 741 | namespace oopse {
741          eta(2, 0) = tokenizer.nextTokenAsDouble();
742          eta(2, 1) = tokenizer.nextTokenAsDouble();
743          eta(2, 2) = tokenizer.nextTokenAsDouble();
744 <        s->setEta(eta);
744 >        s->setBarostat(eta);
745        } else {
746          sprintf(painCave.errMsg,
747                  "DumpReader Error: %s is an invalid property in <FrameData>\n", propertyName.c_str());
# Line 528 | Line 754 | namespace oopse {
754    }
755  
756    
757 < }//end namespace oopse
757 > }//end namespace OpenMD

Comparing trunk/src/io/DumpReader.cpp (property svn:keywords):
Revision 1025 by gezelter, Wed Aug 30 20:33:44 2006 UTC vs.
Revision 2069 by gezelter, Thu Mar 5 16:30:23 2015 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines