Fix 1417

regtest/basic/rt-make-settimestep-rounding/COLVAR.reference

@@ -0,0 +1,4 @@
+#! FIELDS time c
+ 20000.000000 1.000000
+#! FIELDS time c
+ 20000.000000 1.000000

regtest/basic/rt-make-settimestep-rounding/Makefile

@@ -0,0 +1 @@
+include ../../scripts/test.make

regtest/basic/rt-make-settimestep-rounding/config

@@ -0,0 +1 @@
+type=make

regtest/basic/rt-make-settimestep-rounding/main.cpp

@@ -0,0 +1,53 @@
+#include <vector>
+#include <fstream>
+#include "plumed/wrapper/Plumed.h"
+#include <iostream>
+
+using namespace PLMD;
+
+template<typename T>
+void run(){
+
+  auto natoms=10;
+  std::vector<T> masses;
+  std::vector<std::array<T,3>> positions;
+  std::vector<std::array<T,3>> forces;
+  T box[3][3];
+  T virial[3][3];
+
+  Plumed p;
+
+  p.cmd("setRealPrecision",int(sizeof(T)));
+  p.cmd("setNatoms",natoms);
+  p.cmd("setTimestep",(T)0.002);
+  p.cmd("init");
+
+  p.cmd("readInputLines",
+    "c: CONSTANT VALUE=1.0 \n"
+    "PRINT ARG=c FILE=COLVAR RESTART=YES\n"
+  );
+
+  // dummy settings, not really used
+  positions.resize(natoms);
+  masses.resize(natoms);
+  forces.resize(natoms);
+
+  for(unsigned i=0;i<natoms;i++) for(unsigned j=0;j<3;j++) positions[i][j]=0.0;
+  for(unsigned i=0;i<natoms;i++) for(unsigned j=0;j<3;j++) forces[i][j]=0.0;
+  for(unsigned i=0;i<natoms;i++) masses[i]=1.0;
+
+  p.cmd("setStep",10000000);
+  p.cmd("setMasses",&masses[0]);
+  p.cmd("setPositions",&positions[0][0]);
+  p.cmd("setForces",&forces[0][0]);
+
+  p.cmd("setBox",&box[0][0]);
+  p.cmd("setVirial",&virial[0][0]);
+  p.cmd("calc");
+
+}
+
+int main(){
+  run<double>();
+  run<float>();
+}

src/core/Action.cpp

@@ -137,10 +137,7 @@ Action::Action(const ActionOptions&ao):
 }
 
 void Action::resetStoredTimestep() {
-  ActionWithValue* ts = plumed.getActionSet().selectWithLabel<ActionWithValue*>("timestep");
-  if( ts ) {
-    timestep = (ts->copyOutput(0))->get();
-  }
+  timestep = plumed.getTimestep();
 }
 
 Action::~Action() {

src/core/PlumedMain.cpp

@@ -737,14 +737,20 @@ void PlumedMain::cmd(std::string_view word,const TypesafePtr & val) {
         if( !ts->setValuePointer("timestep", val ) ) {
           plumed_error();
         }
+        // store a cached double precision value
+        if( ts ) {
+          timestep = (ts->copyOutput(0))->get();
+        }
+
         // The following is to avoid extra digits in case the MD code uses floats
         // e.g.: float f=0.002 when converted to double becomes 0.002000000094995
         // To avoid this, we keep only up to 6 significant digits after first one
         if( getRealPrecision()<=4 ) {
           Value* tstepv = ts->copyOutput(0);
           double magnitude=std::pow(10,std::floor(std::log10(tstepv->get())));
-          tstepv->set( std::round(tstepv->get()/magnitude*1e6)/1e6*magnitude );
+          timestep = std::round(tstepv->get()/magnitude*1e6)/1e6*magnitude;
         }
+
         ts->updateUnits( passtools.get() );
       }
       break;

src/core/PlumedMain.h

@@ -231,6 +231,8 @@ class PlumedMain:
 /// Flag for parse only mode -- basically just forces restart to turn off
   bool doParseOnly=false;
 
+/// Timestep saved as double
+  double timestep=0.0;
 private:
 /// Forward declaration.
   ForwardDecl<TypesafePtr> stopFlag_fwd;
@@ -533,6 +535,8 @@ class PlumedMain:
   double MDQuantityToPLUMED( const std::string& unit, const TypesafePtr & m) const ;
 /// Get the keywords for a particular action
   void getKeywordsForAction( const std::string& action, Keywords& keys ) const ;
+/// Return timestep
+  double getTimestep() const ;
 };
 
 /////
@@ -618,6 +622,10 @@ bool PlumedMain::callErrorHandler(int code,const char* msg)const {
   }
 }
 
+inline
+double PlumedMain::getTimestep() const {
+  return timestep;
+}
 
 }