--- trunk/src/io/DumpWriter.cpp 2012/08/31 21:16:10 1793 +++ trunk/src/io/DumpWriter.cpp 2014/04/29 17:32:31 1993 @@ -35,10 +35,16 @@ * * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). - * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). + * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ + +#include "config.h" + +#ifdef IS_MPI +#include +#endif #include "io/DumpWriter.hpp" #include "primitives/Molecule.hpp" @@ -54,15 +60,12 @@ #define isinf(x) (!_finite(x) && !_isnan(x)) #endif -#ifdef IS_MPI -#include -#endif - using namespace std; namespace OpenMD { DumpWriter::DumpWriter(SimInfo* info) - : info_(info), filename_(info->getDumpFileName()), eorFilename_(info->getFinalConfigFileName()){ + : info_(info), filename_(info->getDumpFileName()), + eorFilename_(info->getFinalConfigFileName()){ Globals* simParams = info->getSimParams(); needCompression_ = simParams->getCompressDumpFile(); @@ -70,8 +73,10 @@ namespace OpenMD { needParticlePot_ = simParams->getOutputParticlePotential(); needFlucQ_ = simParams->getOutputFluctuatingCharges(); needElectricField_ = simParams->getOutputElectricField(); + needSitePotential_ = simParams->getOutputSitePotential(); - if (needParticlePot_ || needFlucQ_ || needElectricField_) { + if (needParticlePot_ || needFlucQ_ || needElectricField_ || + needSitePotential_) { doSiteData_ = true; } else { doSiteData_ = false; @@ -119,8 +124,10 @@ namespace OpenMD { needParticlePot_ = simParams->getOutputParticlePotential(); needFlucQ_ = simParams->getOutputFluctuatingCharges(); needElectricField_ = simParams->getOutputElectricField(); + needSitePotential_ = simParams->getOutputSitePotential(); - if (needParticlePot_ || needFlucQ_ || needElectricField_) { + if (needParticlePot_ || needFlucQ_ || needElectricField_ || + needSitePotential_) { doSiteData_ = true; } else { doSiteData_ = false; @@ -168,8 +175,10 @@ namespace OpenMD { needParticlePot_ = simParams->getOutputParticlePotential(); needFlucQ_ = simParams->getOutputFluctuatingCharges(); needElectricField_ = simParams->getOutputElectricField(); + needSitePotential_ = simParams->getOutputSitePotential(); - if (needParticlePot_ || needFlucQ_ || needElectricField_) { + if (needParticlePot_ || needFlucQ_ || needElectricField_ || + needSitePotential_) { doSiteData_ = true; } else { doSiteData_ = false; @@ -310,7 +319,6 @@ namespace OpenMD { SimInfo::MoleculeIterator mi; Molecule::IntegrableObjectIterator ii; RigidBody::AtomIterator ai; - Atom* atom; #ifndef IS_MPI os << " \n"; @@ -318,8 +326,8 @@ namespace OpenMD { writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot()); os << " \n"; - for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) { - + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { for (sd = mol->beginIntegrableObject(ii); sd != NULL; sd = mol->nextIntegrableObject(ii)) { @@ -331,11 +339,12 @@ namespace OpenMD { if (doSiteData_) { os << " \n"; - for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) { + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { for (sd = mol->beginIntegrableObject(ii); sd != NULL; - sd = mol->nextIntegrableObject(ii)) { - + sd = mol->nextIntegrableObject(ii)) { + int ioIndex = sd->getGlobalIntegrableObjectIndex(); // do one for the IO itself os << prepareSiteLine(sd, ioIndex, 0); @@ -344,7 +353,7 @@ namespace OpenMD { RigidBody* rb = static_cast(sd); int siteIndex = 0; - for (atom = rb->beginAtom(ai); atom != NULL; + for (Atom* atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { os << prepareSiteLine(atom, ioIndex, siteIndex); siteIndex++; @@ -358,61 +367,163 @@ namespace OpenMD { os.flush(); #else - //every node prepares the dump lines for integrable objects belong to itself - std::string buffer; - for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) { + const int masterNode = 0; + int worldRank; + int nProc; + MPI_Comm_size( MPI_COMM_WORLD, &nProc); + MPI_Comm_rank( MPI_COMM_WORLD, &worldRank); + + + if (worldRank == masterNode) { + os << " \n"; + writeFrameProperties(os, + info_->getSnapshotManager()->getCurrentSnapshot()); + os << " \n"; + } + + //every node prepares the dump lines for integrable objects belong to itself + std::string buffer; + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { for (sd = mol->beginIntegrableObject(ii); sd != NULL; sd = mol->nextIntegrableObject(ii)) { - buffer += prepareDumpLine(sd); + buffer += prepareDumpLine(sd); } } - const int masterNode = 0; - int nProc; - MPI_Comm_size(MPI_COMM_WORLD, &nProc); if (worldRank == masterNode) { - os << " \n"; - writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot()); - os << " \n"; - os << buffer; - + for (int i = 1; i < nProc; ++i) { + // tell processor i to start sending us data: + MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD); // receive the length of the string buffer that was - // prepared by processor i - - MPI_Bcast(&i, 1, MPI_INT,masterNode,MPI_COMM_WORLD); + // prepared by processor i: int recvLength; - MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, &istatus); + MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD, + &istatus); + + // create a buffer to receive the data char* recvBuffer = new char[recvLength]; if (recvBuffer == NULL) { } else { - MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, &istatus); + // receive the data: + MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, + MPI_ANY_TAG, MPI_COMM_WORLD, &istatus); + // send it to the file: os << recvBuffer; + // get rid of the receive buffer: delete [] recvBuffer; } } - os << " \n"; - - os << " \n"; - os.flush(); } else { int sendBufferLength = buffer.size() + 1; int myturn = 0; for (int i = 1; i < nProc; ++i){ - MPI_Bcast(&myturn,1, MPI_INT,masterNode,MPI_COMM_WORLD); + // wait for the master node to call our number: + MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD); if (myturn == worldRank){ + // send the length of our buffer: MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD); - MPI_Send((void *)buffer.c_str(), sendBufferLength, MPI_CHAR, masterNode, 0, MPI_COMM_WORLD); + + // send our buffer: + MPI_Send((void *)buffer.c_str(), sendBufferLength, + MPI_CHAR, masterNode, 0, MPI_COMM_WORLD); + } } } + + if (worldRank == masterNode) { + os << " \n"; + } -#endif // is_mpi + if (doSiteData_) { + if (worldRank == masterNode) { + os << " \n"; + } + buffer.clear(); + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { + + for (sd = mol->beginIntegrableObject(ii); sd != NULL; + sd = mol->nextIntegrableObject(ii)) { + + int ioIndex = sd->getGlobalIntegrableObjectIndex(); + // do one for the IO itself + buffer += prepareSiteLine(sd, ioIndex, 0); + if (sd->isRigidBody()) { + + RigidBody* rb = static_cast(sd); + int siteIndex = 0; + for (Atom* atom = rb->beginAtom(ai); atom != NULL; + atom = rb->nextAtom(ai)) { + buffer += prepareSiteLine(atom, ioIndex, siteIndex); + siteIndex++; + } + } + } + } + + if (worldRank == masterNode) { + os << buffer; + + for (int i = 1; i < nProc; ++i) { + + // tell processor i to start sending us data: + MPI_Bcast(&i, 1, MPI_INT, masterNode, MPI_COMM_WORLD); + + // receive the length of the string buffer that was + // prepared by processor i: + int recvLength; + MPI_Recv(&recvLength, 1, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD, + &istatus); + + // create a buffer to receive the data + char* recvBuffer = new char[recvLength]; + if (recvBuffer == NULL) { + } else { + // receive the data: + MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, + MPI_ANY_TAG, MPI_COMM_WORLD, &istatus); + // send it to the file: + os << recvBuffer; + // get rid of the receive buffer: + delete [] recvBuffer; + } + } + } else { + int sendBufferLength = buffer.size() + 1; + int myturn = 0; + for (int i = 1; i < nProc; ++i){ + // wait for the master node to call our number: + MPI_Bcast(&myturn, 1, MPI_INT, masterNode, MPI_COMM_WORLD); + if (myturn == worldRank){ + // send the length of our buffer: + MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD); + // send our buffer: + MPI_Send((void *)buffer.c_str(), sendBufferLength, + MPI_CHAR, masterNode, 0, MPI_COMM_WORLD); + } + } + } + + if (worldRank == masterNode) { + os << " \n"; + } + } + + if (worldRank == masterNode) { + os << " \n"; + os.flush(); + } + +#endif // is_mpi + } std::string DumpWriter::prepareDumpLine(StuntDouble* sd) { @@ -527,7 +638,7 @@ namespace OpenMD { } std::string DumpWriter::prepareSiteLine(StuntDouble* sd, int ioIndex, int siteIndex) { - + int storageLayout = info_->getSnapshotManager()->getStorageLayout(); std::string id; std::string type; @@ -543,77 +654,102 @@ namespace OpenMD { } if (needFlucQ_) { - type += "cw"; - RealType fqPos = sd->getFlucQPos(); - if (isinf(fqPos) || isnan(fqPos) ) { - sprintf( painCave.errMsg, - "DumpWriter detected a numerical error writing the" - " fluctuating charge for object %s", id.c_str()); - painCave.isFatal = 1; - simError(); - } - sprintf(tempBuffer, " %13e ", fqPos); - line += tempBuffer; - - RealType fqVel = sd->getFlucQVel(); - if (isinf(fqVel) || isnan(fqVel) ) { - sprintf( painCave.errMsg, - "DumpWriter detected a numerical error writing the" - " fluctuating charge velocity for object %s", id.c_str()); - painCave.isFatal = 1; - simError(); - } - sprintf(tempBuffer, " %13e ", fqVel); - line += tempBuffer; + if (storageLayout & DataStorage::dslFlucQPosition) { + type += "c"; + RealType fqPos = sd->getFlucQPos(); + if (isinf(fqPos) || isnan(fqPos) ) { + sprintf( painCave.errMsg, + "DumpWriter detected a numerical error writing the" + " fluctuating charge for object %s", id.c_str()); + painCave.isFatal = 1; + simError(); + } + sprintf(tempBuffer, " %13e ", fqPos); + line += tempBuffer; + } - if (needForceVector_) { - type += "g"; - RealType fqFrc = sd->getFlucQFrc(); - if (isinf(fqFrc) || isnan(fqFrc) ) { + if (storageLayout & DataStorage::dslFlucQVelocity) { + type += "w"; + RealType fqVel = sd->getFlucQVel(); + if (isinf(fqVel) || isnan(fqVel) ) { sprintf( painCave.errMsg, "DumpWriter detected a numerical error writing the" - " fluctuating charge force for object %s", id.c_str()); + " fluctuating charge velocity for object %s", id.c_str()); painCave.isFatal = 1; simError(); } - sprintf(tempBuffer, " %13e ", fqFrc); + sprintf(tempBuffer, " %13e ", fqVel); line += tempBuffer; } - } + if (needForceVector_) { + if (storageLayout & DataStorage::dslFlucQForce) { + type += "g"; + RealType fqFrc = sd->getFlucQFrc(); + if (isinf(fqFrc) || isnan(fqFrc) ) { + sprintf( painCave.errMsg, + "DumpWriter detected a numerical error writing the" + " fluctuating charge force for object %s", id.c_str()); + painCave.isFatal = 1; + simError(); + } + sprintf(tempBuffer, " %13e ", fqFrc); + line += tempBuffer; + } + } + } + if (needElectricField_) { - type += "e"; - Vector3d eField= sd->getElectricField(); - if (isinf(eField[0]) || isnan(eField[0]) || - isinf(eField[1]) || isnan(eField[1]) || - isinf(eField[2]) || isnan(eField[2]) ) { - sprintf( painCave.errMsg, - "DumpWriter detected a numerical error writing the electric" - " field for object %s", id.c_str()); - painCave.isFatal = 1; - simError(); + if (storageLayout & DataStorage::dslElectricField) { + type += "e"; + Vector3d eField= sd->getElectricField(); + if (isinf(eField[0]) || isnan(eField[0]) || + isinf(eField[1]) || isnan(eField[1]) || + isinf(eField[2]) || isnan(eField[2]) ) { + sprintf( painCave.errMsg, + "DumpWriter detected a numerical error writing the electric" + " field for object %s", id.c_str()); + painCave.isFatal = 1; + simError(); + } + sprintf(tempBuffer, " %13e %13e %13e", + eField[0], eField[1], eField[2]); + line += tempBuffer; } - sprintf(tempBuffer, " %13e %13e %13e", - eField[0], eField[1], eField[2]); - line += tempBuffer; } - + if (needSitePotential_) { + if (storageLayout & DataStorage::dslSitePotential) { + type += "s"; + RealType sPot = sd->getSitePotential(); + if (isinf(sPot) || isnan(sPot) ) { + sprintf( painCave.errMsg, + "DumpWriter detected a numerical error writing the" + " site potential for object %s", id.c_str()); + painCave.isFatal = 1; + simError(); + } + sprintf(tempBuffer, " %13e ", sPot); + line += tempBuffer; + } + } + if (needParticlePot_) { - type += "u"; - RealType particlePot = sd->getParticlePot(); - if (isinf(particlePot) || isnan(particlePot)) { - sprintf( painCave.errMsg, - "DumpWriter detected a numerical error writing the particle " - " potential for object %s", id.c_str()); - painCave.isFatal = 1; - simError(); + if (storageLayout & DataStorage::dslParticlePot) { + type += "u"; + RealType particlePot = sd->getParticlePot(); + if (isinf(particlePot) || isnan(particlePot)) { + sprintf( painCave.errMsg, + "DumpWriter detected a numerical error writing the particle " + " potential for object %s", id.c_str()); + painCave.isFatal = 1; + simError(); + } + sprintf(tempBuffer, " %13e", particlePot); + line += tempBuffer; } - sprintf(tempBuffer, " %13e", particlePot); - line += tempBuffer; } - - + sprintf(tempBuffer, "%s %7s %s\n", id.c_str(), type.c_str(), line.c_str()); return std::string(tempBuffer); } @@ -623,25 +759,28 @@ namespace OpenMD { } void DumpWriter::writeEor() { - std::ostream* eorStream; - + + std::ostream* eorStream = NULL; + #ifdef IS_MPI if (worldRank == 0) { #endif // is_mpi - + eorStream = createOStream(eorFilename_); #ifdef IS_MPI } -#endif // is_mpi - +#endif + writeFrame(*eorStream); - + #ifdef IS_MPI if (worldRank == 0) { -#endif // is_mpi +#endif + writeClosing(*eorStream); delete eorStream; + #ifdef IS_MPI } #endif // is_mpi @@ -655,20 +794,15 @@ namespace OpenMD { #ifdef IS_MPI if (worldRank == 0) { #endif // is_mpi - buffers.push_back(dumpFile_->rdbuf()); - eorStream = createOStream(eorFilename_); - buffers.push_back(eorStream->rdbuf()); - #ifdef IS_MPI } #endif // is_mpi TeeBuf tbuf(buffers.begin(), buffers.end()); std::ostream os(&tbuf); - writeFrame(os); #ifdef IS_MPI @@ -678,8 +812,7 @@ namespace OpenMD { delete eorStream; #ifdef IS_MPI } -#endif // is_mpi - +#endif // is_mpi } std::ostream* DumpWriter::createOStream(const std::string& filename) {