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

Comparing branches/development/src/rnemd/RNEMD.cpp (file contents):
Revision 1867 by gezelter, Mon Apr 29 17:53:48 2013 UTC vs.
Revision 1874 by gezelter, Wed May 15 15:09:35 2013 UTC

# Line 71 | Line 71 | namespace OpenMD {
71    RNEMD::RNEMD(SimInfo* info) : info_(info), evaluator_(info), seleMan_(info),
72                                  evaluatorA_(info), seleManA_(info),
73                                  commonA_(info), evaluatorB_(info),
74 <                                seleManB_(info), commonB_(info),
74 >                                seleManB_(info), commonB_(info),
75 >                                hasData_(false), hasDividingArea_(false),
76                                  usePeriodicBoundaryConditions_(info->getSimParams()->getUsePeriodicBoundaryConditions()) {
77  
78      trialCount_ = 0;
# Line 603 | Line 604 | namespace OpenMD {
604   #ifdef IS_MPI
605      }
606   #endif
607 +
608 +    // delete all of the objects we created:
609 +    delete areaAccumulator_;    
610 +    data_.clear();
611    }
612    
613    void RNEMD::doSwap(SelectionManager& smanA, SelectionManager& smanB) {
# Line 1146 | Line 1151 | namespace OpenMD {
1151            //if w is in the right range, so should be x, y, z.
1152            vector<StuntDouble*>::iterator sdi;
1153            Vector3d vel;
1154 <          for (sdi = coldBin.begin(); sdi != coldBin.end(); sdi++) {
1154 >          for (sdi = coldBin.begin(); sdi != coldBin.end(); ++sdi) {
1155              if (rnemdFluxType_ == rnemdFullKE) {
1156                vel = (*sdi)->getVel() * c;
1157                (*sdi)->setVel(vel);
# Line 1157 | Line 1162 | namespace OpenMD {
1162              }
1163            }
1164            w = sqrt(w);
1165 <          for (sdi = hotBin.begin(); sdi != hotBin.end(); sdi++) {
1165 >          for (sdi = hotBin.begin(); sdi != hotBin.end(); ++sdi) {
1166              if (rnemdFluxType_ == rnemdFullKE) {
1167                vel = (*sdi)->getVel();
1168                vel.x() *= x;
# Line 1276 | Line 1281 | namespace OpenMD {
1281        vector<RealType>::iterator ri;
1282        RealType r1, r2, alpha0;
1283        vector<pair<RealType,RealType> > rps;
1284 <      for (ri = realRoots.begin(); ri !=realRoots.end(); ri++) {
1284 >      for (ri = realRoots.begin(); ri !=realRoots.end(); ++ri) {
1285          r2 = *ri;
1286          //check if FindRealRoots() give the right answer
1287          if ( fabs(u0 + r2 * (u1 + r2 * (u2 + r2 * (u3 + r2 * u4)))) > 1e-6 ) {
# Line 1308 | Line 1313 | namespace OpenMD {
1313          RealType diff;
1314          pair<RealType,RealType> bestPair = make_pair(1.0, 1.0);
1315          vector<pair<RealType,RealType> >::iterator rpi;
1316 <        for (rpi = rps.begin(); rpi != rps.end(); rpi++) {
1316 >        for (rpi = rps.begin(); rpi != rps.end(); ++rpi) {
1317            r1 = (*rpi).first;
1318            r2 = (*rpi).second;
1319            switch(rnemdFluxType_) {
# Line 1375 | Line 1380 | namespace OpenMD {
1380          }
1381          vector<StuntDouble*>::iterator sdi;
1382          Vector3d vel;
1383 <        for (sdi = coldBin.begin(); sdi != coldBin.end(); sdi++) {
1383 >        for (sdi = coldBin.begin(); sdi != coldBin.end(); ++sdi) {
1384            vel = (*sdi)->getVel();
1385            vel.x() *= x;
1386            vel.y() *= y;
# Line 1386 | Line 1391 | namespace OpenMD {
1391          x = 1.0 + px * (1.0 - x);
1392          y = 1.0 + py * (1.0 - y);
1393          z = 1.0 + pz * (1.0 - z);
1394 <        for (sdi = hotBin.begin(); sdi != hotBin.end(); sdi++) {
1394 >        for (sdi = hotBin.begin(); sdi != hotBin.end(); ++sdi) {
1395            vel = (*sdi)->getVel();
1396            vel.x() *= x;
1397            vel.y() *= y;
# Line 1664 | Line 1669 | namespace OpenMD {
1669                    Vector3d vel;
1670                    Vector3d rPos;
1671                    
1672 <                  for (sdi = coldBin.begin(); sdi != coldBin.end(); sdi++) {
1672 >                  for (sdi = coldBin.begin(); sdi != coldBin.end(); ++sdi) {
1673                      //vel = (*sdi)->getVel();
1674                      rPos = (*sdi)->getPos() - coordinateOrigin_;
1675                      if (doLinearPart)
# Line 1680 | Line 1685 | namespace OpenMD {
1685                        }
1686                      }
1687                    }
1688 <                  for (sdi = hotBin.begin(); sdi != hotBin.end(); sdi++) {
1688 >                  for (sdi = hotBin.begin(); sdi != hotBin.end(); ++sdi) {
1689                      //vel = (*sdi)->getVel();
1690                      rPos = (*sdi)->getPos() - coordinateOrigin_;
1691                      if (doLinearPart)
# Line 2008 | Line 2013 | namespace OpenMD {
2013          }
2014        }
2015      }
2016 +    hasData_ = true;
2017    }
2018  
2019    void RNEMD::getStarted() {
# Line 2040 | Line 2046 | namespace OpenMD {
2046    
2047    void RNEMD::writeOutputFile() {
2048      if (!doRNEMD_) return;
2049 +    if (!hasData_) return;
2050      
2051   #ifdef IS_MPI
2052      // If we're the root node, should we print out the results
# Line 2061 | Line 2068 | namespace OpenMD {
2068        RealType time = currentSnap_->getTime();
2069        RealType avgArea;
2070        areaAccumulator_->getAverage(avgArea);
2064 –      RealType Jz = kineticExchange_ / (time * avgArea)
2065 –        / PhysicalConstants::energyConvert;
2066 –      Vector3d JzP = momentumExchange_ / (time * avgArea);      
2067 –      Vector3d JzL = angularMomentumExchange_ / (time * avgArea);      
2071  
2072 +      RealType Jz(0.0);
2073 +      Vector3d JzP(V3Zero);
2074 +      Vector3d JzL(V3Zero);
2075 +      if (time >= info_->getSimParams()->getDt()) {
2076 +        Jz = kineticExchange_ / (time * avgArea)
2077 +          / PhysicalConstants::energyConvert;
2078 +        JzP = momentumExchange_ / (time * avgArea);
2079 +        JzL = angularMomentumExchange_ / (time * avgArea);
2080 +      }
2081 +
2082        rnemdFile_ << "#######################################################\n";
2083        rnemdFile_ << "# RNEMD {\n";
2084  
# Line 2218 | Line 2231 | namespace OpenMD {
2231        rnemdFile_ << "\t" << s;
2232      } else{
2233        sprintf( painCave.errMsg,
2234 <               "RNEMD detected a numerical error writing: %s for bin %d",
2234 >               "RNEMD detected a numerical error writing: %s for bin %u",
2235                 data_[index].title.c_str(), bin);
2236        painCave.isFatal = 1;
2237        simError();
# Line 2241 | Line 2254 | namespace OpenMD {
2254          isinf(s[1]) || isnan(s[1]) ||
2255          isinf(s[2]) || isnan(s[2]) ) {      
2256        sprintf( painCave.errMsg,
2257 <               "RNEMD detected a numerical error writing: %s for bin %d",
2257 >               "RNEMD detected a numerical error writing: %s for bin %u",
2258                 data_[index].title.c_str(), bin);
2259        painCave.isFatal = 1;
2260        simError();
# Line 2266 | Line 2279 | namespace OpenMD {
2279        rnemdFile_ << "\t" << s;
2280      } else{
2281        sprintf( painCave.errMsg,
2282 <               "RNEMD detected a numerical error writing: %s std. dev. for bin %d",
2282 >               "RNEMD detected a numerical error writing: %s std. dev. for bin %u",
2283                 data_[index].title.c_str(), bin);
2284        painCave.isFatal = 1;
2285        simError();
# Line 2288 | Line 2301 | namespace OpenMD {
2301          isinf(s[1]) || isnan(s[1]) ||
2302          isinf(s[2]) || isnan(s[2]) ) {      
2303        sprintf( painCave.errMsg,
2304 <               "RNEMD detected a numerical error writing: %s std. dev. for bin %d",
2304 >               "RNEMD detected a numerical error writing: %s std. dev. for bin %u",
2305                 data_[index].title.c_str(), bin);
2306        painCave.isFatal = 1;
2307        simError();

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines