--- trunk/src/io/RestWriter.cpp 2009/10/07 20:49:50 1364 +++ trunk/src/io/RestWriter.cpp 2014/09/23 15:25:08 2021 @@ -6,19 +6,10 @@ * redistribute this software in source and binary code form, provided * that the following conditions are met: * - * 1. Acknowledgement of the program authors must be made in any - * publication of scientific results based in part on use of the - * program. An acceptable form of acknowledgement is citation of - * the article in which the program was described (Matthew - * A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher - * J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented - * Parallel Simulation Engine for Molecular Dynamics," - * J. Comput. Chem. 26, pp. 252-271 (2005)) - * - * 2. Redistributions of source code must retain the above copyright + * 1. Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * - * 3. Redistributions in binary form must reproduce the above copyright + * 2. Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the * distribution. @@ -37,95 +28,273 @@ * arising out of the use of or inability to use software, even if the * University of Notre Dame has been advised of the possibility of * such damages. + * + * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your + * research, please cite the appropriate papers when you publish your + * work. Good starting points are: + * + * [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, 234107 (2008). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ +#ifdef IS_MPI +#include +#endif +#include +#include #include #include "io/RestWriter.hpp" #include "utils/simError.h" #include "brains/SnapshotManager.hpp" -#ifdef IS_MPI -#include -#endif -namespace oopse { +namespace OpenMD { RestWriter::RestWriter(SimInfo* info, const std::string& filename, std::vector restraints ) : info_(info){ - //use master - slave mode, only master node writes to disk + std::vector::const_iterator resti; + + createRestFile_ = false; + +#ifdef IS_MPI + MPI_Status* istatus; +#endif + + int printAny = 0; + for(resti=restraints.begin(); resti != restraints.end(); ++resti){ + if ((*resti)->getPrintRestraint()) { + printAny = 1; + } + } + #ifdef IS_MPI + MPI_Allreduce(MPI_IN_PLACE, &printAny, 1, MPI_INT, MPI_SUM, + MPI_COMM_WORLD); +#endif + + if (printAny) createRestFile_ = true; + +#ifdef IS_MPI if(worldRank == 0){ #endif - - output_.open(filename.c_str()); - - if(!output_){ - sprintf( painCave.errMsg, - "Could not open %s for restraint output.\n", - filename.c_str()); - painCave.isFatal = 1; - simError(); + + if (createRestFile_) { + output_ = new std::ofstream(filename.c_str()); + + if(!output_){ + sprintf( painCave.errMsg, + "Could not open %s for restraint output.\n", + filename.c_str()); + painCave.isFatal = 1; + simError(); + } } + +#ifdef IS_MPI + } +#endif // is_mpi - output_ << "#time\t"; - - // TODO: get Restraint info from slave nodes: - std::vector::const_iterator resti; - for(resti=restraints.begin(); resti != restraints.end(); ++resti){ - if ((*resti)->getPrintRestraint()) { - std::string myName = (*resti)->getRestraintName(); - int myType = (*resti)->getRestraintType(); +#ifndef IS_MPI + + if (createRestFile_) (*output_) << "#time\t"; + + for(resti=restraints.begin(); resti != restraints.end(); ++resti){ + if ((*resti)->getPrintRestraint()) { + std::string myName = (*resti)->getRestraintName(); + int myType = (*resti)->getRestraintType(); + + (*output_) << myName << ":"; + + if (myType & Restraint::rtDisplacement) + (*output_) << "\tPosition(angstroms)\tEnergy(kcal/mol)"; + + if (myType & Restraint::rtTwist) + (*output_) << "\tTwistAngle(radians)\tEnergy(kcal/mol)"; + + if (myType & Restraint::rtSwingX) + (*output_) << "\tSwingXAngle(radians)\tEnergy(kcal/mol)"; - output_ << myName << ":"; - - if (myType & Restraint::rtDisplacement) - output_ << "\tPosition(angstroms)\tEnergy(kcal/mol)"; - - if (myType & Restraint::rtTwist) - output_ << "\tTwistAngle(radians)\tEnergy(kcal/mol)"; - - if (myType & Restraint::rtSwingX) - output_ << "\tSwingXAngle(radians)\tEnergy(kcal/mol)"; - - if (myType & Restraint::rtSwingY) - output_ << "\tSwingYAngle(radians)\tEnergy(kcal/mol)"; - - } + if (myType & Restraint::rtSwingY) + (*output_) << "\tSwingYAngle(radians)\tEnergy(kcal/mol)"; } - output_ << "\n"; -#ifdef IS_MPI } -#endif - } + + if (createRestFile_) (*output_) << "\n"; + if (createRestFile_) (*output_).flush(); + +#else + + std::string buffer; + + for(resti=restraints.begin(); resti != restraints.end(); ++resti){ + if ((*resti)->getPrintRestraint()) { + std::string myName = (*resti)->getRestraintName(); + int myType = (*resti)->getRestraintType(); + + buffer += (myName + ":"); + + if (myType & Restraint::rtDisplacement) + buffer += "\tPosition(angstroms)\tEnergy(kcal/mol)"; + + if (myType & Restraint::rtTwist) + buffer += "\tTwistAngle(radians)\tEnergy(kcal/mol)"; + + if (myType & Restraint::rtSwingX) + buffer += "\tSwingXAngle(radians)\tEnergy(kcal/mol)"; + + if (myType & Restraint::rtSwingY) + buffer += "\tSwingYAngle(radians)\tEnergy(kcal/mol)"; + + buffer += "\n"; + } + } + + const int masterNode = 0; + + if (worldRank == masterNode) { + if (createRestFile_) (*output_) << "#time\t"; + if (createRestFile_) (*output_) << buffer; + + int nProc; + MPI_Comm_size( MPI_COMM_WORLD, &nProc); + + for (int i = 1; i < nProc; ++i) { + + // receive the length of the string buffer that was + // prepared by processor i + + int recvLength; + MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, istatus); + char* recvBuffer = new char[recvLength]; + if (recvBuffer == NULL) { + } else { + MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, + istatus); + if (createRestFile_) (*output_) << recvBuffer; + delete [] recvBuffer; + } + } + if (createRestFile_) (*output_).flush(); + } else { + int sendBufferLength = buffer.size() + 1; + 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); + } + +#endif // is_mpi + + } - RestWriter::~RestWriter() { + void RestWriter::writeRest(std::vector > restInfo) { + #ifdef IS_MPI - if(worldRank == 0 ){ -#endif - output_.close(); -#ifdef IS_MPI - } + MPI_Status* istatus; #endif - } - - void RestWriter::writeRest(std::vector > restInfo){ - - output_ << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); - +#ifndef IS_MPI + if (createRestFile_) (*output_) << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); + // output some information about the molecules std::vector >::const_iterator i; std::map::const_iterator j; - for( i = restInfo.begin(); i != restInfo.end(); ++i){ - for(j = (*i).begin(); j != (*i).end(); ++j){ - output_ << "\t" << (j->second).first << "\t" << (j->second).second; + + if ( createRestFile_ ) { + + for( i = restInfo.begin(); i != restInfo.end(); ++i){ + for(j = (*i).begin(); j != (*i).end(); ++j){ + (*output_) << "\t" << (j->second).first << "\t" << (j->second).second; + } + (*output_) << std::endl; + } + (*output_).flush(); + } +#else + std::string buffer, first, second; + std::stringstream ss; + + std::vector >::const_iterator i; + std::map::const_iterator j; + + if ( createRestFile_ ) { + for( i = restInfo.begin(); i != restInfo.end(); ++i){ + + for(j = (*i).begin(); j != (*i).end(); ++j){ + ss.clear(); + ss << (j->second).first; + ss >> first; + ss.clear(); + ss << (j->second).second; + ss >> second; + buffer += ("\t" + first + "\t" + second); + } + buffer += "\n"; } - output_ << std::endl; } + + const int masterNode = 0; + + if (createRestFile_) { + if (worldRank == masterNode) { + + (*output_) << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); + (*output_) << buffer; + + int nProc; + MPI_Comm_size( MPI_COMM_WORLD, &nProc); + for (int i = 1; i < nProc; ++i) { + + // receive the length of the string buffer that was + // prepared by processor i + + int recvLength; + MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, istatus); + char* recvBuffer = new char[recvLength]; + if (recvBuffer == NULL) { + } else { + MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, + istatus); + if (createRestFile_) (*output_) << recvBuffer; + + delete [] recvBuffer; + } + } + (*output_).flush(); + } else { + int sendBufferLength = buffer.size() + 1; + 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); + } + } +#endif // is_mpi } -}// end oopse + RestWriter::~RestWriter() { + +#ifdef IS_MPI + + if (worldRank == 0) { +#endif // is_mpi + if (createRestFile_){ + writeClosing(*output_); + delete output_; + } +#ifdef IS_MPI + } +#endif // is_mpi + } + + void RestWriter::writeClosing(std::ostream& os) { + os.flush(); + } + +}// end namespace OpenMD +