ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/io/DumpWriter.cpp
(Generate patch)

Comparing:
trunk/src/io/DumpWriter.cpp (file contents), Revision 1437 by gezelter, Wed Apr 21 14:59:18 2010 UTC vs.
branches/development/src/io/DumpWriter.cpp (file contents), Revision 1764 by gezelter, Tue Jul 3 18:32:27 2012 UTC

# Line 36 | Line 36
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, 24107 (2008).          
39 < * [4]  Vardeman & Gezelter, in progress (2009).                        
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   #include "io/DumpWriter.hpp"
# Line 49 | Line 50
50  
51   #ifdef IS_MPI
52   #include <mpi.h>
53 < #endif //is_mpi
53 > #endif
54  
55   using namespace std;
56   namespace OpenMD {
# Line 58 | Line 59 | namespace OpenMD {
59      : info_(info), filename_(info->getDumpFileName()), eorFilename_(info->getFinalConfigFileName()){
60  
61      Globals* simParams = info->getSimParams();
62 <    needCompression_ = simParams->getCompressDumpFile();
63 <    needForceVector_ = simParams->getOutputForceVector();
62 >    needCompression_   = simParams->getCompressDumpFile();
63 >    needForceVector_   = simParams->getOutputForceVector();
64 >    needParticlePot_   = simParams->getOutputParticlePotential();
65 >    needFlucQ_         = simParams->getOutputFluctuatingCharges();
66 >    needElectricField_ = simParams->getOutputElectricField();
67 >
68 >    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
69 >      doSiteData_ = true;
70 >    } else {
71 >      doSiteData_ = false;
72 >    }
73 >
74      createDumpFile_ = true;
75   #ifdef HAVE_LIBZ
76      if (needCompression_) {
# Line 97 | Line 108 | namespace OpenMD {
108      Globals* simParams = info->getSimParams();
109      eorFilename_ = filename_.substr(0, filename_.rfind(".")) + ".eor";    
110  
111 <    needCompression_ = simParams->getCompressDumpFile();
112 <    needForceVector_ = simParams->getOutputForceVector();
111 >    needCompression_   = simParams->getCompressDumpFile();
112 >    needForceVector_   = simParams->getOutputForceVector();
113 >    needParticlePot_   = simParams->getOutputParticlePotential();
114 >    needFlucQ_         = simParams->getOutputFluctuatingCharges();
115 >    needElectricField_ = simParams->getOutputElectricField();
116 >
117 >    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
118 >      doSiteData_ = true;
119 >    } else {
120 >      doSiteData_ = false;
121 >    }
122 >
123      createDumpFile_ = true;
124   #ifdef HAVE_LIBZ
125      if (needCompression_) {
# Line 136 | Line 157 | namespace OpenMD {
157      Globals* simParams = info->getSimParams();
158      eorFilename_ = filename_.substr(0, filename_.rfind(".")) + ".eor";    
159      
160 <    needCompression_ = simParams->getCompressDumpFile();
161 <    needForceVector_ = simParams->getOutputForceVector();
162 <    
160 >    needCompression_   = simParams->getCompressDumpFile();
161 >    needForceVector_   = simParams->getOutputForceVector();
162 >    needParticlePot_   = simParams->getOutputParticlePotential();
163 >    needFlucQ_         = simParams->getOutputFluctuatingCharges();
164 >    needElectricField_ = simParams->getOutputElectricField();
165 >
166 >    if (needParticlePot_ || needFlucQ_ || needElectricField_) {
167 >      doSiteData_ = true;
168 >    } else {
169 >      doSiteData_ = false;
170 >    }
171 >
172   #ifdef HAVE_LIBZ
173      if (needCompression_) {
174        filename_ += ".gz";
# Line 227 | Line 257 | namespace OpenMD {
257              hmat(0, 2), hmat(1, 2), hmat(2, 2));
258      os << buffer;
259  
260 <    RealType chi = s->getChi();
261 <    RealType integralOfChiDt = s->getIntegralOfChiDt();
262 <    if (isinf(chi) || isnan(chi) ||
263 <        isinf(integralOfChiDt) || isnan(integralOfChiDt)) {      
260 >    pair<RealType, RealType> thermostat = s->getThermostat();
261 >
262 >    if (isinf(thermostat.first)  || isnan(thermostat.first) ||
263 >        isinf(thermostat.second) || isnan(thermostat.second)) {      
264        sprintf( painCave.errMsg,
265                 "DumpWriter detected a numerical error writing the thermostat");
266        painCave.isFatal = 1;
267        simError();
268      }
269 <    sprintf(buffer, "  Thermostat: %.10g , %.10g\n", chi, integralOfChiDt);
269 >    sprintf(buffer, "  Thermostat: %.10g , %.10g\n", thermostat.first,
270 >            thermostat.second);
271      os << buffer;
272  
273      Mat3x3d eta;
274 <    eta = s->getEta();
274 >    eta = s->getBarostat();
275  
276      for (unsigned int i = 0; i < 3; i++) {
277        for (unsigned int j = 0; j < 3; j++) {
# Line 272 | Line 303 | namespace OpenMD {
303      StuntDouble* integrableObject;
304      SimInfo::MoleculeIterator mi;
305      Molecule::IntegrableObjectIterator ii;
306 +    RigidBody::AtomIterator ai;
307 +    Atom* atom;
308  
309   #ifndef IS_MPI
310      os << "  <Snapshot>\n";
# Line 289 | Line 322 | namespace OpenMD {
322        }
323      }    
324      os << "    </StuntDoubles>\n";
325 <    
325 >
326 >    if (doSiteData_) {
327 >      os << "    <SiteData>\n";
328 >      for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) {
329 >              
330 >        for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL;  
331 >           integrableObject = mol->nextIntegrableObject(ii)) {  
332 >
333 >          int ioIndex = integrableObject->getGlobalIntegrableObjectIndex();
334 >          // do one for the IO itself
335 >          os << prepareSiteLine(integrableObject, ioIndex, 0);
336 >
337 >          if (integrableObject->isRigidBody()) {
338 >            
339 >            RigidBody* rb = static_cast<RigidBody*>(integrableObject);
340 >            int siteIndex = 0;
341 >            for (atom = rb->beginAtom(ai); atom != NULL;  
342 >                 atom = rb->nextAtom(ai)) {                                            
343 >              os << prepareSiteLine(atom, ioIndex, siteIndex);
344 >              siteIndex++;
345 >            }
346 >          }
347 >        }
348 >      }    
349 >      os << "    </SiteData>\n";
350 >    }
351      os << "  </Snapshot>\n";
352  
353      os.flush();
# Line 306 | Line 364 | namespace OpenMD {
364      }
365      
366      const int masterNode = 0;
367 <
367 >    int nProc;
368 >    MPI_Comm_size(MPI_COMM_WORLD, &nProc);
369      if (worldRank == masterNode) {      
370        os << "  <Snapshot>\n";  
371        writeFrameProperties(os, info_->getSnapshotManager()->getCurrentSnapshot());
# Line 314 | Line 373 | namespace OpenMD {
373          
374        os << buffer;
375  
317      int nProc;
318      MPI_Comm_size(MPI_COMM_WORLD, &nProc);
376        for (int i = 1; i < nProc; ++i) {
377  
378          // receive the length of the string buffer that was
379          // prepared by processor i
380  
381 +        MPI_Bcast(&i, 1, MPI_INT,masterNode,MPI_COMM_WORLD);
382          int recvLength;
383          MPI_Recv(&recvLength, 1, MPI_INT, i, 0, MPI_COMM_WORLD, &istatus);
384          char* recvBuffer = new char[recvLength];
# Line 337 | Line 395 | namespace OpenMD {
395        os.flush();
396      } else {
397        int sendBufferLength = buffer.size() + 1;
398 <      MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
399 <      MPI_Send((void *)buffer.c_str(), sendBufferLength, MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
398 >      int myturn = 0;
399 >      for (int i = 1; i < nProc; ++i){
400 >        MPI_Bcast(&myturn,1, MPI_INT,masterNode,MPI_COMM_WORLD);
401 >        if (myturn == worldRank){
402 >          MPI_Send(&sendBufferLength, 1, MPI_INT, masterNode, 0, MPI_COMM_WORLD);
403 >          MPI_Send((void *)buffer.c_str(), sendBufferLength, MPI_CHAR, masterNode, 0, MPI_COMM_WORLD);
404 >        }
405 >      }
406      }
407  
408   #endif // is_mpi
# Line 420 | Line 484 | namespace OpenMD {
484  
485      if (needForceVector_) {
486        type += "f";
487 <      Vector3d frc;
424 <
425 <      frc = integrableObject->getFrc();
426 <
487 >      Vector3d frc = integrableObject->getFrc();
488        if (isinf(frc[0]) || isnan(frc[0]) ||
489            isinf(frc[1]) || isnan(frc[1]) ||
490            isinf(frc[2]) || isnan(frc[2]) ) {      
# Line 439 | Line 500 | namespace OpenMD {
500        
501        if (integrableObject->isDirectional()) {
502          type += "t";
503 <        Vector3d trq;
443 <        
444 <        trq = integrableObject->getTrq();
445 <        
503 >        Vector3d trq = integrableObject->getTrq();        
504          if (isinf(trq[0]) || isnan(trq[0]) ||
505              isinf(trq[1]) || isnan(trq[1]) ||
506              isinf(trq[2]) || isnan(trq[2]) ) {      
# Line 451 | Line 509 | namespace OpenMD {
509                     " for object %d", index);      
510            painCave.isFatal = 1;
511            simError();
512 <        }
455 <        
512 >        }        
513          sprintf(tempBuffer, " %13e %13e %13e",
514                  trq[0], trq[1], trq[2]);
515          line += tempBuffer;
516        }      
517      }
518 <    
518 >
519      sprintf(tempBuffer, "%10d %7s %s\n", index, type.c_str(), line.c_str());
520      return std::string(tempBuffer);
521    }
522  
523 +  std::string DumpWriter::prepareSiteLine(StuntDouble* integrableObject, int ioIndex, int siteIndex) {
524 +        
525 +
526 +    std::string id;
527 +    std::string type;
528 +    std::string line;
529 +    char tempBuffer[4096];
530 +
531 +    if (integrableObject->isRigidBody()) {
532 +      sprintf(tempBuffer, "%10d           ", ioIndex);
533 +      id = std::string(tempBuffer);
534 +    } else {
535 +      sprintf(tempBuffer, "%10d %10d", ioIndex, siteIndex);
536 +      id = std::string(tempBuffer);
537 +    }
538 +              
539 +    if (needFlucQ_) {
540 +      type += "cw";
541 +      RealType fqPos = integrableObject->getFlucQPos();
542 +      if (isinf(fqPos) || isnan(fqPos) ) {      
543 +        sprintf( painCave.errMsg,
544 +                 "DumpWriter detected a numerical error writing the"
545 +                 " fluctuating charge for object %s", id.c_str());      
546 +        painCave.isFatal = 1;
547 +        simError();
548 +      }
549 +      sprintf(tempBuffer, " %13e ", fqPos);
550 +      line += tempBuffer;
551 +    
552 +      RealType fqVel = integrableObject->getFlucQVel();
553 +      if (isinf(fqVel) || isnan(fqVel) ) {      
554 +        sprintf( painCave.errMsg,
555 +                 "DumpWriter detected a numerical error writing the"
556 +                 " fluctuating charge velocity for object %s", id.c_str());      
557 +        painCave.isFatal = 1;
558 +        simError();
559 +      }
560 +      sprintf(tempBuffer, " %13e ", fqVel);
561 +      line += tempBuffer;
562 +
563 +      if (needForceVector_) {
564 +        type += "g";
565 +        RealType fqFrc = integrableObject->getFlucQFrc();        
566 +        if (isinf(fqFrc) || isnan(fqFrc) ) {      
567 +          sprintf( painCave.errMsg,
568 +                   "DumpWriter detected a numerical error writing the"
569 +                   " fluctuating charge force for object %s", id.c_str());      
570 +          painCave.isFatal = 1;
571 +          simError();
572 +        }
573 +        sprintf(tempBuffer, " %13e ", fqFrc);        
574 +        line += tempBuffer;
575 +      }
576 +    }
577 +
578 +    if (needElectricField_) {
579 +      type += "e";
580 +      Vector3d eField= integrableObject->getElectricField();
581 +      if (isinf(eField[0]) || isnan(eField[0]) ||
582 +          isinf(eField[1]) || isnan(eField[1]) ||
583 +          isinf(eField[2]) || isnan(eField[2]) ) {      
584 +        sprintf( painCave.errMsg,
585 +                 "DumpWriter detected a numerical error writing the electric"
586 +                 " field for object %s", id.c_str());      
587 +        painCave.isFatal = 1;
588 +        simError();
589 +      }
590 +      sprintf(tempBuffer, " %13e %13e %13e",
591 +              eField[0], eField[1], eField[2]);
592 +      line += tempBuffer;
593 +    }
594 +
595 +
596 +    if (needParticlePot_) {
597 +      type += "u";
598 +      RealType particlePot = integrableObject->getParticlePot();
599 +      if (isinf(particlePot) || isnan(particlePot)) {      
600 +        sprintf( painCave.errMsg,
601 +                 "DumpWriter detected a numerical error writing the particle "
602 +                 " potential for object %s", id.c_str());      
603 +        painCave.isFatal = 1;
604 +        simError();
605 +      }
606 +      sprintf(tempBuffer, " %13e", particlePot);
607 +      line += tempBuffer;
608 +    }
609 +    
610 +
611 +    sprintf(tempBuffer, "%s %7s %s\n", id.c_str(), type.c_str(), line.c_str());
612 +    return std::string(tempBuffer);
613 +  }
614 +
615    void DumpWriter::writeDump() {
616      writeFrame(*dumpFile_);
617    }
# Line 540 | Line 689 | namespace OpenMD {
689      newOStream = new std::ofstream(filename.c_str());
690   #endif
691      //write out MetaData first
692 <    (*newOStream) << "<OpenMD version=1>" << std::endl;
692 >    (*newOStream) << "<OpenMD version=2>" << std::endl;
693      (*newOStream) << "  <MetaData>" << std::endl;
694      (*newOStream) << info_->getRawMetaData();
695      (*newOStream) << "  </MetaData>" << std::endl;

Comparing:
trunk/src/io/DumpWriter.cpp (property svn:keywords), Revision 1437 by gezelter, Wed Apr 21 14:59:18 2010 UTC vs.
branches/development/src/io/DumpWriter.cpp (property svn:keywords), Revision 1764 by gezelter, Tue Jul 3 18:32:27 2012 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines