| 1 | /* | 
| 2 | * Copyright (c) 2005 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. Acknowledgement of the program authors must be made in any | 
| 10 | *    publication of scientific results based in part on use of the | 
| 11 | *    program.  An acceptable form of acknowledgement is citation of | 
| 12 | *    the article in which the program was described (Matthew | 
| 13 | *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher | 
| 14 | *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented | 
| 15 | *    Parallel Simulation Engine for Molecular Dynamics," | 
| 16 | *    J. Comput. Chem. 26, pp. 252-271 (2005)) | 
| 17 | * | 
| 18 | * 2. Redistributions of source code must retain the above copyright | 
| 19 | *    notice, this list of conditions and the following disclaimer. | 
| 20 | * | 
| 21 | * 3. Redistributions in binary form must reproduce the above copyright | 
| 22 | *    notice, this list of conditions and the following disclaimer in the | 
| 23 | *    documentation and/or other materials provided with the | 
| 24 | *    distribution. | 
| 25 | * | 
| 26 | * This software is provided "AS IS," without a warranty of any | 
| 27 | * kind. All express or implied conditions, representations and | 
| 28 | * warranties, including any implied warranty of merchantability, | 
| 29 | * fitness for a particular purpose or non-infringement, are hereby | 
| 30 | * excluded.  The University of Notre Dame and its licensors shall not | 
| 31 | * be liable for any damages suffered by licensee as a result of | 
| 32 | * using, modifying or distributing the software or its | 
| 33 | * derivatives. In no event will the University of Notre Dame or its | 
| 34 | * licensors be liable for any lost revenue, profit or data, or for | 
| 35 | * direct, indirect, special, consequential, incidental or punitive | 
| 36 | * damages, however caused and regardless of the theory of liability, | 
| 37 | * arising out of the use of or inability to use software, even if the | 
| 38 | * University of Notre Dame has been advised of the possibility of | 
| 39 | * such damages. | 
| 40 | */ | 
| 41 |  | 
| 42 | #define _LARGEFILE_SOURCE64 | 
| 43 | #define _FILE_OFFSET_BITS 64 | 
| 44 |  | 
| 45 | #include <sys/types.h> | 
| 46 | #include <sys/stat.h> | 
| 47 |  | 
| 48 | #include <iostream> | 
| 49 | #include <math.h> | 
| 50 |  | 
| 51 | #include <stdio.h> | 
| 52 | #include <stdlib.h> | 
| 53 | #include <string.h> | 
| 54 |  | 
| 55 | #include "primitives/Molecule.hpp" | 
| 56 | #include "utils/MemoryUtils.hpp" | 
| 57 | #include "utils/StringTokenizer.hpp" | 
| 58 | #include "io/RestReader.hpp" | 
| 59 | #include "utils/simError.h" | 
| 60 |  | 
| 61 | #ifdef IS_MPI | 
| 62 | #include <mpi.h> | 
| 63 | #define TAKE_THIS_TAG_CHAR 0 | 
| 64 | #define TAKE_THIS_TAG_INT 1 | 
| 65 | #define TAKE_THIS_TAG_DOUBLE 2 | 
| 66 | #endif // is_mpi | 
| 67 |  | 
| 68 | namespace oopse { | 
| 69 |  | 
| 70 | RestReader::RestReader( SimInfo* info ) : info_(info){ | 
| 71 |  | 
| 72 | idealName = "idealCrystal.in"; | 
| 73 |  | 
| 74 | #ifdef IS_MPI | 
| 75 | if (worldRank == 0) { | 
| 76 | #endif | 
| 77 |  | 
| 78 | inIdealFile = new std::ifstream(idealName.c_str()); | 
| 79 |  | 
| 80 | if(inIdealFile->fail()){ | 
| 81 | sprintf(painCave.errMsg, | 
| 82 | "RestReader: Cannot open file: %s\n", | 
| 83 | idealName.c_str()); | 
| 84 | painCave.isFatal = 1; | 
| 85 | simError(); | 
| 86 | } | 
| 87 |  | 
| 88 | #ifdef IS_MPI | 
| 89 | } | 
| 90 | strcpy( checkPointMsg, | 
| 91 | "File \"idealCrystal.in\" opened successfully for reading." ); | 
| 92 | MPIcheckPoint(); | 
| 93 | #endif | 
| 94 |  | 
| 95 | return; | 
| 96 | } | 
| 97 |  | 
| 98 | RestReader :: ~RestReader( ){ | 
| 99 | #ifdef IS_MPI | 
| 100 | if (worldRank == 0) { | 
| 101 | #endif | 
| 102 |  | 
| 103 | delete inIdealFile; | 
| 104 | delete inAngFile; | 
| 105 |  | 
| 106 | #ifdef IS_MPI | 
| 107 | } | 
| 108 | strcpy( checkPointMsg, | 
| 109 | "File idealCrystal.in (and .zang0 if present) closed successfully." ); | 
| 110 | MPIcheckPoint(); | 
| 111 | #endif | 
| 112 |  | 
| 113 | return; | 
| 114 | } | 
| 115 |  | 
| 116 |  | 
| 117 | void RestReader :: readIdealCrystal(){ | 
| 118 |  | 
| 119 | #ifdef IS_MPI | 
| 120 | int which_node; | 
| 121 | int i, j; | 
| 122 | #endif //is_mpi | 
| 123 |  | 
| 124 | const int BUFFERSIZE = 2000; // size of the read buffer | 
| 125 | int nTotObjs; // the number of atoms | 
| 126 | char read_buffer[BUFFERSIZE]; //the line buffer for reading | 
| 127 |  | 
| 128 | char *parseErr; | 
| 129 |  | 
| 130 | std::vector<StuntDouble*> integrableObjects; | 
| 131 |  | 
| 132 | Molecule* mol; | 
| 133 | StuntDouble* integrableObject; | 
| 134 | SimInfo::MoleculeIterator mi; | 
| 135 | Molecule::IntegrableObjectIterator ii; | 
| 136 |  | 
| 137 | #ifndef IS_MPI | 
| 138 |  | 
| 139 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 140 |  | 
| 141 | if( inIdealFile->eof() ){ | 
| 142 | sprintf( painCave.errMsg, | 
| 143 | "RestraintReader error: error reading 1st line of \"%s\"\n", | 
| 144 | idealName.c_str() ); | 
| 145 | painCave.isFatal = 1; | 
| 146 | simError(); | 
| 147 | } | 
| 148 |  | 
| 149 | nTotObjs = atoi( read_buffer ); | 
| 150 |  | 
| 151 | if( nTotObjs != info_->getNGlobalIntegrableObjects() ){ | 
| 152 | sprintf( painCave.errMsg, | 
| 153 | "RestraintReader error. %s nIntegrable, %d, " | 
| 154 | "does not match the meta-data file's nIntegrable, %d.\n", | 
| 155 | idealName.c_str(), | 
| 156 | nTotObjs, | 
| 157 | info_->getNGlobalIntegrableObjects()); | 
| 158 | painCave.isFatal = 1; | 
| 159 | simError(); | 
| 160 | } | 
| 161 |  | 
| 162 | // skip over the comment line | 
| 163 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 164 |  | 
| 165 | if( inIdealFile->eof() ){ | 
| 166 | sprintf( painCave.errMsg, | 
| 167 | "error in reading commment in %s\n", | 
| 168 | idealName.c_str()); | 
| 169 | painCave.isFatal = 1; | 
| 170 | simError(); | 
| 171 | } | 
| 172 |  | 
| 173 | // parse the ideal crystal lines | 
| 174 | /* | 
| 175 | * Note: we assume that there is a one-to-one correspondence between | 
| 176 | * integrable objects and lines in the idealCrystal.in file.  Thermodynamic | 
| 177 | * integration is only supported for simple rigid bodies. | 
| 178 | */ | 
| 179 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 180 | mol = info_->nextMolecule(mi)) { | 
| 181 |  | 
| 182 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 183 | integrableObject != NULL; | 
| 184 | integrableObject = mol->nextIntegrableObject(ii)) { | 
| 185 |  | 
| 186 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 187 |  | 
| 188 | if( inIdealFile->eof() ){ | 
| 189 | sprintf(painCave.errMsg, | 
| 190 | "RestReader Error: error in reading file %s\n" | 
| 191 | "natoms  = %d; index = %d\n" | 
| 192 | "error reading the line from the file.\n", | 
| 193 | idealName.c_str(), nTotObjs, | 
| 194 | integrableObject->getGlobalIndex() ); | 
| 195 | painCave.isFatal = 1; | 
| 196 | simError(); | 
| 197 | } | 
| 198 |  | 
| 199 | parseErr = parseIdealLine( read_buffer, integrableObject); | 
| 200 | if( parseErr != NULL ){ | 
| 201 | strcpy( painCave.errMsg, parseErr ); | 
| 202 | painCave.isFatal = 1; | 
| 203 | simError(); | 
| 204 | } | 
| 205 | } | 
| 206 | } | 
| 207 |  | 
| 208 | // MPI Section of code.......... | 
| 209 | #else //IS_MPI | 
| 210 |  | 
| 211 | // first thing first, suspend fatalities. | 
| 212 | painCave.isEventLoop = 1; | 
| 213 |  | 
| 214 | int masterNode = 0; | 
| 215 |  | 
| 216 | MPI_Status istatus; | 
| 217 | int nCurObj; | 
| 218 | int nitems; | 
| 219 | int haveError; | 
| 220 |  | 
| 221 | nTotObjs = info_->getNGlobalIntegrableObjects(); | 
| 222 | haveError = 0; | 
| 223 |  | 
| 224 | if (worldRank == masterNode) { | 
| 225 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 226 |  | 
| 227 | if( inIdealFile->eof() ){ | 
| 228 | sprintf( painCave.errMsg, | 
| 229 | "Error reading 1st line of %s \n ",idealName.c_str()); | 
| 230 | painCave.isFatal = 1; | 
| 231 | simError(); | 
| 232 | } | 
| 233 |  | 
| 234 | nitems = atoi( read_buffer ); | 
| 235 |  | 
| 236 | // Check to see that the number of integrable objects in the | 
| 237 | // intial configuration file is the same as derived from the | 
| 238 | // meta-data file. | 
| 239 | if( nTotObjs != nitems){ | 
| 240 | sprintf( painCave.errMsg, | 
| 241 | "RestraintReader Error. %s nIntegrable, %d, " | 
| 242 | "does not match the meta-data file's nIntegrable, %d.\n", | 
| 243 | idealName.c_str(), nTotObjs, | 
| 244 | info_->getNGlobalIntegrableObjects()); | 
| 245 | painCave.isFatal = 1; | 
| 246 | simError(); | 
| 247 | } | 
| 248 |  | 
| 249 | // skip over the comment line | 
| 250 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 251 |  | 
| 252 | if( inIdealFile->eof() ){ | 
| 253 | sprintf( painCave.errMsg, | 
| 254 | "error in reading commment in %s\n", idealName.c_str()); | 
| 255 | painCave.isFatal = 1; | 
| 256 | simError(); | 
| 257 | } | 
| 258 |  | 
| 259 | MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, masterNode, MPI_COMM_WORLD); | 
| 260 |  | 
| 261 | for (i=0 ; i < info_->getNGlobalMolecules(); i++) { | 
| 262 | which_node = info_->getMolToProc(i); | 
| 263 |  | 
| 264 | if(which_node == masterNode){ | 
| 265 | //molecules belong to master node | 
| 266 |  | 
| 267 | mol = info_->getMoleculeByGlobalIndex(i); | 
| 268 |  | 
| 269 | if(mol == NULL) { | 
| 270 | sprintf(painCave.errMsg, | 
| 271 | "RestReader Error: Molecule not found on node %d!\n", | 
| 272 | worldRank); | 
| 273 | painCave.isFatal = 1; | 
| 274 | simError(); | 
| 275 | } | 
| 276 |  | 
| 277 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 278 | integrableObject != NULL; | 
| 279 | integrableObject = mol->nextIntegrableObject(ii)){ | 
| 280 |  | 
| 281 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 282 |  | 
| 283 | if( inIdealFile->eof() ){ | 
| 284 | sprintf(painCave.errMsg, | 
| 285 | "RestReader Error: error in reading file %s\n" | 
| 286 | "natoms  = %d; index = %d\n" | 
| 287 | "error reading the line from the file.\n", | 
| 288 | idealName.c_str(), nTotObjs, i ); | 
| 289 | painCave.isFatal = 1; | 
| 290 | simError(); | 
| 291 | } | 
| 292 |  | 
| 293 | parseIdealLine(read_buffer, integrableObject); | 
| 294 |  | 
| 295 | } | 
| 296 |  | 
| 297 | } else { | 
| 298 | //molecule belongs to slave nodes | 
| 299 |  | 
| 300 | MPI_Recv(&nCurObj, 1, MPI_INT, which_node, | 
| 301 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus); | 
| 302 |  | 
| 303 | for(j = 0; j < nCurObj; j++){ | 
| 304 | inIdealFile->getline(read_buffer, sizeof(read_buffer)); | 
| 305 |  | 
| 306 | if( inIdealFile->eof() ){ | 
| 307 | sprintf(painCave.errMsg, | 
| 308 | "RestReader Error: error in reading file %s\n" | 
| 309 | "natoms  = %d; index = %d\n" | 
| 310 | "error reading the line from the file.\n", | 
| 311 | idealName.c_str(), nTotObjs, i ); | 
| 312 | painCave.isFatal = 1; | 
| 313 | simError(); | 
| 314 | } | 
| 315 |  | 
| 316 | MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node, | 
| 317 | TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD); | 
| 318 | } | 
| 319 | } | 
| 320 | } | 
| 321 | } else { | 
| 322 | //actions taken at slave nodes | 
| 323 | MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, masterNode, MPI_COMM_WORLD); | 
| 324 |  | 
| 325 | for (i=0 ; i < info_->getNGlobalMolecules(); i++) { | 
| 326 | int which_node = info_->getMolToProc(i); | 
| 327 |  | 
| 328 | if(which_node == worldRank){ | 
| 329 | //molecule with global index i belongs to this processor | 
| 330 |  | 
| 331 | mol = info_->getMoleculeByGlobalIndex(i); | 
| 332 |  | 
| 333 | if(mol == NULL) { | 
| 334 | sprintf(painCave.errMsg, | 
| 335 | "RestReader Error: molecule not found on node %d\n", | 
| 336 | worldRank); | 
| 337 | painCave.isFatal = 1; | 
| 338 | simError(); | 
| 339 | } | 
| 340 |  | 
| 341 | nCurObj = mol->getNIntegrableObjects(); | 
| 342 |  | 
| 343 | MPI_Send(&nCurObj, 1, MPI_INT, masterNode, | 
| 344 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD); | 
| 345 |  | 
| 346 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 347 | integrableObject != NULL; | 
| 348 | integrableObject = mol->nextIntegrableObject(ii)){ | 
| 349 |  | 
| 350 | MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, masterNode, | 
| 351 | TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus); | 
| 352 |  | 
| 353 | parseErr = parseIdealLine(read_buffer, integrableObject); | 
| 354 |  | 
| 355 | if( parseErr != NULL ){ | 
| 356 | strcpy( painCave.errMsg, parseErr ); | 
| 357 | simError(); | 
| 358 | } | 
| 359 | } | 
| 360 | } | 
| 361 | } | 
| 362 | } | 
| 363 | #endif | 
| 364 | } | 
| 365 |  | 
| 366 | char* RestReader::parseIdealLine(char* readLine, StuntDouble* sd){ | 
| 367 |  | 
| 368 | RealType pos[3];        // position place holders | 
| 369 | RealType q[4];          // the quaternions | 
| 370 | RealType RfromQ[3][3];  // the rotation matrix | 
| 371 | RealType normalize;     // to normalize the reference unit vector | 
| 372 | RealType uX, uY, uZ;    // reference unit vector place holders | 
| 373 | RealType uselessToken; | 
| 374 | StringTokenizer tokenizer(readLine); | 
| 375 | int nTokens; | 
| 376 |  | 
| 377 | nTokens = tokenizer.countTokens(); | 
| 378 |  | 
| 379 | if (nTokens < 14) { | 
| 380 | sprintf(painCave.errMsg, | 
| 381 | "RestReader Error: Not enough Tokens.\n"); | 
| 382 | painCave.isFatal = 1; | 
| 383 | simError(); | 
| 384 | } | 
| 385 |  | 
| 386 | std::string name = tokenizer.nextToken(); | 
| 387 |  | 
| 388 | if (name != sd->getType()) { | 
| 389 |  | 
| 390 | sprintf(painCave.errMsg, | 
| 391 | "RestReader Error: Atom type [%s] in %s does not " | 
| 392 | "match Atom Type [%s] in .md file.\n", | 
| 393 | name.c_str(), idealName.c_str(), | 
| 394 | sd->getType().c_str()); | 
| 395 | painCave.isFatal = 1; | 
| 396 | simError(); | 
| 397 | } | 
| 398 |  | 
| 399 | pos[0] = tokenizer.nextTokenAsDouble(); | 
| 400 | pos[1] = tokenizer.nextTokenAsDouble(); | 
| 401 | pos[2] = tokenizer.nextTokenAsDouble(); | 
| 402 |  | 
| 403 | // store the positions in the stuntdouble as generic data doubles | 
| 404 | DoubleGenericData* refPosX = new DoubleGenericData(); | 
| 405 | refPosX->setID("refPosX"); | 
| 406 | refPosX->setData(pos[0]); | 
| 407 | sd->addProperty(refPosX); | 
| 408 |  | 
| 409 | DoubleGenericData* refPosY = new DoubleGenericData(); | 
| 410 | refPosY->setID("refPosY"); | 
| 411 | refPosY->setData(pos[1]); | 
| 412 | sd->addProperty(refPosY); | 
| 413 |  | 
| 414 | DoubleGenericData* refPosZ = new DoubleGenericData(); | 
| 415 | refPosZ->setID("refPosZ"); | 
| 416 | refPosZ->setData(pos[2]); | 
| 417 | sd->addProperty(refPosZ); | 
| 418 |  | 
| 419 | // we don't need the velocities | 
| 420 | uselessToken = tokenizer.nextTokenAsDouble(); | 
| 421 | uselessToken = tokenizer.nextTokenAsDouble(); | 
| 422 | uselessToken = tokenizer.nextTokenAsDouble(); | 
| 423 |  | 
| 424 | if (sd->isDirectional()) { | 
| 425 |  | 
| 426 | q[0] = tokenizer.nextTokenAsDouble(); | 
| 427 | q[1] = tokenizer.nextTokenAsDouble(); | 
| 428 | q[2] = tokenizer.nextTokenAsDouble(); | 
| 429 | q[3] = tokenizer.nextTokenAsDouble(); | 
| 430 |  | 
| 431 | // now build the rotation matrix and find the unit vectors | 
| 432 | RfromQ[0][0] = q[0]*q[0] + q[1]*q[1] - q[2]*q[2] - q[3]*q[3]; | 
| 433 | RfromQ[0][1] = 2*(q[1]*q[2] + q[0]*q[3]); | 
| 434 | RfromQ[0][2] = 2*(q[1]*q[3] - q[0]*q[2]); | 
| 435 | RfromQ[1][0] = 2*(q[1]*q[2] - q[0]*q[3]); | 
| 436 | RfromQ[1][1] = q[0]*q[0] - q[1]*q[1] + q[2]*q[2] - q[3]*q[3]; | 
| 437 | RfromQ[1][2] = 2*(q[2]*q[3] + q[0]*q[1]); | 
| 438 | RfromQ[2][0] = 2*(q[1]*q[3] + q[0]*q[2]); | 
| 439 | RfromQ[2][1] = 2*(q[2]*q[3] - q[0]*q[1]); | 
| 440 | RfromQ[2][2] = q[0]*q[0] - q[1]*q[1] - q[2]*q[2] + q[3]*q[3]; | 
| 441 |  | 
| 442 | normalize = sqrt(RfromQ[2][0]*RfromQ[2][0] + RfromQ[2][1]*RfromQ[2][1] | 
| 443 | + RfromQ[2][2]*RfromQ[2][2]); | 
| 444 | uX = RfromQ[2][0]/normalize; | 
| 445 | uY = RfromQ[2][1]/normalize; | 
| 446 | uZ = RfromQ[2][2]/normalize; | 
| 447 |  | 
| 448 | // store reference unit vectors as generic data in the stuntdouble | 
| 449 | DoubleGenericData* refVectorX = new DoubleGenericData(); | 
| 450 | refVectorX->setID("refVectorX"); | 
| 451 | refVectorX->setData(uX); | 
| 452 | sd->addProperty(refVectorX); | 
| 453 |  | 
| 454 | DoubleGenericData* refVectorY = new DoubleGenericData(); | 
| 455 | refVectorY->setID("refVectorY"); | 
| 456 | refVectorY->setData(uY); | 
| 457 | sd->addProperty(refVectorY); | 
| 458 |  | 
| 459 | DoubleGenericData* refVectorZ = new DoubleGenericData(); | 
| 460 | refVectorZ->setID("refVectorZ"); | 
| 461 | refVectorZ->setData(uZ); | 
| 462 | sd->addProperty(refVectorZ); | 
| 463 | } | 
| 464 |  | 
| 465 | // we don't need the angular velocities, so let's exit the line | 
| 466 | return NULL; | 
| 467 | } | 
| 468 |  | 
| 469 | void RestReader::readZangle(){ | 
| 470 |  | 
| 471 | int i; | 
| 472 | int isPresent; | 
| 473 |  | 
| 474 | Molecule* mol; | 
| 475 | StuntDouble* integrableObject; | 
| 476 | SimInfo::MoleculeIterator mi; | 
| 477 | Molecule::IntegrableObjectIterator ii; | 
| 478 |  | 
| 479 | #ifdef IS_MPI | 
| 480 | int which_node; | 
| 481 | MPI_Status istatus; | 
| 482 | #endif //is_mpi | 
| 483 |  | 
| 484 | const int BUFFERSIZE = 2000; // size of the read buffer | 
| 485 | unsigned int nTotObjs; // the number of atoms | 
| 486 | char read_buffer[BUFFERSIZE]; //the line buffer for reading | 
| 487 |  | 
| 488 | std::vector<StuntDouble*> vecParticles; | 
| 489 | std::vector<RealType> tempZangs; | 
| 490 |  | 
| 491 | angFile = info_->getRestFileName(); | 
| 492 |  | 
| 493 | angFile += "0"; | 
| 494 |  | 
| 495 | // open the omega value file for reading | 
| 496 | #ifdef IS_MPI | 
| 497 | if (worldRank == 0) { | 
| 498 | #endif | 
| 499 | isPresent = 1; | 
| 500 |  | 
| 501 | inAngFile = new std::ifstream(angFile.c_str()); | 
| 502 |  | 
| 503 | if(inAngFile->fail()){ | 
| 504 | sprintf(painCave.errMsg, | 
| 505 | "Restraints Warning: %s file is not present\n" | 
| 506 | "\tAll omega values will be initialized to zero. If the\n" | 
| 507 | "\tsimulation is starting from the idealCrystal.in reference\n" | 
| 508 | "\tconfiguration, this is the desired action. If this is not\n" | 
| 509 | "\tthe case, the energy calculations will be incorrect.\n", | 
| 510 | angFile.c_str()); | 
| 511 | painCave.severity = OOPSE_WARNING; | 
| 512 | painCave.isFatal = 0; | 
| 513 | simError(); | 
| 514 | isPresent = 0; | 
| 515 | } | 
| 516 |  | 
| 517 | #ifdef IS_MPI | 
| 518 | // let the other nodes know the status of the file search | 
| 519 | MPI_Bcast(&isPresent, 1, MPI_INT, 0, MPI_COMM_WORLD); | 
| 520 | #endif // is_mpi | 
| 521 |  | 
| 522 | if (!isPresent) { | 
| 523 | zeroZangle(); | 
| 524 | return; | 
| 525 | } | 
| 526 |  | 
| 527 | #ifdef IS_MPI | 
| 528 | if (!isPresent) { | 
| 529 | // master node zeroes out its zAngles if .zang0 isn't present | 
| 530 | zeroZangle(); | 
| 531 | return; | 
| 532 | } | 
| 533 |  | 
| 534 | } | 
| 535 |  | 
| 536 | // listen to node 0 to see if we should exit this function | 
| 537 | if (worldRank != 0) { | 
| 538 | MPI_Bcast(&isPresent, 1, MPI_INT, 0, MPI_COMM_WORLD); | 
| 539 | if (!isPresent) { | 
| 540 | zeroZangle(); | 
| 541 | return; | 
| 542 | } | 
| 543 | } | 
| 544 |  | 
| 545 | strcpy( checkPointMsg, "zAngle file opened successfully for reading." ); | 
| 546 | MPIcheckPoint(); | 
| 547 | #endif | 
| 548 |  | 
| 549 | #ifndef IS_MPI | 
| 550 |  | 
| 551 | // read the first line and die if there is a failure | 
| 552 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 553 |  | 
| 554 | if( inAngFile->eof() ){ | 
| 555 | sprintf( painCave.errMsg, | 
| 556 | "RestraintReader error: error reading 1st line of \"%s\"\n", | 
| 557 | angFile.c_str() ); | 
| 558 | painCave.isFatal = 1; | 
| 559 | simError(); | 
| 560 | } | 
| 561 |  | 
| 562 | // read the file and load the values into a vector | 
| 563 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 564 |  | 
| 565 | while ( !inAngFile->eof() ) { | 
| 566 | // check for and ignore blank lines | 
| 567 | if ( read_buffer != NULL ) | 
| 568 | tempZangs.push_back( atof(read_buffer) ); | 
| 569 |  | 
| 570 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 571 | } | 
| 572 |  | 
| 573 | nTotObjs = info_->getNGlobalIntegrableObjects(); | 
| 574 |  | 
| 575 | if( nTotObjs != tempZangs.size() ){ | 
| 576 | sprintf( painCave.errMsg, | 
| 577 | "RestraintReader zAngle reading error. %s nIntegrable, %d, " | 
| 578 | "does not match the meta-data file's nIntegrable, %i.\n", | 
| 579 | angFile.c_str(), | 
| 580 | tempZangs.size(), | 
| 581 | nTotObjs ); | 
| 582 | painCave.isFatal = 1; | 
| 583 | simError(); | 
| 584 | } | 
| 585 |  | 
| 586 | // load the zAngles into the integrable objects | 
| 587 | i = 0; | 
| 588 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 589 | mol = info_->nextMolecule(mi)) { | 
| 590 |  | 
| 591 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 592 | integrableObject != NULL; | 
| 593 | integrableObject = mol->nextIntegrableObject(ii)) { | 
| 594 |  | 
| 595 | integrableObject->setZangle(tempZangs[i]); | 
| 596 | i++; | 
| 597 | } | 
| 598 | } | 
| 599 |  | 
| 600 | // MPI Section of code.......... | 
| 601 | #else //IS_MPI | 
| 602 |  | 
| 603 | // first thing first, suspend fatalities. | 
| 604 | painCave.isEventLoop = 1; | 
| 605 |  | 
| 606 | int masterNode = 0; | 
| 607 |  | 
| 608 | int haveError; | 
| 609 | int intObjIndex; | 
| 610 | int intObjIndexTransfer; | 
| 611 |  | 
| 612 | int j; | 
| 613 | int nCurObj; | 
| 614 | RealType angleTransfer; | 
| 615 |  | 
| 616 | nTotObjs = info_->getNGlobalIntegrableObjects(); | 
| 617 | haveError = 0; | 
| 618 |  | 
| 619 | if (worldRank == masterNode) { | 
| 620 |  | 
| 621 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 622 |  | 
| 623 | if( inAngFile->eof() ){ | 
| 624 | sprintf( painCave.errMsg, | 
| 625 | "Error reading 1st line of %s \n ",angFile.c_str()); | 
| 626 | haveError = 1; | 
| 627 | simError(); | 
| 628 | } | 
| 629 |  | 
| 630 | // let the master node read the file and load the temporary angle vector | 
| 631 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 632 |  | 
| 633 | while ( !inAngFile->eof() ) { | 
| 634 | // check for and ignore blank lines | 
| 635 | if ( read_buffer != NULL ) | 
| 636 | tempZangs.push_back( atof(read_buffer) ); | 
| 637 |  | 
| 638 | inAngFile->getline(read_buffer, sizeof(read_buffer)); | 
| 639 | } | 
| 640 |  | 
| 641 | // Check to see that the number of integrable objects in the | 
| 642 | // intial configuration file is the same as derived from the | 
| 643 | // meta-data file. | 
| 644 | if( nTotObjs != tempZangs.size() ){ | 
| 645 | sprintf( painCave.errMsg, | 
| 646 | "RestraintReader zAngle reading Error. %s nIntegrable, %d, " | 
| 647 | "does not match the meta-data file's nIntegrable, %d.\n", | 
| 648 | angFile.c_str(), | 
| 649 | tempZangs.size(), | 
| 650 | nTotObjs); | 
| 651 | haveError= 1; | 
| 652 | simError(); | 
| 653 | } | 
| 654 |  | 
| 655 | // At this point, node 0 has a tempZangs vector completed, and | 
| 656 | // everyone else has nada | 
| 657 |  | 
| 658 | for (i=0 ; i < info_->getNGlobalMolecules(); i++) { | 
| 659 | // Get the Node number which has this atom | 
| 660 | which_node = info_->getMolToProc(i); | 
| 661 |  | 
| 662 | if (which_node == masterNode) { | 
| 663 | mol = info_->getMoleculeByGlobalIndex(i); | 
| 664 |  | 
| 665 | if(mol == NULL) { | 
| 666 | strcpy(painCave.errMsg, "Molecule not found on node 0!"); | 
| 667 | haveError = 1; | 
| 668 | simError(); | 
| 669 | } | 
| 670 |  | 
| 671 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 672 | integrableObject != NULL; | 
| 673 | integrableObject = mol->nextIntegrableObject(ii)){ | 
| 674 | intObjIndex = integrableObject->getGlobalIntegrableObjectIndex(); | 
| 675 | integrableObject->setZangle(tempZangs[intObjIndex]); | 
| 676 | } | 
| 677 |  | 
| 678 | } else { | 
| 679 | // I am MASTER OF THE UNIVERSE, but I don't own this molecule | 
| 680 | // listen for the number of integrableObjects in the molecule | 
| 681 | MPI_Recv(&nCurObj, 1, MPI_INT, which_node, | 
| 682 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus); | 
| 683 |  | 
| 684 | for(j=0; j < nCurObj; j++){ | 
| 685 | // listen for which integrableObject we need to send the value for | 
| 686 | MPI_Recv(&intObjIndexTransfer, 1, MPI_INT, which_node, | 
| 687 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus); | 
| 688 | angleTransfer = tempZangs[intObjIndexTransfer]; | 
| 689 | // send the value to the node so it can initialize the object | 
| 690 | MPI_Send(&angleTransfer, 1, MPI_REALTYPE, which_node, | 
| 691 | TAKE_THIS_TAG_DOUBLE, MPI_COMM_WORLD); | 
| 692 | } | 
| 693 | } | 
| 694 | } | 
| 695 | } else { | 
| 696 | // I am SLAVE TO THE MASTER | 
| 697 | for (i=0 ; i < info_->getNGlobalMolecules(); i++) { | 
| 698 | which_node = info_->getMolToProc(i); | 
| 699 |  | 
| 700 | if (which_node == worldRank) { | 
| 701 |  | 
| 702 | // BUT I OWN THIS MOLECULE!!! | 
| 703 |  | 
| 704 | mol = info_->getMoleculeByGlobalIndex(i); | 
| 705 |  | 
| 706 | if(mol == NULL) { | 
| 707 | sprintf(painCave.errMsg, | 
| 708 | "RestReader Error: molecule not found on node %d\n", | 
| 709 | worldRank); | 
| 710 | painCave.isFatal = 1; | 
| 711 | simError(); | 
| 712 | } | 
| 713 |  | 
| 714 | nCurObj = mol->getNIntegrableObjects(); | 
| 715 | // send the number of integrableObjects in the molecule | 
| 716 | MPI_Send(&nCurObj, 1, MPI_INT, 0, | 
| 717 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD); | 
| 718 |  | 
| 719 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 720 | integrableObject != NULL; | 
| 721 | integrableObject = mol->nextIntegrableObject(ii)){ | 
| 722 | intObjIndexTransfer = integrableObject->getGlobalIntegrableObjectIndex(); | 
| 723 | // send the global index of the integrableObject | 
| 724 | MPI_Send(&intObjIndexTransfer, 1, MPI_INT, 0, | 
| 725 | TAKE_THIS_TAG_INT, MPI_COMM_WORLD); | 
| 726 | // listen for the value we want to set locally | 
| 727 | MPI_Recv(&angleTransfer, 1, MPI_REALTYPE, 0, | 
| 728 | TAKE_THIS_TAG_DOUBLE, MPI_COMM_WORLD, &istatus); | 
| 729 |  | 
| 730 | integrableObject->setZangle(angleTransfer); | 
| 731 | } | 
| 732 | } | 
| 733 | } | 
| 734 | } | 
| 735 | #endif | 
| 736 | } | 
| 737 |  | 
| 738 | void RestReader :: zeroZangle(){ | 
| 739 |  | 
| 740 | Molecule* mol; | 
| 741 | StuntDouble* integrableObject; | 
| 742 | SimInfo::MoleculeIterator mi; | 
| 743 | Molecule::IntegrableObjectIterator ii; | 
| 744 |  | 
| 745 | #ifndef IS_MPI | 
| 746 | // set all zAngles to 0.0 | 
| 747 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 748 | mol = info_->nextMolecule(mi)) | 
| 749 |  | 
| 750 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 751 | integrableObject != NULL; | 
| 752 | integrableObject = mol->nextIntegrableObject(ii)) | 
| 753 | integrableObject->setZangle( 0.0 ); | 
| 754 |  | 
| 755 |  | 
| 756 | // MPI Section of code.......... | 
| 757 | #else //IS_MPI | 
| 758 |  | 
| 759 | // first thing first, suspend fatalities. | 
| 760 | painCave.isEventLoop = 1; | 
| 761 |  | 
| 762 | int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone | 
| 763 | int haveError; | 
| 764 | int which_node; | 
| 765 |  | 
| 766 | haveError = 0; | 
| 767 |  | 
| 768 | for (int i=0 ; i < info_->getNGlobalMolecules(); i++) { | 
| 769 |  | 
| 770 | // Get the Node number which has this atom | 
| 771 | which_node = info_->getMolToProc(i); | 
| 772 |  | 
| 773 | // each processor zeroes its own integrable objects | 
| 774 | if (which_node == worldRank) { | 
| 775 | mol = info_->getMoleculeByGlobalIndex(i); | 
| 776 |  | 
| 777 | if(mol == NULL) { | 
| 778 | sprintf( painCave.errMsg, | 
| 779 | "Molecule not found on node %i!", | 
| 780 | which_node ); | 
| 781 | haveError = 1; | 
| 782 | simError(); | 
| 783 | } | 
| 784 |  | 
| 785 | for (integrableObject = mol->beginIntegrableObject(ii); | 
| 786 | integrableObject != NULL; | 
| 787 | integrableObject = mol->nextIntegrableObject(ii)){ | 
| 788 |  | 
| 789 | integrableObject->setZangle( 0.0 ); | 
| 790 |  | 
| 791 | } | 
| 792 | } | 
| 793 | } | 
| 794 |  | 
| 795 | #endif | 
| 796 | } | 
| 797 |  | 
| 798 | } // end namespace oopse |