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

Comparing trunk/src/io/DumpWriter.cpp (file contents):
Revision 1796 by gezelter, Mon Sep 10 18:38:44 2012 UTC vs.
Revision 1971 by gezelter, Fri Feb 28 13:25:13 2014 UTC

# Line 35 | Line 35
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).          
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 + #include "config.h"
44 +
45 + #ifdef IS_MPI
46 + #include <mpi.h>
47 + #endif
48  
49   #include "io/DumpWriter.hpp"
50   #include "primitives/Molecule.hpp"
# Line 54 | Line 60
60   #define isinf(x) (!_finite(x) && !_isnan(x))
61   #endif
62  
57 – #ifdef IS_MPI
58 – #include <mpi.h>
59 – #endif
60 –
63   using namespace std;
64   namespace OpenMD {
65  
# Line 302 | Line 304 | namespace OpenMD {
304    void DumpWriter::writeFrame(std::ostream& os) {
305  
306   #ifdef IS_MPI
307 <    MPI::Status istatus;
307 >    MPI_Status istatus;
308   #endif
309  
310      Molecule* mol;
# Line 310 | Line 312 | namespace OpenMD {
312      SimInfo::MoleculeIterator mi;
313      Molecule::IntegrableObjectIterator ii;
314      RigidBody::AtomIterator ai;
313 –    Atom* atom;
315  
316   #ifndef IS_MPI
317      os << "  <Snapshot>\n";
# Line 318 | Line 319 | namespace OpenMD {
319      writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot());
320  
321      os << "    <StuntDoubles>\n";
322 <    for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
323 <
322 >    for (mol = info_->beginMolecule(mi); mol != NULL;
323 >         mol = info_->nextMolecule(mi)) {
324        
325        for (sd = mol->beginIntegrableObject(ii); sd != NULL;  
326             sd = mol->nextIntegrableObject(ii)) {        
# Line 335 | Line 336 | namespace OpenMD {
336             mol = info_->nextMolecule(mi)) {
337                
338          for (sd = mol->beginIntegrableObject(ii); sd != NULL;  
339 <             sd = mol->nextIntegrableObject(ii)) {      
340 <          
339 >           sd = mol->nextIntegrableObject(ii)) {        
340 >
341            int ioIndex = sd->getGlobalIntegrableObjectIndex();
342            // do one for the IO itself
343            os << prepareSiteLine(sd, ioIndex, 0);
# Line 345 | Line 346 | namespace OpenMD {
346              
347              RigidBody* rb = static_cast<RigidBody*>(sd);
348              int siteIndex = 0;
349 <            for (atom = rb->beginAtom(ai); atom != NULL;  
349 >            for (Atom* atom = rb->beginAtom(ai); atom != NULL;  
350                   atom = rb->nextAtom(ai)) {                                            
351                os << prepareSiteLine(atom, ioIndex, siteIndex);
352                siteIndex++;
# Line 361 | Line 362 | namespace OpenMD {
362   #else
363  
364      const int masterNode = 0;
365 <    int worldRank = MPI::COMM_WORLD.Get_rank();
366 <    int nProc = MPI::COMM_WORLD.Get_size();
365 >    int worldRank;
366 >    int nProc;
367 >
368 >    MPI_Comm_size( MPI_COMM_WORLD, &nProc);
369 >    MPI_Comm_rank( MPI_COMM_WORLD, &worldRank);
370  
371 +
372      if (worldRank == masterNode) {      
373        os << "  <Snapshot>\n";  
374        writeFrameProperties(os,
# Line 386 | Line 391 | namespace OpenMD {
391        
392        for (int i = 1; i < nProc; ++i) {
393          // tell processor i to start sending us data:
394 <        MPI::COMM_WORLD.Bcast(&i, 1, MPI::INT, masterNode);
394 >        MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
395  
396          // receive the length of the string buffer that was
397          // prepared by processor i:        
398          int recvLength;
399 <        MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, MPI::ANY_TAG,
400 <                             istatus);
399 >        MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD,
400 >                 &istatus);
401  
402          // create a buffer to receive the data
403          char* recvBuffer = new char[recvLength];
404          if (recvBuffer == NULL) {
405          } else {
406            // receive the data:
407 <          MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i,
408 <                               MPI::ANY_TAG, istatus);
407 >          MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i,
408 >                               MPI_ANY_TAG, MPI_COMM_WORLD, &istatus);
409            // send it to the file:
410            os << recvBuffer;
411            // get rid of the receive buffer:
# Line 412 | Line 417 | namespace OpenMD {
417        int myturn = 0;
418        for (int i = 1; i < nProc; ++i){
419          // wait for the master node to call our number:
420 <        MPI::COMM_WORLD.Bcast(&myturn, 1, MPI::INT, masterNode);
420 >        MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
421          if (myturn == worldRank){
422            // send the length of our buffer:
423 <          MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0);
423 >          MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
424  
425            // send our buffer:
426 <          MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength,
427 <                               MPI::CHAR, masterNode, 0);
426 >          MPI_Send((void *)buffer.c_str(), sendBufferLength,
427 >                   MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
428 >
429          }
430        }
431      }
# Line 447 | Line 453 | namespace OpenMD {
453              
454              RigidBody* rb = static_cast<RigidBody*>(sd);
455              int siteIndex = 0;
456 <            for (atom = rb->beginAtom(ai); atom != NULL;  
456 >            for (Atom* atom = rb->beginAtom(ai); atom != NULL;  
457                   atom = rb->nextAtom(ai)) {                                            
458                buffer += prepareSiteLine(atom, ioIndex, siteIndex);
459                siteIndex++;
# Line 462 | Line 468 | namespace OpenMD {
468          for (int i = 1; i < nProc; ++i) {
469            
470            // tell processor i to start sending us data:
471 <          MPI::COMM_WORLD.Bcast(&i, 1, MPI::INT, masterNode);
471 >          MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
472            
473            // receive the length of the string buffer that was
474            // prepared by processor i:        
475            int recvLength;
476 <          MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, MPI::ANY_TAG,
477 <                               istatus);
476 >          MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD,
477 >                   &istatus);
478            
479            // create a buffer to receive the data
480            char* recvBuffer = new char[recvLength];
481            if (recvBuffer == NULL) {
482            } else {
483              // receive the data:
484 <            MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i,
485 <                                 MPI::ANY_TAG, istatus);
484 >            MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i,
485 >                     MPI_ANY_TAG, MPI_COMM_WORLD, &istatus);
486              // send it to the file:
487              os << recvBuffer;
488              // get rid of the receive buffer:
# Line 488 | Line 494 | namespace OpenMD {
494          int myturn = 0;
495          for (int i = 1; i < nProc; ++i){
496            // wait for the master node to call our number:
497 <          MPI::COMM_WORLD.Bcast(&myturn, 1, MPI::INT, masterNode);
497 >          MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
498            if (myturn == worldRank){
499              // send the length of our buffer:
500 <            MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0);
500 >            MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
501              // send our buffer:
502 <            MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength,
503 <                                 MPI::CHAR, masterNode, 0);
502 >            MPI_Send((void *)buffer.c_str(), sendBufferLength,
503 >                     MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
504            }
505          }
506        }
# Line 625 | Line 631 | namespace OpenMD {
631    }
632  
633    std::string DumpWriter::prepareSiteLine(StuntDouble* sd, int ioIndex, int siteIndex) {
634 <        
634 >    int storageLayout = info_->getSnapshotManager()->getStorageLayout();
635  
636      std::string id;
637      std::string type;
# Line 641 | Line 647 | namespace OpenMD {
647      }
648                
649      if (needFlucQ_) {
650 <      type += "cw";
651 <      RealType fqPos = sd->getFlucQPos();
652 <      if (isinf(fqPos) || isnan(fqPos) ) {      
653 <        sprintf( painCave.errMsg,
654 <                 "DumpWriter detected a numerical error writing the"
655 <                 " fluctuating charge for object %s", id.c_str());      
656 <        painCave.isFatal = 1;
657 <        simError();
658 <      }
659 <      sprintf(tempBuffer, " %13e ", fqPos);
660 <      line += tempBuffer;
661 <    
662 <      RealType fqVel = sd->getFlucQVel();
657 <      if (isinf(fqVel) || isnan(fqVel) ) {      
658 <        sprintf( painCave.errMsg,
659 <                 "DumpWriter detected a numerical error writing the"
660 <                 " fluctuating charge velocity for object %s", id.c_str());      
661 <        painCave.isFatal = 1;
662 <        simError();
663 <      }
664 <      sprintf(tempBuffer, " %13e ", fqVel);
665 <      line += tempBuffer;
650 >      if (storageLayout & DataStorage::dslFlucQPosition) {
651 >        type += "c";
652 >        RealType fqPos = sd->getFlucQPos();
653 >        if (isinf(fqPos) || isnan(fqPos) ) {      
654 >          sprintf( painCave.errMsg,
655 >                   "DumpWriter detected a numerical error writing the"
656 >                   " fluctuating charge for object %s", id.c_str());      
657 >          painCave.isFatal = 1;
658 >          simError();
659 >        }
660 >        sprintf(tempBuffer, " %13e ", fqPos);
661 >        line += tempBuffer;
662 >      }
663  
664 <      if (needForceVector_) {
665 <        type += "g";
666 <        RealType fqFrc = sd->getFlucQFrc();        
667 <        if (isinf(fqFrc) || isnan(fqFrc) ) {      
664 >      if (storageLayout & DataStorage::dslFlucQVelocity) {
665 >        type += "w";    
666 >        RealType fqVel = sd->getFlucQVel();
667 >        if (isinf(fqVel) || isnan(fqVel) ) {      
668            sprintf( painCave.errMsg,
669                     "DumpWriter detected a numerical error writing the"
670 <                   " fluctuating charge force for object %s", id.c_str());      
670 >                   " fluctuating charge velocity for object %s", id.c_str());      
671            painCave.isFatal = 1;
672            simError();
673          }
674 <        sprintf(tempBuffer, " %13e ", fqFrc);        
674 >        sprintf(tempBuffer, " %13e ", fqVel);
675          line += tempBuffer;
676        }
680 –    }
677  
678 +      if (needForceVector_) {
679 +        if (storageLayout & DataStorage::dslFlucQForce) {          
680 +          type += "g";
681 +          RealType fqFrc = sd->getFlucQFrc();        
682 +          if (isinf(fqFrc) || isnan(fqFrc) ) {      
683 +            sprintf( painCave.errMsg,
684 +                     "DumpWriter detected a numerical error writing the"
685 +                     " fluctuating charge force for object %s", id.c_str());      
686 +            painCave.isFatal = 1;
687 +            simError();
688 +          }
689 +          sprintf(tempBuffer, " %13e ", fqFrc);        
690 +          line += tempBuffer;
691 +        }
692 +      }
693 +    }
694 +    
695      if (needElectricField_) {
696 <      type += "e";
697 <      Vector3d eField= sd->getElectricField();
698 <      if (isinf(eField[0]) || isnan(eField[0]) ||
699 <          isinf(eField[1]) || isnan(eField[1]) ||
700 <          isinf(eField[2]) || isnan(eField[2]) ) {      
701 <        sprintf( painCave.errMsg,
702 <                 "DumpWriter detected a numerical error writing the electric"
703 <                 " field for object %s", id.c_str());      
704 <        painCave.isFatal = 1;
705 <        simError();
696 >      if (storageLayout & DataStorage::dslElectricField) {
697 >        type += "e";
698 >        Vector3d eField= sd->getElectricField();
699 >        if (isinf(eField[0]) || isnan(eField[0]) ||
700 >            isinf(eField[1]) || isnan(eField[1]) ||
701 >            isinf(eField[2]) || isnan(eField[2]) ) {      
702 >          sprintf( painCave.errMsg,
703 >                   "DumpWriter detected a numerical error writing the electric"
704 >                   " field for object %s", id.c_str());      
705 >          painCave.isFatal = 1;
706 >          simError();
707 >        }
708 >        sprintf(tempBuffer, " %13e %13e %13e",
709 >                eField[0], eField[1], eField[2]);
710 >        line += tempBuffer;
711        }
694 –      sprintf(tempBuffer, " %13e %13e %13e",
695 –              eField[0], eField[1], eField[2]);
696 –      line += tempBuffer;
712      }
713  
714  
715      if (needParticlePot_) {
716 <      type += "u";
717 <      RealType particlePot = sd->getParticlePot();
718 <      if (isinf(particlePot) || isnan(particlePot)) {      
719 <        sprintf( painCave.errMsg,
720 <                 "DumpWriter detected a numerical error writing the particle "
721 <                 " potential for object %s", id.c_str());      
722 <        painCave.isFatal = 1;
723 <        simError();
716 >      if (storageLayout & DataStorage::dslParticlePot) {
717 >        type += "u";
718 >        RealType particlePot = sd->getParticlePot();
719 >        if (isinf(particlePot) || isnan(particlePot)) {      
720 >          sprintf( painCave.errMsg,
721 >                   "DumpWriter detected a numerical error writing the particle "
722 >                   " potential for object %s", id.c_str());      
723 >          painCave.isFatal = 1;
724 >          simError();
725 >        }
726 >        sprintf(tempBuffer, " %13e", particlePot);
727 >        line += tempBuffer;
728        }
710 –      sprintf(tempBuffer, " %13e", particlePot);
711 –      line += tempBuffer;
729      }
730 <    
714 <
730 >  
731      sprintf(tempBuffer, "%s %7s %s\n", id.c_str(), type.c_str(), line.c_str());
732      return std::string(tempBuffer);
733    }
# Line 726 | Line 742 | namespace OpenMD {
742   #ifdef IS_MPI
743      if (worldRank == 0) {
744   #endif // is_mpi
745 <
745 >      
746        eorStream = createOStream(eorFilename_);
747  
748   #ifdef IS_MPI
749      }
750 < #endif // is_mpi    
751 <
750 > #endif
751 >    
752      writeFrame(*eorStream);
753 <
753 >      
754   #ifdef IS_MPI
755      if (worldRank == 0) {
756 < #endif // is_mpi
756 > #endif
757 >      
758        writeClosing(*eorStream);
759        delete eorStream;
760 +      
761   #ifdef IS_MPI
762      }
763   #endif // is_mpi  
# Line 753 | Line 771 | namespace OpenMD {
771   #ifdef IS_MPI
772      if (worldRank == 0) {
773   #endif // is_mpi
756 –
774        buffers.push_back(dumpFile_->rdbuf());
758 –
775        eorStream = createOStream(eorFilename_);
760 –
776        buffers.push_back(eorStream->rdbuf());
762 –        
777   #ifdef IS_MPI
778      }
779   #endif // is_mpi    
780  
781      TeeBuf tbuf(buffers.begin(), buffers.end());
782      std::ostream os(&tbuf);
769 –
783      writeFrame(os);
784  
785   #ifdef IS_MPI
# Line 776 | Line 789 | namespace OpenMD {
789        delete eorStream;
790   #ifdef IS_MPI
791      }
792 < #endif // is_mpi  
780 <    
792 > #endif // is_mpi      
793    }
794  
795    std::ostream* DumpWriter::createOStream(const std::string& filename) {

Diff Legend

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