Line data    Source code

       1             : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       2             :    Copyright (c) 2012-2023 The plumed team
       3             :    (see the PEOPLE file at the root of the distribution for a list of names)
       4             : 
       5             :    See http://www.plumed.org for more information.
       6             : 
       7             :    This file is part of plumed, version 2.
       8             : 
       9             :    plumed is free software: you can redistribute it and/or modify
      10             :    it under the terms of the GNU Lesser General Public License as published by
      11             :    the Free Software Foundation, either version 3 of the License, or
      12             :    (at your option) any later version.
      13             : 
      14             :    plumed is distributed in the hope that it will be useful,
      15             :    but WITHOUT ANY WARRANTY; without even the implied warranty of
      16             :    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
      17             :    GNU Lesser General Public License for more details.
      18             : 
      19             :    You should have received a copy of the GNU Lesser General Public License
      20             :    along with plumed.  If not, see <http://www.gnu.org/licenses/>.
      21             : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
      22             : #include "Bias.h"
      23             : #include "core/ActionRegister.h"
      24             : #include "tools/Grid.h"
      25             : #include "tools/Exception.h"
      26             : #include "tools/File.h"
      27             : 
      28             : namespace PLMD {
      29             : namespace bias {
      30             : 
      31             : //+PLUMEDOC BIAS EXTERNAL
      32             : /*
      33             : Calculate a restraint that is defined on a grid that is read during start up
      34             : 
      35             : In this action you define an external potential that you would like to use as a bias in a file.
      36             : This file contains the value of your bias potential on a grid of points.  The instantaneous value
      37             : of the bias potential is then computed via interpolation.
      38             : 
      39             : 
      40             : The following input illustrates how this works in practice.  Here the external potential is
      41             : defined in the file `extras/1d_bias.dat`. This potential acts on the distance between atoms 3 and 5.
      42             : 
      43             : ```plumed
      44             : #SETTINGS INPUTFILES=extras/1d_bias.grid
      45             : d1: DISTANCE ATOMS=3,5
      46             : m: EXTERNAL ARG=d1 FILE=extras/1d_bias.grid
      47             : ```
      48             : 
      49             : As you can see from the example above, the file that is input to the EXTERNAL command
      50             : contains the value of the function and its derivative on a grid of points.  This file
      51             : has the grid format that is discussed in the documentation for [DUMPGRID](DUMPGRID.md).
      52             : 
      53             : You can also include grids that are a function of more than one collective
      54             : variable.  For instance the following would be the input for an external
      55             : potential acting on two torsional angles:
      56             : 
      57             : ```plumed
      58             : #SETTINGS INPUTFILES=extras/2d_bias.grid
      59             : t1: TORSION ATOMS=4,5,6,7
      60             : t2: TORSION ATOMS=6,7,8,9
      61             : ext: EXTERNAL ARG=t1,t2 FILE=extras/2d_bias.grid
      62             : ```
      63             : 
      64             : Please note the order that the order of arguments in the plumed.dat file must be the same as
      65             : the order of arguments in the header of the grid file.
      66             : */
      67             : //+ENDPLUMEDOC
      68             : 
      69             : class External : public Bias {
      70             : 
      71             : private:
      72             :   std::unique_ptr<GridBase> BiasGrid_;
      73             :   double scale_;
      74             : 
      75             : public:
      76             :   explicit External(const ActionOptions&);
      77             :   void calculate() override;
      78             :   static void registerKeywords(Keywords& keys);
      79             : };
      80             : 
      81             : PLUMED_REGISTER_ACTION(External,"EXTERNAL")
      82             : 
      83           7 : void External::registerKeywords(Keywords& keys) {
      84           7 :   Bias::registerKeywords(keys);
      85           7 :   keys.add("compulsory","FILE","the name of the file containing the external potential.");
      86           7 :   keys.addFlag("NOSPLINE",false,"specifies that no spline interpolation is to be used when calculating the energy and forces due to the external potential");
      87           7 :   keys.addFlag("SPARSE",false,"specifies that the external potential uses a sparse grid");
      88           7 :   keys.add("compulsory","SCALE","1.0","a factor that multiplies the external potential, useful to invert free energies");
      89           7 : }
      90             : 
      91           3 : External::External(const ActionOptions& ao):
      92           3 :   PLUMED_BIAS_INIT(ao) {
      93             :   std::string filename;
      94           6 :   parse("FILE",filename);
      95           3 :   if( filename.length()==0 ) {
      96           0 :     error("No external potential file was specified");
      97             :   }
      98           3 :   bool sparsegrid=false;
      99           3 :   parseFlag("SPARSE",sparsegrid);
     100           3 :   bool nospline=false;
     101           3 :   parseFlag("NOSPLINE",nospline);
     102           3 :   bool spline=!nospline;
     103           5 :   parse("SCALE",scale_);
     104             : 
     105           3 :   checkRead();
     106             : 
     107           3 :   log.printf("  External potential from file %s\n",filename.c_str());
     108           3 :   log.printf("  Multiplied by %lf\n",scale_);
     109           3 :   if(spline) {
     110           3 :     log.printf("  External potential uses spline interpolation\n");
     111             :   }
     112           3 :   if(sparsegrid) {
     113           0 :     log.printf("  External potential uses sparse grid\n");
     114             :   }
     115             : 
     116             : // read grid
     117           3 :   IFile gridfile;
     118           3 :   gridfile.open(filename);
     119           1 :   std::string funcl=getLabel() + ".bias";
     120           2 :   BiasGrid_=GridBase::create(funcl,getArguments(),gridfile,sparsegrid,spline,true);
     121           1 :   if(BiasGrid_->getDimension()!=getNumberOfArguments()) {
     122           0 :     error("mismatch between dimensionality of input grid and number of arguments");
     123             :   }
     124           3 :   for(unsigned i=0; i<getNumberOfArguments(); ++i) {
     125           4 :     if( getPntrToArgument(i)->isPeriodic()!=BiasGrid_->getIsPeriodic()[i] ) {
     126           0 :       error("periodicity mismatch between arguments and input bias");
     127             :     }
     128             :   }
     129          10 : }
     130             : 
     131           5 : void External::calculate() {
     132           5 :   unsigned ncv=getNumberOfArguments();
     133           5 :   std::vector<double> cv(ncv), der(ncv);
     134             : 
     135          15 :   for(unsigned i=0; i<ncv; ++i) {
     136          10 :     cv[i]=getArgument(i);
     137             :   }
     138             : 
     139           5 :   double ene=scale_*BiasGrid_->getValueAndDerivatives(cv,der);
     140             : 
     141           5 :   setBias(ene);
     142             : 
     143          15 :   for(unsigned i=0; i<ncv; ++i) {
     144          10 :     const double f=-scale_*der[i];
     145          10 :     setOutputForce(i,f);
     146             :   }
     147           5 : }
     148             : 
     149             : }
     150             : }