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 1880 by gezelter, Mon Jun 17 18:28:30 2013 UTC vs.
Revision 2020 by gezelter, Mon Sep 22 19:18:35 2014 UTC

# Line 39 | Line 39
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  
66    DumpWriter::DumpWriter(SimInfo* info)
67 <    : info_(info), filename_(info->getDumpFileName()), eorFilename_(info->getFinalConfigFileName()){
67 >    : info_(info), filename_(info->getDumpFileName()),
68 >      eorFilename_(info->getFinalConfigFileName()){
69  
70      Globals* simParams = info->getSimParams();
71      needCompression_   = simParams->getCompressDumpFile();
# Line 70 | Line 73 | namespace OpenMD {
73      needParticlePot_   = simParams->getOutputParticlePotential();
74      needFlucQ_         = simParams->getOutputFluctuatingCharges();
75      needElectricField_ = simParams->getOutputElectricField();
76 +    needSitePotential_ = simParams->getOutputSitePotential();
77  
78 <    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
78 >    if (needParticlePot_ || needFlucQ_ || needElectricField_ ||
79 >        needSitePotential_) {
80        doSiteData_ = true;
81      } else {
82        doSiteData_ = false;
# Line 119 | Line 124 | namespace OpenMD {
124      needParticlePot_   = simParams->getOutputParticlePotential();
125      needFlucQ_         = simParams->getOutputFluctuatingCharges();
126      needElectricField_ = simParams->getOutputElectricField();
127 +    needSitePotential_ = simParams->getOutputSitePotential();
128  
129 <    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
129 >    if (needParticlePot_ || needFlucQ_ || needElectricField_ ||
130 >        needSitePotential_) {
131        doSiteData_ = true;
132      } else {
133        doSiteData_ = false;
# Line 168 | Line 175 | namespace OpenMD {
175      needParticlePot_   = simParams->getOutputParticlePotential();
176      needFlucQ_         = simParams->getOutputFluctuatingCharges();
177      needElectricField_ = simParams->getOutputElectricField();
178 +    needSitePotential_ = simParams->getOutputSitePotential();
179  
180 <    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
180 >    if (needParticlePot_ || needFlucQ_ || needElectricField_ ||
181 >        needSitePotential_) {
182        doSiteData_ = true;
183      } else {
184        doSiteData_ = false;
# Line 302 | Line 311 | namespace OpenMD {
311    void DumpWriter::writeFrame(std::ostream& os) {
312  
313   #ifdef IS_MPI
314 <    MPI::Status istatus;
314 >    MPI_Status istatus;
315   #endif
316  
317      Molecule* mol;
# Line 360 | Line 369 | namespace OpenMD {
369   #else
370  
371      const int masterNode = 0;
372 <    int worldRank = MPI::COMM_WORLD.Get_rank();
373 <    int nProc = MPI::COMM_WORLD.Get_size();
372 >    int worldRank;
373 >    int nProc;
374  
375 +    MPI_Comm_size( MPI_COMM_WORLD, &nProc);
376 +    MPI_Comm_rank( MPI_COMM_WORLD, &worldRank);
377 +
378 +
379      if (worldRank == masterNode) {      
380        os << "  <Snapshot>\n";  
381        writeFrameProperties(os,
# Line 385 | Line 398 | namespace OpenMD {
398        
399        for (int i = 1; i < nProc; ++i) {
400          // tell processor i to start sending us data:
388        MPI::COMM_WORLD.Bcast(&i, 1, MPI::INT, masterNode);
401  
402 +        MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
403 +
404          // receive the length of the string buffer that was
405          // prepared by processor i:        
406          int recvLength;
407 <        MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, MPI::ANY_TAG,
408 <                             istatus);
407 >        MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD,
408 >                 &istatus);
409  
410          // create a buffer to receive the data
411          char* recvBuffer = new char[recvLength];
412          if (recvBuffer == NULL) {
413          } else {
414            // receive the data:
415 <          MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i,
416 <                               MPI::ANY_TAG, istatus);
415 >          MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i,
416 >                               MPI_ANY_TAG, MPI_COMM_WORLD, &istatus);
417            // send it to the file:
418            os << recvBuffer;
419            // get rid of the receive buffer:
# Line 411 | Line 425 | namespace OpenMD {
425        int myturn = 0;
426        for (int i = 1; i < nProc; ++i){
427          // wait for the master node to call our number:
428 <        MPI::COMM_WORLD.Bcast(&myturn, 1, MPI::INT, masterNode);
428 >        MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
429          if (myturn == worldRank){
430            // send the length of our buffer:
417          MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0);
431  
432 +          MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
433 +
434            // send our buffer:
435 <          MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength,
436 <                               MPI::CHAR, masterNode, 0);
435 >          MPI_Send((void *)buffer.c_str(), sendBufferLength,
436 >                   MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
437  
438          }
439        }
# Line 462 | Line 477 | namespace OpenMD {
477          for (int i = 1; i < nProc; ++i) {
478            
479            // tell processor i to start sending us data:
480 <          MPI::COMM_WORLD.Bcast(&i, 1, MPI::INT, masterNode);
480 >          MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
481            
482            // receive the length of the string buffer that was
483            // prepared by processor i:        
484            int recvLength;
485 <          MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, MPI::ANY_TAG,
486 <                               istatus);
485 >          MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD,
486 >                   &istatus);
487            
488            // create a buffer to receive the data
489            char* recvBuffer = new char[recvLength];
490            if (recvBuffer == NULL) {
491            } else {
492              // receive the data:
493 <            MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i,
494 <                                 MPI::ANY_TAG, istatus);
493 >            MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i,
494 >                     MPI_ANY_TAG, MPI_COMM_WORLD, &istatus);
495              // send it to the file:
496              os << recvBuffer;
497              // get rid of the receive buffer:
# Line 488 | Line 503 | namespace OpenMD {
503          int myturn = 0;
504          for (int i = 1; i < nProc; ++i){
505            // wait for the master node to call our number:
506 <          MPI::COMM_WORLD.Bcast(&myturn, 1, MPI::INT, masterNode);
506 >          MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
507            if (myturn == worldRank){
508              // send the length of our buffer:
509 <            MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0);
509 >            MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
510              // send our buffer:
511 <            MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength,
512 <                                 MPI::CHAR, masterNode, 0);
511 >            MPI_Send((void *)buffer.c_str(), sendBufferLength,
512 >                     MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
513            }
514          }
515        }
# Line 625 | Line 640 | namespace OpenMD {
640    }
641  
642    std::string DumpWriter::prepareSiteLine(StuntDouble* sd, int ioIndex, int siteIndex) {
643 <        
643 >    int storageLayout = info_->getSnapshotManager()->getStorageLayout();
644  
645      std::string id;
646      std::string type;
# Line 641 | Line 656 | namespace OpenMD {
656      }
657                
658      if (needFlucQ_) {
659 <      type += "cw";
660 <      RealType fqPos = sd->getFlucQPos();
661 <      if (isinf(fqPos) || isnan(fqPos) ) {      
662 <        sprintf( painCave.errMsg,
663 <                 "DumpWriter detected a numerical error writing the"
664 <                 " fluctuating charge for object %s", id.c_str());      
665 <        painCave.isFatal = 1;
666 <        simError();
667 <      }
668 <      sprintf(tempBuffer, " %13e ", fqPos);
669 <      line += tempBuffer;
670 <    
671 <      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;
659 >      if (storageLayout & DataStorage::dslFlucQPosition) {
660 >        type += "c";
661 >        RealType fqPos = sd->getFlucQPos();
662 >        if (isinf(fqPos) || isnan(fqPos) ) {      
663 >          sprintf( painCave.errMsg,
664 >                   "DumpWriter detected a numerical error writing the"
665 >                   " fluctuating charge for object %s", id.c_str());      
666 >          painCave.isFatal = 1;
667 >          simError();
668 >        }
669 >        sprintf(tempBuffer, " %13e ", fqPos);
670 >        line += tempBuffer;
671 >      }
672  
673 <      if (needForceVector_) {
674 <        type += "g";
675 <        RealType fqFrc = sd->getFlucQFrc();        
676 <        if (isinf(fqFrc) || isnan(fqFrc) ) {      
673 >      if (storageLayout & DataStorage::dslFlucQVelocity) {
674 >        type += "w";    
675 >        RealType fqVel = sd->getFlucQVel();
676 >        if (isinf(fqVel) || isnan(fqVel) ) {      
677            sprintf( painCave.errMsg,
678                     "DumpWriter detected a numerical error writing the"
679 <                   " fluctuating charge force for object %s", id.c_str());      
679 >                   " fluctuating charge velocity for object %s", id.c_str());      
680            painCave.isFatal = 1;
681            simError();
682          }
683 <        sprintf(tempBuffer, " %13e ", fqFrc);        
683 >        sprintf(tempBuffer, " %13e ", fqVel);
684          line += tempBuffer;
685        }
680    }
686  
687 +      if (needForceVector_) {
688 +        if (storageLayout & DataStorage::dslFlucQForce) {          
689 +          type += "g";
690 +          RealType fqFrc = sd->getFlucQFrc();        
691 +          if (isinf(fqFrc) || isnan(fqFrc) ) {      
692 +            sprintf( painCave.errMsg,
693 +                     "DumpWriter detected a numerical error writing the"
694 +                     " fluctuating charge force for object %s", id.c_str());      
695 +            painCave.isFatal = 1;
696 +            simError();
697 +          }
698 +          sprintf(tempBuffer, " %13e ", fqFrc);        
699 +          line += tempBuffer;
700 +        }
701 +      }
702 +    }
703 +    
704      if (needElectricField_) {
705 <      type += "e";
706 <      Vector3d eField= sd->getElectricField();
707 <      if (isinf(eField[0]) || isnan(eField[0]) ||
708 <          isinf(eField[1]) || isnan(eField[1]) ||
709 <          isinf(eField[2]) || isnan(eField[2]) ) {      
710 <        sprintf( painCave.errMsg,
711 <                 "DumpWriter detected a numerical error writing the electric"
712 <                 " field for object %s", id.c_str());      
713 <        painCave.isFatal = 1;
714 <        simError();
705 >      if (storageLayout & DataStorage::dslElectricField) {
706 >        type += "e";
707 >        Vector3d eField= sd->getElectricField();
708 >        if (isinf(eField[0]) || isnan(eField[0]) ||
709 >            isinf(eField[1]) || isnan(eField[1]) ||
710 >            isinf(eField[2]) || isnan(eField[2]) ) {      
711 >          sprintf( painCave.errMsg,
712 >                   "DumpWriter detected a numerical error writing the electric"
713 >                   " field for object %s", id.c_str());      
714 >          painCave.isFatal = 1;
715 >          simError();
716 >        }
717 >        sprintf(tempBuffer, " %13e %13e %13e",
718 >                eField[0], eField[1], eField[2]);
719 >        line += tempBuffer;
720        }
694      sprintf(tempBuffer, " %13e %13e %13e",
695              eField[0], eField[1], eField[2]);
696      line += tempBuffer;
721      }
722  
723 <
723 >    if (needSitePotential_) {
724 >      if (storageLayout & DataStorage::dslSitePotential) {          
725 >        type += "s";
726 >        RealType sPot = sd->getSitePotential();        
727 >        if (isinf(sPot) || isnan(sPot) ) {      
728 >          sprintf( painCave.errMsg,
729 >                   "DumpWriter detected a numerical error writing the"
730 >                   " site potential for object %s", id.c_str());      
731 >          painCave.isFatal = 1;
732 >          simError();
733 >        }
734 >        sprintf(tempBuffer, " %13e ", sPot);        
735 >        line += tempBuffer;
736 >      }
737 >    }    
738 >    
739      if (needParticlePot_) {
740 <      type += "u";
741 <      RealType particlePot = sd->getParticlePot();
742 <      if (isinf(particlePot) || isnan(particlePot)) {      
743 <        sprintf( painCave.errMsg,
744 <                 "DumpWriter detected a numerical error writing the particle "
745 <                 " potential for object %s", id.c_str());      
746 <        painCave.isFatal = 1;
747 <        simError();
740 >      if (storageLayout & DataStorage::dslParticlePot) {
741 >        type += "u";
742 >        RealType particlePot = sd->getParticlePot();
743 >        if (isinf(particlePot) || isnan(particlePot)) {      
744 >          sprintf( painCave.errMsg,
745 >                   "DumpWriter detected a numerical error writing the particle "
746 >                   " potential for object %s", id.c_str());      
747 >          painCave.isFatal = 1;
748 >          simError();
749 >        }
750 >        sprintf(tempBuffer, " %13e", particlePot);
751 >        line += tempBuffer;
752        }
710      sprintf(tempBuffer, " %13e", particlePot);
711      line += tempBuffer;
753      }
754 <    
714 <
754 >  
755      sprintf(tempBuffer, "%s %7s %s\n", id.c_str(), type.c_str(), line.c_str());
756      return std::string(tempBuffer);
757    }
# Line 721 | Line 761 | namespace OpenMD {
761    }
762  
763    void DumpWriter::writeEor() {
764 <    std::ostream* eorStream;
765 <    
764 >
765 >    std::ostream* eorStream = NULL;
766 >
767   #ifdef IS_MPI
768      if (worldRank == 0) {
769   #endif // is_mpi
# Line 755 | Line 796 | namespace OpenMD {
796   #ifdef IS_MPI
797      if (worldRank == 0) {
798   #endif // is_mpi
758
799        buffers.push_back(dumpFile_->rdbuf());
760
800        eorStream = createOStream(eorFilename_);
762
801        buffers.push_back(eorStream->rdbuf());
764        
802   #ifdef IS_MPI
803      }
804   #endif // is_mpi    
805  
806      TeeBuf tbuf(buffers.begin(), buffers.end());
807      std::ostream os(&tbuf);
771
808      writeFrame(os);
809  
810   #ifdef IS_MPI
# Line 778 | Line 814 | namespace OpenMD {
814        delete eorStream;
815   #ifdef IS_MPI
816      }
817 < #endif // is_mpi  
782 <    
817 > #endif // is_mpi      
818    }
819  
820    std::ostream* DumpWriter::createOStream(const std::string& filename) {

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines