Start of transfer of functionality from hack-the-tree to maser

regtest/basic/rt-make-extracv/main.cpp

@@ -30,6 +30,13 @@ int main(){
   plumed->cmd("readInputLine","PRINT ARG=e3 FILE=extra3 STRIDE=3");
 
   std::ofstream ofs("output");
+  long int rank=-1; plumed->cmd("getDataRank d",&rank,1 ); 
+  ofs<<"distance rank "<<rank<<"\n";
+  std::vector<long int> shape(rank); plumed->cmd("getDataShape d", &shape[0],rank);
+  ofs<<"distance shape "<<shape[0]; int nvals = shape[0];
+  for(unsigned i=1;i<shape.size(); ++i) { nvals*=shape[i]; ofs<<" "<<shape[i]; }
+  ofs<<" totalvals "<<nvals<<"\n";
+  std::vector<double> dmem(nvals); plumed->cmd("setMemoryForData d", &dmem[0],nvals);
 
   for(int step=0;step<10;step++){
     double extracv=step;
@@ -101,6 +108,7 @@ int main(){
     for(auto & f:forces) ofs<<" "<<f;
     ofs<<"\n";
     plumed->cmd("update");
+    ofs<<"distance value: "<<dmem[0]<<"\n";
   }
 
   delete plumed;

regtest/basic/rt-make-extracv/output.reference

@@ -1,3 +1,5 @@
+distance rank 1
+distance shape 1 totalvals 1
 extracv3_needed: 1
 bias_pre: -4987.88
 extracvf_pre: 0
@@ -11,6 +13,7 @@ bias: 13.5
 extracvf: 0
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4987.38
 extracvf_pre: -1
@@ -24,6 +27,7 @@ bias: 14
 extracvf: -1
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4985.88
 extracvf_pre: -2
@@ -37,6 +41,7 @@ bias: 15.5
 extracvf: -2
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 1
 bias_pre: -4983.38
 extracvf_pre: -3
@@ -50,6 +55,7 @@ bias: 18
 extracvf: -3
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4979.88
 extracvf_pre: -4
@@ -63,6 +69,7 @@ bias: 21.5
 extracvf: -4
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4975.38
 extracvf_pre: -5
@@ -76,6 +83,7 @@ bias: 26
 extracvf: -5
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 1
 bias_pre: -4969.88
 extracvf_pre: -6
@@ -89,6 +97,7 @@ bias: 31.5
 extracvf: -6
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4963.38
 extracvf_pre: -7
@@ -102,6 +111,7 @@ bias: 38
 extracvf: -7
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 0
 bias_pre: -4955.88
 extracvf_pre: -8
@@ -115,6 +125,7 @@ bias: 45.5
 extracvf: -8
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615
 extracv3_needed: 1
 bias_pre: -4947.38
 extracvf_pre: -9
@@ -128,3 +139,4 @@ bias: 54
 extracvf: -9
 extracvf2: 0
 f: 3 3 3 -3 -3 -3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
+distance value: 5.19615

regtest/basic/rt-simplemd2/COLVAR.reference

@@ -1,4 +1,4 @@
-#! FIELDS time a d ene @3.bias @3.force2
+#! FIELDS time a d ene @2.bias @2.force2
  0.050000  1.049  1.154 -646.446 54.357  3.239
  0.100000  1.081  1.105 -526.972  0.091  1.094
  0.150000  1.105  1.053 -552.967  2.124  0.317

regtest/basic/rt73/act.reference

@@ -1,4 +1,4 @@
-activity at step 0: ++++-
-activity at step 2: +-+--
-activity at step 3: -+-+-
-activity at step 4: +-+--
+activity at step 0: -+++++++++++-
+activity at step 2: -++++++++-+--
+activity at step 3: -+++++++-+-+-
+activity at step 4: -++++++++-+--

src/adjmat/Sprint.cpp

@@ -212,29 +212,17 @@ void Sprint::calculate() {
 }
 
 void Sprint::apply() {
-  std::vector<Vector>&   f(modifyForces());
-  Tensor&           v(modifyVirial());
-  unsigned          nat=getNumberOfAtoms();
-
+  bool hasforce=false;
   std::vector<double> forces( 3*getNumberOfAtoms() + 9 );
+  std::vector<double> fforces( 3*getNumberOfAtoms() + 9, 0 );
   for(int i=0; i<getNumberOfComponents(); ++i) {
     if( getPntrToComponent(i)->applyForce( forces ) ) {
-      for(unsigned j=0; j<nat; ++j) {
-        f[j][0]+=forces[3*j+0];
-        f[j][1]+=forces[3*j+1];
-        f[j][2]+=forces[3*j+2];
-      }
-      v(0,0)+=forces[3*nat+0];
-      v(0,1)+=forces[3*nat+1];
-      v(0,2)+=forces[3*nat+2];
-      v(1,0)+=forces[3*nat+3];
-      v(1,1)+=forces[3*nat+4];
-      v(1,2)+=forces[3*nat+5];
-      v(2,0)+=forces[3*nat+6];
-      v(2,1)+=forces[3*nat+7];
-      v(2,2)+=forces[3*nat+8];
+      hasforce=true; for(unsigned j=0; j<fforces.size(); ++j) fforces[j] += forces[j];
     }
   }
+  if( hasforce ) {
+    unsigned ind=0; setForcesOnAtoms( fforces, ind );
+  }
 }
 
 }

src/analysis/Histogram.cpp

@@ -525,7 +525,7 @@ void Histogram::apply() {
     tbase += myvessels[i]->getNumberOfDerivatives();
   }
   // And set the final forces on the atoms
-  setForcesOnAtoms( finalForces );
+  unsigned ind=0; setForcesOnAtoms( finalForces, ind );
   // Reset everything for next regular loop
   in_apply=false;
 }

src/analysis/OutputPDBFile.cpp

@@ -25,7 +25,6 @@
 #include "core/ActionRegister.h"
 #include "core/PlumedMain.h"
 #include "core/ActionSet.h"
-#include "core/Atoms.h"
 #include "core/GenericMolInfo.h"
 #include "tools/PDB.h"
 
@@ -90,8 +89,8 @@ void OutputPDBFile::performAnalysis() {
     afile.printf("DESCRIPTION: analysis data from calculation done by %s at time %f \n",getLabel().c_str(),getTime() );
     if( dissimilaritiesWereSet() ) afile.printf("REMARK %s \n", getDissimilarityInstruction().c_str() );
     afile.printf(descr.c_str(),getWeight(j) ); getStoredData(j,false).transferDataToPDB( mypdb );
-    if( plumed.getAtoms().usingNaturalUnits() ) mypdb.print( 1.0, mymoldat, afile, fmt );
-    else mypdb.print( plumed.getAtoms().getUnits().getLength()/0.1, mymoldat, afile, fmt );
+    if( usingNaturalUnits() ) mypdb.print( 1.0, mymoldat, afile, fmt );
+    else mypdb.print( getUnits().getLength()/0.1, mymoldat, afile, fmt );
   }
   afile.close();
 }

src/analysis/ReadDissimilarityMatrix.cpp

@@ -27,6 +27,7 @@
 #include "core/ActionSet.h"
 #include "core/ActionRegister.h"
 #include "core/ActionSetup.h"
+#include "core/ActionForInterface.h"
 #include "tools/IFile.h"
 
 //+PLUMEDOC ANALYSIS READ_DISSIMILARITY_MATRIX
@@ -89,8 +90,9 @@ ReadDissimilarityMatrix::ReadDissimilarityMatrix( const ActionOptions& ao ):
 {
   setStride(1); // Set the stride equal to one to ensure we don't get stuck in an infinite loop
   std::vector<ActionSetup*> setupActions=plumed.getActionSet().select<ActionSetup*>();
-  if( my_input_data && (plumed.getActionSet().size()-setupActions.size())!=1 ) error("should only be this action and the READ_ANALYSIS_FRAMES command in the input file");
-  if( !my_input_data && plumed.getActionSet().size()!=0 ) error("read dissimilarity matrix command must be at top of input file");
+  std::vector<ActionForInterface*> interActions=plumed.getActionSet().select<ActionForInterface*>();
+  if( my_input_data && (plumed.getActionSet().size()-setupActions.size()-interActions.size())!=1 ) error("should only be this action and the READ_ANALYSIS_FRAMES command in the input file");
+  if( !my_input_data && plumed.getActionSet().size()!=interActions.size() ) error("read dissimilarity matrix command must be at top of input file");
 
   parse("FILE",fname);
   log.printf("  reading dissimilarity matrix from file %s \n",fname.c_str() );

src/bias/ExtendedLagrangian.cpp

@@ -23,7 +23,6 @@
 #include "ActionRegister.h"
 #include "tools/Random.h"
 #include "core/PlumedMain.h"
-#include "core/Atoms.h"
 
 namespace PLMD {
 namespace bias {
@@ -161,11 +160,7 @@ ExtendedLagrangian::ExtendedLagrangian(const ActionOptions&ao):
   parseVector("TAU",tau);
   parseVector("FRICTION",friction);
   parseVector("KAPPA",kappa);
-  double temp=-1.0;
-  parse("TEMP",temp);
-  if(temp>=0.0) kbt=plumed.getAtoms().getKBoltzmann()*temp;
-  else kbt=plumed.getAtoms().getKbT();
-  checkRead();
+  kbt=getkBT(); checkRead();
 
   log.printf("  with harmonic force constant");
   for(unsigned i=0; i<kappa.size(); i++) log.printf(" %f",kappa[i]);

src/bias/MaxEnt.cpp

@@ -21,7 +21,6 @@
 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
 #include "Bias.h"
 #include "core/PlumedMain.h"
-#include "core/Atoms.h"
 #include "core/ActionRegister.h"
 #include "core/ActionWithValue.h"
 #include "tools/Communicator.h"
@@ -260,10 +259,7 @@ MaxEnt::MaxEnt(const ActionOptions&ao):
   stride_=pace_;  //if no STRIDE is passed, then Lagrangian multipliers willbe printed at each update
   parse("PRINT_STRIDE",stride_);
   if(stride_<=0 ) error("frequency for Lagrangian multipliers printing (STRIDE) is nonsensical");
-  simtemp=0.;
-  parse("TEMP",simtemp);
-  if(simtemp>0) simtemp*=plumed.getAtoms().getKBoltzmann();
-  else simtemp=plumed.getAtoms().getKbT();
+  simtemp=getkBT();
   parseFlag("REWEIGHT",reweight);
   if(simtemp<=0 && reweight) error("Set the temperature (TEMP) if you want to do reweighting.");
 

src/bias/MetaD.cpp

@@ -23,14 +23,14 @@
 #include "ActionRegister.h"
 #include "core/ActionSet.h"
 #include "core/PlumedMain.h"
-#include "core/Atoms.h"
 #include "core/FlexibleBin.h"
 #include "tools/Exception.h"
 #include "tools/Grid.h"
 #include "tools/Matrix.h"
 #include "tools/OpenMP.h"
 #include "tools/Random.h"
 #include "tools/File.h"
+#include "tools/Communicator.h"
 #include <ctime>
 #include <numeric>
 #if defined(__PLUMED_HAS_GETCWD)
@@ -715,11 +715,7 @@ MetaD::MetaD(const ActionOptions& ao):
     parse("BIASFACTOR",biasf_);
   }
   if( biasf_<1.0  && biasf_!=-1.0) error("well tempered bias factor is nonsensical");
-  parse("DAMPFACTOR",dampfactor_);
-  double temp=0.0;
-  parse("TEMP",temp);
-  if(temp>0.0) kbt_=plumed.getAtoms().getKBoltzmann()*temp;
-  else kbt_=plumed.getAtoms().getKbT();
+  parse("DAMPFACTOR",dampfactor_); kbt_=getkBT();
   if(biasf_>=1.0) {
     if(kbt_==0.0) error("Unless the MD engine passes the temperature to plumed, with well-tempered metad you must specify it using TEMP");
     welltemp_=true;

src/bias/PBMetaD.cpp

@@ -23,14 +23,14 @@
 #include "ActionRegister.h"
 #include "core/ActionSet.h"
 #include "core/PlumedMain.h"
-#include "core/Atoms.h"
 #include "core/FlexibleBin.h"
 #include "tools/Exception.h"
 #include "tools/Grid.h"
 #include "tools/Matrix.h"
 #include "tools/OpenMP.h"
 #include "tools/Random.h"
 #include "tools/File.h"
+#include "tools/Communicator.h"
 #include <ctime>
 #include <numeric>
 #if defined(__PLUMED_HAS_GETCWD)
@@ -516,10 +516,7 @@ PBMetaD::PBMetaD(const ActionOptions& ao):
 
   parse("BIASFACTOR",biasf_);
   if( biasf_<1.0 ) error("well tempered bias factor is nonsensical");
-  double temp=0.0;
-  parse("TEMP",temp);
-  if(temp>0.0) kbt_=plumed.getAtoms().getKBoltzmann()*temp;
-  else kbt_=plumed.getAtoms().getKbT();
+  kbt_=getkBT();
   if(biasf_>1.0) {
     if(kbt_==0.0) error("Unless the MD engine passes the temperature to plumed, with well-tempered metad you must specify it using TEMP");
     welltemp_=true;

src/bias/ReweightBase.cpp

@@ -20,8 +20,6 @@
    along with plumed.  If not, see <http://www.gnu.org/licenses/>.
 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
 #include "ReweightBase.h"
-#include "core/PlumedMain.h"
-#include "core/Atoms.h"
 
 namespace PLMD {
 namespace bias {
@@ -40,9 +38,7 @@ ReweightBase::ReweightBase(const ActionOptions&ao):
   ActionWithValue(ao),
   ActionWithArguments(ao)
 {
-  simtemp=0.; parse("TEMP",simtemp);
-  if(simtemp>0) simtemp*=plumed.getAtoms().getKBoltzmann();
-  else simtemp=plumed.getAtoms().getKbT();
+  simtemp=getkBT();
   if(simtemp==0) error("The MD engine does not pass the temperature to plumed so you have to specify it using TEMP");
   // Create something to hold the weight
   addValue(); setNotPeriodic();

src/bias/ReweightTemperaturePressure.cpp

@@ -20,8 +20,6 @@
    along with plumed.  If not, see <http://www.gnu.org/licenses/>.
 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
 #include "core/ActionRegister.h"
-#include "core/PlumedMain.h"
-#include "core/Atoms.h"
 #include "ReweightBase.h"
 
 //+PLUMEDOC REWEIGHTING REWEIGHT_TEMP_PRESS
@@ -210,7 +208,7 @@ ReweightTemperaturePressure::ReweightTemperaturePressure(const ActionOptions&ao)
   parse("REWEIGHT_PRESSURE",rpress_);
   parse("PRESSURE",press_);
   parse("REWEIGHT_TEMP",rtemp_);
-  rtemp_*=plumed.getAtoms().getKBoltzmann();
+  rtemp_*=getKBoltzmann();
 
   parseArgumentList("ENERGY",myenergy);
   if(!myenergy.empty()) {
@@ -235,15 +233,15 @@ ReweightTemperaturePressure::ReweightTemperaturePressure(const ActionOptions&ao)
   // 4 possible cases
   // Case 1) Reweight from T to T' with V=const (canonical)
   if (rtemp_>=0 && press_<0 && rpress_<0 && !myenergy.empty() && myvol.empty() ) {
-    log.printf("  reweighting simulation from temperature %f to temperature %f at constant volume \n",simtemp/plumed.getAtoms().getKBoltzmann(),rtemp_/plumed.getAtoms().getKBoltzmann() );
+    log.printf("  reweighting simulation from temperature %f to temperature %f at constant volume \n",simtemp/getKBoltzmann(),rtemp_/getKBoltzmann() );
     log.printf("  WARNING: If the simulation is performed at constant pressure add the keywords PRESSURE and VOLUME \n" );
   }
   // Case 2) Reweight from T to T' with P=const (isothermal-isobaric)
-  else if (rtemp_>=0 && press_>=0 && rpress_<0 && !myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from temperature %f to temperature %f at constant pressure %f \n",simtemp/plumed.getAtoms().getKBoltzmann(),rtemp_/plumed.getAtoms().getKBoltzmann(), press_ );
+  else if (rtemp_>=0 && press_>=0 && rpress_<0 && !myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from temperature %f to temperature %f at constant pressure %f \n",simtemp/getKBoltzmann(),rtemp_/getKBoltzmann(), press_ );
   // Case 3) Reweight from P to P' with T=const (isothermal-isobaric)
-  else if (rtemp_<0 && press_>=0 && rpress_>=0 && myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from pressure %f to pressure %f at constant temperature %f\n",press_,rpress_,simtemp/plumed.getAtoms().getKBoltzmann() );
+  else if (rtemp_<0 && press_>=0 && rpress_>=0 && myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from pressure %f to pressure %f at constant temperature %f\n",press_,rpress_,simtemp/getKBoltzmann() );
   // Case 4) Reweight from T,P to T',P' (isothermal-isobaric)
-  else if (rtemp_>0 && press_>=0 && rpress_>=0 && !myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from temperature %f and pressure %f to temperature %f and pressure %f \n",simtemp/plumed.getAtoms().getKBoltzmann(), press_, rtemp_/plumed.getAtoms().getKBoltzmann(), rpress_);
+  else if (rtemp_>0 && press_>=0 && rpress_>=0 && !myenergy.empty() && !myvol.empty() ) log.printf("  reweighting simulation from temperature %f and pressure %f to temperature %f and pressure %f \n",simtemp/getKBoltzmann(), press_, rtemp_/getKBoltzmann(), rpress_);
   else error("Combination of ENERGY, VOLUME, REWEIGHT_PRESSURE, PRESSURE and REWEIGHT_TEMP not supported. Please refer to the manual for supported combinations.");
 }
 

src/cltools/Driver.cpp

@@ -727,7 +727,7 @@ int Driver<real>::main(FILE* in,FILE*out,Communicator& pc) {
         std::map<std::string,std::vector<std::string> > data;
         IFile ifile; ifile.open(plumedFile); std::vector<std::string> words;
         while( Tools::getParsedLine(ifile,words) && !p.getEndPlumed() ) {
-          p.readInputWords(words); Action* aa=p.getActionSet()[p.getActionSet().size()-1].get();
+          p.readInputWords(words,false); Action* aa=p.getActionSet()[p.getActionSet().size()-1].get();
           ActionWithValue* av=dynamic_cast<ActionWithValue*>(aa);
           if( av && aa->getDefaultString().length()>0 ) {
             std::vector<std::string> def; def.push_back( "defaults " + aa->getDefaultString() );

src/cltools/GenExample.cpp

@@ -29,6 +29,7 @@
 #include "core/ActionShortcut.h"
 #include "core/ActionSet.h"
 #include "core/PlumedMain.h"
+#include "tools/Communicator.h"
 #include "tools/IFile.h"
 #include <cstdio>
 #include <string>

src/cltools/SimpleMD.cpp

@@ -493,8 +493,8 @@ class SimpleMD:
     randomize_velocities(natoms,ndim,temperature,masses,velocities,random);
 
     if(plumed) {
-      plumed->cmd("setNoVirial");
       plumed->cmd("setNatoms",natoms);
+      plumed->cmd("setNoVirial");
       plumed->cmd("setMDEngine","simpleMD");
       plumed->cmd("setTimestep",tstep);
       plumed->cmd("setPlumedDat","plumed.dat");