| 1 | /* | 
| 2 | * Copyright (c) 2009 The University of Notre Dame. All Rights Reserved. | 
| 3 | * | 
| 4 | * The University of Notre Dame grants you ("Licensee") a | 
| 5 | * non-exclusive, royalty free, license to use, modify and | 
| 6 | * redistribute this software in source and binary code form, provided | 
| 7 | * that the following conditions are met: | 
| 8 | * | 
| 9 | * 1. Redistributions of source code must retain the above copyright | 
| 10 | *    notice, this list of conditions and the following disclaimer. | 
| 11 | * | 
| 12 | * 2. Redistributions in binary form must reproduce the above copyright | 
| 13 | *    notice, this list of conditions and the following disclaimer in the | 
| 14 | *    documentation and/or other materials provided with the | 
| 15 | *    distribution. | 
| 16 | * | 
| 17 | * This software is provided "AS IS," without a warranty of any | 
| 18 | * kind. All express or implied conditions, representations and | 
| 19 | * warranties, including any implied warranty of merchantability, | 
| 20 | * fitness for a particular purpose or non-infringement, are hereby | 
| 21 | * excluded.  The University of Notre Dame and its licensors shall not | 
| 22 | * be liable for any damages suffered by licensee as a result of | 
| 23 | * using, modifying or distributing the software or its | 
| 24 | * derivatives. In no event will the University of Notre Dame or its | 
| 25 | * licensors be liable for any lost revenue, profit or data, or for | 
| 26 | * direct, indirect, special, consequential, incidental or punitive | 
| 27 | * damages, however caused and regardless of the theory of liability, | 
| 28 | * arising out of the use of or inability to use software, even if the | 
| 29 | * University of Notre Dame has been advised of the possibility of | 
| 30 | * such damages. | 
| 31 | * | 
| 32 | * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your | 
| 33 | * research, please cite the appropriate papers when you publish your | 
| 34 | * work.  Good starting points are: | 
| 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, 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 | #ifdef IS_MPI | 
| 44 | #include <mpi.h> | 
| 45 | #endif | 
| 46 |  | 
| 47 | #include <string> | 
| 48 | #include <sstream> | 
| 49 | #include <iostream> | 
| 50 |  | 
| 51 | #include "io/RestWriter.hpp" | 
| 52 | #include "utils/simError.h" | 
| 53 | #include "brains/SnapshotManager.hpp" | 
| 54 |  | 
| 55 | namespace OpenMD { | 
| 56 | RestWriter::RestWriter(SimInfo* info, const std::string& filename, | 
| 57 | std::vector<Restraint*> restraints ) : | 
| 58 | info_(info){ | 
| 59 | createRestFile_ = true; | 
| 60 |  | 
| 61 | #ifdef IS_MPI | 
| 62 | if(worldRank == 0){ | 
| 63 | #endif | 
| 64 |  | 
| 65 | output_ = new std::ofstream(filename.c_str()); | 
| 66 |  | 
| 67 | if(!output_){ | 
| 68 | sprintf( painCave.errMsg, | 
| 69 | "Could not open %s for restraint output.\n", | 
| 70 | filename.c_str()); | 
| 71 | painCave.isFatal = 1; | 
| 72 | simError(); | 
| 73 | } | 
| 74 |  | 
| 75 | #ifdef IS_MPI | 
| 76 | } | 
| 77 | #endif // is_mpi | 
| 78 |  | 
| 79 |  | 
| 80 | #ifdef IS_MPI | 
| 81 | MPI_Status* istatus; | 
| 82 | #endif | 
| 83 |  | 
| 84 | #ifndef IS_MPI | 
| 85 |  | 
| 86 | (*output_) << "#time\t"; | 
| 87 |  | 
| 88 | std::vector<Restraint*>::const_iterator resti; | 
| 89 |  | 
| 90 | for(resti=restraints.begin(); resti != restraints.end(); ++resti){ | 
| 91 | if ((*resti)->getPrintRestraint()) { | 
| 92 |  | 
| 93 | std::string myName = (*resti)->getRestraintName(); | 
| 94 | int myType = (*resti)->getRestraintType(); | 
| 95 |  | 
| 96 | (*output_) << myName << ":"; | 
| 97 |  | 
| 98 | if (myType & Restraint::rtDisplacement) | 
| 99 | (*output_) << "\tPosition(angstroms)\tEnergy(kcal/mol)"; | 
| 100 |  | 
| 101 | if (myType & Restraint::rtTwist) | 
| 102 | (*output_) << "\tTwistAngle(radians)\tEnergy(kcal/mol)"; | 
| 103 |  | 
| 104 | if (myType & Restraint::rtSwingX) | 
| 105 | (*output_) << "\tSwingXAngle(radians)\tEnergy(kcal/mol)"; | 
| 106 |  | 
| 107 | if (myType & Restraint::rtSwingY) | 
| 108 | (*output_) << "\tSwingYAngle(radians)\tEnergy(kcal/mol)"; | 
| 109 | } | 
| 110 | } | 
| 111 |  | 
| 112 | (*output_) << "\n"; | 
| 113 | (*output_).flush(); | 
| 114 |  | 
| 115 | #else | 
| 116 |  | 
| 117 | std::string buffer; | 
| 118 |  | 
| 119 | std::vector<Restraint*>::const_iterator resti; | 
| 120 |  | 
| 121 | for(resti=restraints.begin(); resti != restraints.end(); ++resti){ | 
| 122 | if ((*resti)->getPrintRestraint()) { | 
| 123 |  | 
| 124 | std::string myName = (*resti)->getRestraintName(); | 
| 125 | int myType = (*resti)->getRestraintType(); | 
| 126 |  | 
| 127 | buffer += (myName + ":"); | 
| 128 |  | 
| 129 | if (myType & Restraint::rtDisplacement) | 
| 130 | buffer += "\tPosition(angstroms)\tEnergy(kcal/mol)"; | 
| 131 |  | 
| 132 | if (myType & Restraint::rtTwist) | 
| 133 | buffer += "\tTwistAngle(radians)\tEnergy(kcal/mol)"; | 
| 134 |  | 
| 135 | if (myType & Restraint::rtSwingX) | 
| 136 | buffer += "\tSwingXAngle(radians)\tEnergy(kcal/mol)"; | 
| 137 |  | 
| 138 | if (myType & Restraint::rtSwingY) | 
| 139 | buffer += "\tSwingYAngle(radians)\tEnergy(kcal/mol)"; | 
| 140 |  | 
| 141 | buffer += "\n"; | 
| 142 | } | 
| 143 | } | 
| 144 |  | 
| 145 | const int masterNode = 0; | 
| 146 |  | 
| 147 | if (worldRank == masterNode) { | 
| 148 | (*output_) << "#time\t"; | 
| 149 | (*output_) << buffer; | 
| 150 |  | 
| 151 | int nProc; | 
| 152 | MPI_Comm_size( MPI_COMM_WORLD, &nProc); | 
| 153 |  | 
| 154 | for (int i = 1; i < nProc; ++i) { | 
| 155 |  | 
| 156 | // receive the length of the string buffer that was | 
| 157 | // prepared by processor i | 
| 158 |  | 
| 159 | int recvLength; | 
| 160 | MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, istatus); | 
| 161 | char* recvBuffer = new char[recvLength]; | 
| 162 | if (recvBuffer == NULL) { | 
| 163 | } else { | 
| 164 | MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, | 
| 165 | istatus); | 
| 166 | (*output_) << recvBuffer; | 
| 167 | delete [] recvBuffer; | 
| 168 | } | 
| 169 | } | 
| 170 | (*output_).flush(); | 
| 171 | } else { | 
| 172 | int sendBufferLength = buffer.size() + 1; | 
| 173 | MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD); | 
| 174 | MPI_Send((void *)buffer.c_str(), sendBufferLength, MPI_CHAR, | 
| 175 | masterNode, 0, MPI_COMM_WORLD); | 
| 176 | } | 
| 177 |  | 
| 178 | #endif // is_mpi | 
| 179 |  | 
| 180 | } | 
| 181 |  | 
| 182 | void RestWriter::writeRest(std::vector<std::map<int, Restraint::RealPair> > restInfo) { | 
| 183 |  | 
| 184 | #ifdef IS_MPI | 
| 185 | MPI_Status* istatus; | 
| 186 | #endif | 
| 187 |  | 
| 188 | #ifndef IS_MPI | 
| 189 | (*output_) << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); | 
| 190 |  | 
| 191 | // output some information about the molecules | 
| 192 | std::vector<std::map<int, Restraint::RealPair> >::const_iterator i; | 
| 193 | std::map<int, Restraint::RealPair>::const_iterator j; | 
| 194 |  | 
| 195 | for( i = restInfo.begin(); i != restInfo.end(); ++i){ | 
| 196 | for(j = (*i).begin(); j != (*i).end(); ++j){ | 
| 197 | (*output_) << "\t" << (j->second).first << "\t" << (j->second).second; | 
| 198 | } | 
| 199 | (*output_) << std::endl; | 
| 200 | } | 
| 201 | (*output_).flush(); | 
| 202 | #else | 
| 203 | std::string buffer, first, second; | 
| 204 | std::stringstream ss; | 
| 205 |  | 
| 206 | std::vector<std::map<int, Restraint::RealPair> >::const_iterator i; | 
| 207 | std::map<int, Restraint::RealPair>::const_iterator j; | 
| 208 |  | 
| 209 | for( i = restInfo.begin(); i != restInfo.end(); ++i){ | 
| 210 | for(j = (*i).begin(); j != (*i).end(); ++j){ | 
| 211 | ss.clear(); | 
| 212 | ss << (j->second).first; | 
| 213 | ss >> first; | 
| 214 | ss.clear(); | 
| 215 | ss << (j->second).second; | 
| 216 | ss >> second; | 
| 217 | buffer += ("\t" + first + "\t" + second); | 
| 218 | } | 
| 219 | buffer += "\n"; | 
| 220 | } | 
| 221 |  | 
| 222 | const int masterNode = 0; | 
| 223 |  | 
| 224 | if (worldRank == masterNode) { | 
| 225 | (*output_) << info_->getSnapshotManager()->getCurrentSnapshot()->getTime(); | 
| 226 | (*output_) << buffer; | 
| 227 |  | 
| 228 | int nProc; | 
| 229 | MPI_Comm_size( MPI_COMM_WORLD, &nProc); | 
| 230 | for (int i = 1; i < nProc; ++i) { | 
| 231 |  | 
| 232 | // receive the length of the string buffer that was | 
| 233 | // prepared by processor i | 
| 234 |  | 
| 235 | int recvLength; | 
| 236 | MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, istatus); | 
| 237 | char* recvBuffer = new char[recvLength]; | 
| 238 | if (recvBuffer == NULL) { | 
| 239 | } else { | 
| 240 | MPI_Recv(recvBuffer, recvLength, MPI_CHAR, i, 0, MPI_COMM_WORLD, | 
| 241 | istatus); | 
| 242 | (*output_) << recvBuffer; | 
| 243 |  | 
| 244 | delete [] recvBuffer; | 
| 245 | } | 
| 246 | } | 
| 247 | (*output_).flush(); | 
| 248 | } else { | 
| 249 | int sendBufferLength = buffer.size() + 1; | 
| 250 | MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD); | 
| 251 | MPI_Send((void *)buffer.c_str(), sendBufferLength, | 
| 252 | MPI_CHAR, masterNode, 0, MPI_COMM_WORLD); | 
| 253 | } | 
| 254 | #endif // is_mpi | 
| 255 | } | 
| 256 |  | 
| 257 |  | 
| 258 | RestWriter::~RestWriter() { | 
| 259 |  | 
| 260 | #ifdef IS_MPI | 
| 261 |  | 
| 262 | if (worldRank == 0) { | 
| 263 | #endif // is_mpi | 
| 264 | if (createRestFile_){ | 
| 265 | writeClosing(*output_); | 
| 266 | delete output_; | 
| 267 | } | 
| 268 | #ifdef IS_MPI | 
| 269 | } | 
| 270 | #endif // is_mpi | 
| 271 | } | 
| 272 |  | 
| 273 | void RestWriter::writeClosing(std::ostream& os) { | 
| 274 | os.flush(); | 
| 275 | } | 
| 276 |  | 
| 277 | }// end namespace OpenMD | 
| 278 |  |