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 1790 by gezelter, Thu Aug 30 17:18:22 2012 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, 24107 (2008).          
39 > * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 > * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41 > */
42    
43   #define _LARGEFILE_SOURCE64
44   #define _FILE_OFFSET_BITS 64
# Line 57 | Line 58
58   #include "utils/simError.h"
59   #include "utils/MemoryUtils.hpp"
60   #include "utils/StringTokenizer.hpp"
61 + #include "brains/Thermo.hpp"
62  
63   #ifdef IS_MPI
62
64   #include <mpi.h>
65 < #define TAKE_THIS_TAG_CHAR 0
65 < #define TAKE_THIS_TAG_INT 1
65 > #endif
66  
67 #endif // is_mpi
67  
68 <
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 113 | Line 112 | namespace oopse {
112      }
113      
114      strcpy(checkPointMsg, "Dump file closed successfully.");
115 <    MPIcheckPoint();
115 >    errorCheckPoint();
116      
117   #endif
118      
# Line 142 | Line 141 | namespace oopse {
141        prevPos = currPos;
142        bool foundOpenSnapshotTag = false;
143        bool foundClosedSnapshotTag = false;
144 +      bool foundOpenSiteDataTag = false;
145        while(inFile_->getline(buffer, bufferSize)) {
146          ++lineNo;
147          
# Line 240 | Line 240 | namespace oopse {
240      }
241      
242      readSet(whichFrame);
243 +
244 +    if (needCOMprops_) {
245 +      Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
246 +      Thermo thermo(info_);
247 +      Vector3d com;
248 +
249 +      if (needPos_ && needVel_) {
250 +        Vector3d comvel;
251 +        Vector3d comw;
252 +        thermo.getComAll(com, comvel);
253 +        comw = thermo.getAngularMomentum();
254 +      } else {
255 +        com = thermo.getCom();
256 +      }                    
257 +    }
258    }
259    
260    void DumpReader::readSet(int whichFrame) {    
261      std::string line;
262  
263   #ifndef IS_MPI
249
264      inFile_->clear();  
265      inFile_->seekg(framePos_[whichFrame]);
266  
# Line 280 | Line 294 | namespace oopse {
294        int sendBufferSize;
295        MPI_Bcast(&sendBufferSize, 1, MPI_INT, masterNode, MPI_COMM_WORLD);    
296        char * recvBuffer = new char[sendBufferSize+1];
297 +      assert(recvBuffer);
298 +      recvBuffer[sendBufferSize] = '\0';
299        MPI_Bcast(recvBuffer, sendBufferSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);    
300        sstream.str(recvBuffer);
301 +      delete [] recvBuffer;
302      }      
303  
304      std::istream& inputStream = sstream;  
# Line 305 | Line 322 | namespace oopse {
322  
323      inputStream.getline(buffer, bufferSize);
324      line = buffer;
325 <    if (line.find("</Snapshot>") == std::string::npos) {
326 <      sprintf(painCave.errMsg,
327 <              "DumpReader Error: can not find </Snapshot>\n");
328 <      painCave.isFatal = 1;
329 <      simError();
330 <    }        
331 <  
325 >
326 >    if (line.find("<SiteData>") != std::string::npos) {
327 >      //read SiteData
328 >      readSiteData(inputStream);        
329 >    } else {
330 >      if (line.find("</Snapshot>") == std::string::npos) {
331 >        sprintf(painCave.errMsg,
332 >                "DumpReader Error: can not find </Snapshot>\n");
333 >        painCave.isFatal = 1;
334 >        simError();
335 >      }        
336 >    }
337    }
338    
339    void DumpReader::parseDumpLine(const std::string& line) {
# Line 322 | Line 344 | namespace oopse {
344      
345      nTokens = tokenizer.countTokens();
346      
347 <    if (nTokens < 2) {
347 >    if (nTokens < 2) {  
348        sprintf(painCave.errMsg,
349                "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
350        painCave.isFatal = 1;
# Line 331 | Line 353 | namespace oopse {
353  
354      int index = tokenizer.nextTokenAsInt();
355  
356 <    StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index);
356 >    StuntDouble* sd = info_->getIOIndexToIntegrableObject(index);
357  
358 <    if (integrableObject == NULL) {
358 >    if (sd == NULL) {
359        return;
360      }
361      std::string type = tokenizer.nextToken();
362      int size = type.size();
363 +
364 +    size_t found;
365 +    
366 +    if (needPos_) {
367 +      found = type.find("p");      
368 +      if (found == std::string::npos) {
369 +        sprintf(painCave.errMsg,
370 +                "DumpReader Error: StuntDouble %d has no Position\n"
371 +                "\tField (\"p\") specified.\n%s\n", index,
372 +                line.c_str());  
373 +        painCave.isFatal = 1;
374 +        simError();
375 +      }
376 +    }
377 +    
378 +    if (sd->isDirectional()) {
379 +      if (needQuaternion_) {
380 +        found = type.find("q");      
381 +        if (found == std::string::npos) {
382 +          sprintf(painCave.errMsg,
383 +                  "DumpReader Error: Directional StuntDouble %d has no\n"
384 +                  "\tQuaternion Field (\"q\") specified.\n%s\n", index,
385 +                  line.c_str());  
386 +          painCave.isFatal = 1;
387 +          simError();
388 +        }
389 +      }      
390 +    }
391 +
392      for(int i = 0; i < size; ++i) {
393        switch(type[i]) {
394          
# Line 347 | Line 398 | namespace oopse {
398              pos[1] = tokenizer.nextTokenAsDouble();
399              pos[2] = tokenizer.nextTokenAsDouble();
400              if (needPos_) {
401 <              integrableObject->setPos(pos);
401 >              sd->setPos(pos);
402              }            
403              break;
404          }
# Line 357 | Line 408 | namespace oopse {
408              vel[1] = tokenizer.nextTokenAsDouble();
409              vel[2] = tokenizer.nextTokenAsDouble();
410              if (needVel_) {
411 <              integrableObject->setVel(vel);
411 >              sd->setVel(vel);
412              }
413              break;
414          }
415  
416          case 'q' : {
417             Quat4d q;
418 <           if (integrableObject->isDirectional()) {
418 >           if (sd->isDirectional()) {
419                
420               q[0] = tokenizer.nextTokenAsDouble();
421               q[1] = tokenizer.nextTokenAsDouble();
# Line 372 | Line 423 | namespace oopse {
423               q[3] = tokenizer.nextTokenAsDouble();
424                
425               RealType qlen = q.length();
426 <             if (qlen < oopse::epsilon) { //check quaternion is not equal to 0
426 >             if (qlen < OpenMD::epsilon) { //check quaternion is not equal to 0
427                  
428                 sprintf(painCave.errMsg,
429                         "DumpReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n");
# Line 383 | Line 434 | namespace oopse {
434                
435               q.normalize();
436               if (needQuaternion_) {            
437 <               integrableObject->setQ(q);
437 >               sd->setQ(q);
438               }              
439             }            
440             break;
441          }  
442          case 'j' : {
443            Vector3d ji;
444 <          if (integrableObject->isDirectional()) {
444 >          if (sd->isDirectional()) {
445               ji[0] = tokenizer.nextTokenAsDouble();
446               ji[1] = tokenizer.nextTokenAsDouble();
447               ji[2] = tokenizer.nextTokenAsDouble();
448               if (needAngMom_) {
449 <               integrableObject->setJ(ji);
449 >               sd->setJ(ji);
450               }
451            }
452 +          break;
453          }  
454          case 'f': {
455  
# Line 405 | Line 457 | namespace oopse {
457            force[0] = tokenizer.nextTokenAsDouble();
458            force[1] = tokenizer.nextTokenAsDouble();
459            force[2] = tokenizer.nextTokenAsDouble();          
460 <          integrableObject->setFrc(force);
460 >          sd->setFrc(force);
461            break;
462          }
463          case 't' : {
# Line 414 | Line 466 | namespace oopse {
466             torque[0] = tokenizer.nextTokenAsDouble();
467             torque[1] = tokenizer.nextTokenAsDouble();
468             torque[2] = tokenizer.nextTokenAsDouble();          
469 <           integrableObject->setTrq(torque);          
469 >           sd->setTrq(torque);          
470             break;
471          }
472 +        case 'u' : {
473 +
474 +           RealType particlePot;
475 +           particlePot = tokenizer.nextTokenAsDouble();
476 +           sd->setParticlePot(particlePot);          
477 +           break;
478 +        }
479 +        case 'c' : {
480 +
481 +           RealType flucQPos;
482 +           flucQPos = tokenizer.nextTokenAsDouble();
483 +           sd->setFlucQPos(flucQPos);          
484 +           break;
485 +        }
486 +        case 'w' : {
487 +
488 +           RealType flucQVel;
489 +           flucQVel = tokenizer.nextTokenAsDouble();
490 +           sd->setFlucQVel(flucQVel);          
491 +           break;
492 +        }
493 +        case 'g' : {
494 +
495 +           RealType flucQFrc;
496 +           flucQFrc = tokenizer.nextTokenAsDouble();
497 +           sd->setFlucQFrc(flucQFrc);          
498 +           break;
499 +        }
500 +        case 'e' : {
501 +
502 +           Vector3d eField;
503 +           eField[0] = tokenizer.nextTokenAsDouble();
504 +           eField[1] = tokenizer.nextTokenAsDouble();
505 +           eField[2] = tokenizer.nextTokenAsDouble();          
506 +           sd->setElectricField(eField);          
507 +           break;
508 +        }
509          default: {
510                 sprintf(painCave.errMsg,
511                         "DumpReader Error: %s is an unrecognized type\n", type.c_str());
# Line 431 | Line 520 | namespace oopse {
520    }
521    
522  
523 <  void  DumpReader::readStuntDoubles(std::istream& inputStream) {
523 >  void DumpReader::parseSiteLine(const std::string& line) {
524  
525 +    StringTokenizer tokenizer(line);
526 +    int nTokens;
527 +    
528 +    nTokens = tokenizer.countTokens();
529 +    
530 +    if (nTokens < 2) {  
531 +      sprintf(painCave.errMsg,
532 +              "DumpReader Error: Not enough Tokens.\n%s\n", line.c_str());
533 +      painCave.isFatal = 1;
534 +      simError();
535 +    }
536 +
537 +    /**
538 +     * The first token is the global integrable object index.
539 +     */
540 +
541 +    int index = tokenizer.nextTokenAsInt();
542 +    StuntDouble* sd = info_->getIOIndexToIntegrableObject(index);
543 +    if (sd == NULL) {
544 +      return;
545 +    }
546 +
547 +    /**
548 +     * Test to see if the next token is an integer or not.  If not,
549 +     * we've got data on the integrable object itself.  If there is an
550 +     * integer, we're parsing data for a site on a rigid body.
551 +     */
552 +
553 +    std::string indexTest = tokenizer.peekNextToken();
554 +    std::istringstream i(indexTest);
555 +    int siteIndex;
556 +    if (i >> siteIndex) {
557 +      // chew up this token and parse as an int:
558 +      siteIndex = tokenizer.nextTokenAsInt();
559 +      RigidBody* rb = static_cast<RigidBody*>(sd);
560 +      sd = rb->getAtoms()[siteIndex];
561 +    }
562 +
563 +    /**
564 +     * The next token contains information on what follows.
565 +     */
566 +    std::string type = tokenizer.nextToken();
567 +    int size = type.size();
568 +    
569 +    for(int i = 0; i < size; ++i) {
570 +      switch(type[i]) {
571 +        
572 +      case 'u' : {
573 +        
574 +        RealType particlePot;
575 +        particlePot = tokenizer.nextTokenAsDouble();
576 +        sd->setParticlePot(particlePot);
577 +        break;
578 +      }
579 +      case 'c' : {
580 +        
581 +        RealType flucQPos;
582 +        flucQPos = tokenizer.nextTokenAsDouble();
583 +        sd->setFlucQPos(flucQPos);
584 +        break;
585 +      }
586 +      case 'w' : {
587 +        
588 +        RealType flucQVel;
589 +        flucQVel = tokenizer.nextTokenAsDouble();
590 +        sd->setFlucQVel(flucQVel);
591 +        break;
592 +      }
593 +      case 'g' : {
594 +        
595 +        RealType flucQFrc;
596 +        flucQFrc = tokenizer.nextTokenAsDouble();
597 +        sd->setFlucQFrc(flucQFrc);
598 +        break;
599 +      }
600 +      case 'e' : {
601 +        
602 +        Vector3d eField;
603 +        eField[0] = tokenizer.nextTokenAsDouble();
604 +        eField[1] = tokenizer.nextTokenAsDouble();
605 +        eField[2] = tokenizer.nextTokenAsDouble();  
606 +        sd->setElectricField(eField);          
607 +        break;
608 +      }
609 +      default: {
610 +        sprintf(painCave.errMsg,
611 +                "DumpReader Error: %s is an unrecognized type\n", type.c_str());
612 +        painCave.isFatal = 1;
613 +        simError();
614 +        break;  
615 +      }
616 +      }
617 +    }    
618 +  }
619 +  
620 +  
621 +  void  DumpReader::readStuntDoubles(std::istream& inputStream) {
622 +    
623      inputStream.getline(buffer, bufferSize);
624      std::string line(buffer);
625      
# Line 455 | Line 642 | namespace oopse {
642    
643    }
644  
645 +  void  DumpReader::readSiteData(std::istream& inputStream) {
646 +
647 +    inputStream.getline(buffer, bufferSize);
648 +    std::string line(buffer);
649 +    
650 +    if (line.find("<SiteData>") == std::string::npos) {
651 +      // site data isn't required for a simulation, so skip
652 +      return;
653 +    }
654 +
655 +    while(inputStream.getline(buffer, bufferSize)) {
656 +      line = buffer;
657 +      
658 +      if(line.find("</SiteData>") != std::string::npos) {
659 +        break;
660 +      }
661 +
662 +      parseSiteLine(line);
663 +    }
664 +  
665 +  }
666 +
667    void DumpReader::readFrameProperties(std::istream& inputStream) {
668  
669      Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot();
# Line 500 | Line 709 | namespace oopse {
709          hmat(2, 2) = tokenizer.nextTokenAsDouble();
710          s->setHmat(hmat);      
711        } else if (propertyName == "Thermostat") {
712 <        RealType chi = tokenizer.nextTokenAsDouble();
713 <        RealType integralOfChiDt = tokenizer.nextTokenAsDouble();
714 <        s->setChi(chi);
715 <        s->setIntegralOfChiDt(integralOfChiDt);        
712 >        pair<RealType, RealType> thermostat;
713 >        thermostat.first = tokenizer.nextTokenAsDouble();
714 >        thermostat.second = tokenizer.nextTokenAsDouble();
715 >        s->setThermostat(thermostat);
716       } else if (propertyName == "Barostat") {
717          Mat3x3d eta;
718          eta(0, 0) = tokenizer.nextTokenAsDouble();
# Line 515 | Line 724 | namespace oopse {
724          eta(2, 0) = tokenizer.nextTokenAsDouble();
725          eta(2, 1) = tokenizer.nextTokenAsDouble();
726          eta(2, 2) = tokenizer.nextTokenAsDouble();
727 <        s->setEta(eta);
727 >        s->setBarostat(eta);
728        } else {
729          sprintf(painCave.errMsg,
730                  "DumpReader Error: %s is an invalid property in <FrameData>\n", propertyName.c_str());
# Line 528 | Line 737 | namespace oopse {
737    }
738  
739    
740 < }//end namespace oopse
740 > }//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 1790 by gezelter, Thu Aug 30 17:18:22 2012 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines