--- trunk/src/io/RestWriter.cpp 2005/03/10 15:10:24 417 +++ branches/development/src/io/RestWriter.cpp 2012/09/13 14:10:11 1798 @@ -1,24 +1,15 @@ - /* - * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. +/* + * Copyright (c) 2009 The University of Notre Dame. All Rights Reserved. * * The University of Notre Dame grants you ("Licensee") a * non-exclusive, royalty free, license to use, modify and * 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,151 +28,249 @@ * 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, 24107 (2008). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ - -#include + + +#include +#include #include -#include -#include "primitives/Molecule.hpp" #include "io/RestWriter.hpp" #include "utils/simError.h" +#include "brains/SnapshotManager.hpp" +#ifdef IS_MPI +#include +#endif - -namespace oopse { - RestWriter::RestWriter(SimInfo* info) : - info_(info) { +namespace OpenMD { + RestWriter::RestWriter(SimInfo* info, const std::string& filename, + std::vector restraints ) : + info_(info){ + createRestFile_ = true; - //we use master - slave mode, only master node writes to disk - outName = info_->getRestFileName(); - } - - RestWriter::~RestWriter() {} - - void RestWriter::writeZangle(){ - const int BUFFERSIZE = 2000; - char tempBuffer[BUFFERSIZE]; - char writeLine[BUFFERSIZE]; - - std::ofstream finalOut; - - Molecule* mol; - StuntDouble* integrableObject; - SimInfo::MoleculeIterator mi; - Molecule::IntegrableObjectIterator ii; - #ifdef IS_MPI - if(worldRank == 0 ){ -#endif - finalOut.open( outName.c_str(), std::ios::out | std::ios::trunc ); - if( !finalOut ){ + if(worldRank == 0){ +#endif + + output_ = new std::ofstream(filename.c_str()); + + if(!output_){ sprintf( painCave.errMsg, - "Could not open \"%s\" for zAngle output.\n", - outName ); + "Could not open %s for restraint output.\n", + filename.c_str()); painCave.isFatal = 1; simError(); } + #ifdef IS_MPI } #endif // is_mpi + + +#ifdef IS_MPI + MPI::Status istatus; +#endif #ifndef IS_MPI - // first we do output for the single processor version - finalOut - << info_->getSnapshotManager()->getCurrentSnapshot()->getTime() - << " : omega values at this time\n"; - - for (mol = info_->beginMolecule(mi); mol != NULL; - mol = info_->nextMolecule(mi)) { - - for (integrableObject = mol->beginIntegrableObject(ii); - integrableObject != NULL; - integrableObject = mol->nextIntegrableObject(ii)) { + + (*output_) << "#time\t"; + + std::vector::const_iterator resti; + + for(resti=restraints.begin(); resti != restraints.end(); ++resti){ + if ((*resti)->getPrintRestraint()) { - sprintf( tempBuffer, - "%14.10lf\n", - integrableObject->getZangle()); - strcpy( writeLine, tempBuffer ); + std::string myName = (*resti)->getRestraintName(); + int myType = (*resti)->getRestraintType(); - finalOut << writeLine; + (*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)"; } - } + + (*output_) << "\n"; + (*output_).flush(); #else - int nproc; - MPI_Comm_size(MPI_COMM_WORLD, &nproc); - const int masterNode = 0; - int myNode = worldRank; - std::vector tmpNIntObjects(nproc, 0); - std::vector nIntObjectsInProc(nproc, 0); - tmpNIntObjects[myNode] = info_->getNGlobalIntegrableObjects(); - //do MPI_ALLREDUCE to exchange the total number of atoms, rigidbodies and cutoff groups - MPI_Allreduce(&tmpNIntObjects[0], &nIntObjectsInProc[0], nproc, MPI_INT, - MPI_SUM, MPI_COMM_WORLD); + std::string buffer; + + 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(); + + 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"; + } + } - MPI_Status ierr; - int intObIndex; - double zAngle; + const int masterNode = 0; - if (masterNode == 0) { - std::map zAngData; - for(int i = 0 ; i < nproc; ++i) { - if (i == masterNode) { - for (mol = info_->beginMolecule(mi); mol != NULL; - mol = info_->nextMolecule(mi)) { - - for (integrableObject = mol->beginIntegrableObject(ii); - integrableObject != NULL; - integrableObject = mol->nextIntegrableObject(ii)) { - - intObIndex = integrableObject->getGlobalIndex() ; - zAngle = integrableObject->getZangle(); - zAngData.insert(pair(intObIndex, zAngle)); - } - } - + if (worldRank == masterNode) { + (*output_) << "#time\t"; + (*output_) << buffer; + + int nProc = MPI::COMM_WORLD.Get_size(); + for (int i = 1; i < nProc; ++i) { + + // receive the length of the string buffer that was + // prepared by processor i + + int recvLength; + MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, 0, istatus); + char* recvBuffer = new char[recvLength]; + if (recvBuffer == NULL) { } else { - for(int k = 0; k < nIntObjectsInProc[i]; ++k) { - MPI_Recv(&intObIndex, 1, MPI_INT, i, 0, MPI_COMM_WORLD,&ierr); - MPI_Recv(&zAngle, 1, MPI_DOUBLE, i, 0, MPI_COMM_WORLD,&ierr); - zAngData.insert(pair(intObIndex, zAngle)); - } + MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i, 0, + istatus); + (*output_) << recvBuffer; + delete [] recvBuffer; } - + } + (*output_).flush(); + } else { + int sendBufferLength = buffer.size() + 1; + MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0); + MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength, MPI::CHAR, + masterNode, 0); + } + +#endif // is_mpi + + } + + void RestWriter::writeRest(std::vector > restInfo) { + +#ifdef IS_MPI + MPI::Status istatus; +#endif + +#ifndef IS_MPI + (*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; } - - finalOut - << info_->getSnapshotManager()->getCurrentSnapshot()->getTime() - << " : omega values at this time\n"; - - std::map::iterator l; - for (l = zAngData.begin(); l != zAngData.end(); ++l) { - finalOut << l->second << "\n"; + (*output_) << std::endl; + } + (*output_).flush(); +#else + std::string buffer, first, second; + std::stringstream ss; + + 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){ + ss.clear(); + ss << (j->second).first; + ss >> first; + ss.clear(); + ss << (j->second).second; + ss >> second; + buffer += ("\t" + first + "\t" + second); } + buffer += "\n"; + } + + const int masterNode = 0; + + if (worldRank == masterNode) { + (*output_) << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); + (*output_) << buffer; - } else { - - for (mol = info_->beginMolecule(mi); mol != NULL; - mol = info_->nextMolecule(mi)) { + int nProc; + nProc = MPI::COMM_WORLD.Get_size(); + for (int i = 1; i < nProc; ++i) { - for (integrableObject = mol->beginIntegrableObject(ii); - integrableObject != NULL; - integrableObject = mol->nextIntegrableObject(ii)) { - intObIndex = integrableObject->getGlobalIndex(); - zAngle = integrableObject->getZangle(); - MPI_Send(&intObIndex, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD); - MPI_Send(&zAngle, 1, MPI_DOUBLE, masterNode, 0, MPI_COMM_WORLD); + // receive the length of the string buffer that was + // prepared by processor i + + int recvLength; + MPI::COMM_WORLD.Recv(&recvLength, 1, MPI::INT, i, 0, istatus); + char* recvBuffer = new char[recvLength]; + if (recvBuffer == NULL) { + } else { + MPI::COMM_WORLD.Recv(recvBuffer, recvLength, MPI::CHAR, i, 0, + istatus); + (*output_) << recvBuffer; + + delete [] recvBuffer; } - } + } + (*output_).flush(); + } else { + int sendBufferLength = buffer.size() + 1; + MPI::COMM_WORLD.Send(&sendBufferLength, 1, MPI::INT, masterNode, 0); + MPI::COMM_WORLD.Send((void *)buffer.c_str(), sendBufferLength, + MPI::CHAR, masterNode, 0); } -#endif +#endif // is_mpi + } + + + RestWriter::~RestWriter() { #ifdef IS_MPI - finalOut.close(); -#endif + 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 +