Line data    Source code

       1             : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       2             :    Copyright (c) 2013-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 "PathBase.h"
      23             : #include "core/ActionRegister.h"
      24             : 
      25             : //+PLUMEDOC COLVAR PATH
      26             : /*
      27             : Path collective variables with a more flexible framework for the distance metric being used.
      28             : 
      29             : The Path Collective Variables developed by Branduardi and co-workers \cite brand07 allow one
      30             : to compute the progress along a high-dimensional path and the distance from the high-dimensional
      31             : path.  The progress along the path (s) is computed using:
      32             : 
      33             : \f[
      34             : s = \frac{ \sum_{i=1}^N i \exp( -\lambda R[X - X_i] ) }{ \sum_{i=1}^N \exp( -\lambda R[X - X_i] ) }
      35             : \f]
      36             : 
      37             : while the distance from the path (z) is measured using:
      38             : 
      39             : \f[
      40             : z = -\frac{1}{\lambda} \ln\left[ \sum_{i=1}^N \exp( -\lambda R[X - X_i] ) \right]
      41             : \f]
      42             : 
      43             : In these expressions \f$N\f$ high-dimensional frames (\f$X_i\f$) are used to describe the path in the high-dimensional
      44             : space. The two expressions above are then functions of the distances from each of the high-dimensional frames \f$R[X - X_i]\f$.
      45             : Within PLUMED there are multiple ways to define the distance from a high-dimensional configuration.  You could calculate
      46             : the RMSD distance or you could calculate the amount by which a set of collective variables change.  As such this implementation
      47             : of the path CV allows one to use all the difference distance metrics that are discussed in \ref dists. This is as opposed to
      48             : the alternative implementation of path (\ref PATHMSD) which is a bit faster but which only allows one to use the RMSD distance.
      49             : 
      50             : The \f$s\f$ and \f$z\f$ variables are calculated using the above formulas by default.  However, there is an alternative method
      51             : of calculating these collective variables, which is detailed in \cite bernd-path.  This alternative method uses the tools of
      52             : geometry (as opposed to algebra, which is used in the equations above).  In this alternative formula the progress along the path
      53             : \f$s\f$ is calculated using:
      54             : 
      55             : \f[
      56             : s = i_2 + \textrm{sign}(i_2-i_1) \frac{ \sqrt{( \mathbf{v}_1\cdot\mathbf{v}_2 )^2 - |\mathbf{v}_3|^2(|\mathbf{v}_1|^2 - |\mathbf{v}_2|^2) } }{2|\mathbf{v}_3|^2} - \frac{\mathbf{v}_1\cdot\mathbf{v}_3 - |\mathbf{v}_3|^2}{2|\mathbf{v}_3|^2}
      57             : \f]
      58             : 
      59             : where \f$\mathbf{v}_1\f$ and \f$\mathbf{v}_3\f$ are the vectors connecting the current position to the closest and second closest node of the path,
      60             : respectfully and \f$i_1\f$ and \f$i_2\f$ are the projections of the closest and second closest frames of the path. \f$\mathbf{v}_2\f$, meanwhile, is the
      61             : vector connecting the closest frame to the second closest frame.  The distance from the path, \f$z\f$ is calculated using:
      62             : 
      63             : \f[
      64             : z = \sqrt{ \left[ |\mathbf{v}_1|^2 - |\mathbf{v}_2| \left( \frac{ \sqrt{( \mathbf{v}_1\cdot\mathbf{v}_2 )^2 - |\mathbf{v}_3|^2(|\mathbf{v}_1|^2 - |\mathbf{v}_2|^2) } }{2|\mathbf{v}_3|^2} - \frac{\mathbf{v}_1\cdot\mathbf{v}_3 - |\mathbf{v}_3|^2}{2|\mathbf{v}_3|^2} \right) \right]^2 }
      65             : \f]
      66             : 
      67             : The symbols here are as they were for \f$s\f$.  If you would like to use these equations to calculate \f$s\f$ and \f$z\f$ then you should use the GPATH flag.
      68             : The values of \f$s\f$ and \f$z\f$ can then be referenced using the gspath and gzpath labels.
      69             : 
      70             : \par Examples
      71             : 
      72             : In the example below the path is defined using RMSD distance from frames.
      73             : 
      74             : \plumedfile
      75             : p1: PATH REFERENCE=file.pdb TYPE=OPTIMAL LAMBDA=500.0
      76             : PRINT ARG=p1.spath,p1.zpath STRIDE=1 FILE=colvar FMT=%8.4f
      77             : \endplumedfile
      78             : 
      79             : The reference frames in the path are defined in the pdb file shown below.  In this frame
      80             : each configuration in the path is separated by a line containing just the word END.
      81             : 
      82             : \auxfile{file.pdb}
      83             : ATOM      1  CL  ALA     1      -3.171   0.295   2.045  1.00  1.00
      84             : ATOM      5  CLP ALA     1      -1.819  -0.143   1.679  1.00  1.00
      85             : ATOM      6  OL  ALA     1      -1.177  -0.889   2.401  1.00  1.00
      86             : ATOM      7  NL  ALA     1      -1.313   0.341   0.529  1.00  1.00
      87             : END
      88             : ATOM      1  CL  ALA     1      -3.175   0.365   2.024  1.00  1.00
      89             : ATOM      5  CLP ALA     1      -1.814  -0.106   1.685  1.00  1.00
      90             : ATOM      6  OL  ALA     1      -1.201  -0.849   2.425  1.00  1.00
      91             : ATOM      7  NL  ALA     1      -1.296   0.337   0.534  1.00  1.00
      92             : END
      93             : ATOM      1  CL  ALA     1      -2.990   0.383   2.277  1.00  1.00
      94             : ATOM      5  CLP ALA     1      -1.664  -0.085   1.831  1.00  1.00
      95             : ATOM      6  OL  ALA     1      -0.987  -0.835   2.533  1.00  1.00
      96             : ATOM      7  NL  ALA     1      -1.227   0.364   0.646  1.00  1.00
      97             : END
      98             : \endauxfile
      99             : 
     100             : In the example below the path is defined using the values of two torsional angles (t1 and t2).
     101             : In addition, the \f$s\f$ and \f$z\f$ are calculated using the geometric expressions described
     102             : above rather than the algebraic expressions that are used by default.
     103             : 
     104             : \plumedfile
     105             : t1: TORSION ATOMS=5,7,9,15
     106             : t2: TORSION ATOMS=7,9,15,17
     107             : pp: PATH TYPE=EUCLIDEAN REFERENCE=epath.pdb GPATH NOSPATH NOZPATH
     108             : PRINT ARG=pp.* FILE=colvar
     109             : \endplumedfile
     110             : 
     111             : Notice that the LAMBDA parameter is not required here as we are not calculating \f$s\f$ and \f$s\f$
     112             : using the algebraic formulas defined earlier.  The positions of the frames in the path are defined
     113             : in the file epath.pdb.  An extract from this file looks as shown below.
     114             : 
     115             : \auxfile{epath.pdb}
     116             : REMARK ARG=t1,t2 t1=-4.25053  t2=3.88053
     117             : END
     118             : REMARK ARG=t1,t2 t1=-4.11     t2=3.75
     119             : END
     120             : REMARK ARG=t1,t2 t1=-3.96947  t2=3.61947
     121             : END
     122             : \endauxfile
     123             : 
     124             : The remarks in this pdb file tell PLUMED the labels that are being used to define the position in the
     125             : high dimensional space and the values that these arguments have at each point on the path.
     126             : 
     127             : */
     128             : //+ENDPLUMEDOC
     129             : 
     130             : namespace PLMD {
     131             : namespace mapping {
     132             : 
     133             : class Path : public PathBase {
     134             : public:
     135             :   static void registerKeywords( Keywords& keys );
     136             :   explicit Path(const ActionOptions&);
     137             : };
     138             : 
     139       10431 : PLUMED_REGISTER_ACTION(Path,"PATH")
     140             : 
     141           7 : void Path::registerKeywords( Keywords& keys ) {
     142           7 :   PathBase::registerKeywords( keys ); keys.remove("PROPERTY");
     143          14 :   keys.addFlag("NOSPATH",false,"do not calculate the spath position");
     144          14 :   keys.remove("LOWMEM"); keys.use("GPATH");
     145           7 : }
     146             : 
     147           6 : Path::Path(const ActionOptions& ao):
     148             :   Action(ao),
     149           6 :   PathBase(ao)
     150             : {
     151             :   setLowMemOption( true );
     152          12 :   bool nos; parseFlag("NOSPATH",nos);
     153             : 
     154             :   std::string empty;
     155           6 :   if(!nos) {
     156           4 :     if( getLambda()==0 ) error("you must set LAMBDA parameter in order to calculate spath position.  Use LAMBDA/NOSPATH keyword");
     157             :     empty="LABEL=spath";
     158           8 :     addVessel("SPATH",empty,0);
     159             :   }
     160           6 :   readVesselKeywords();
     161           6 :   checkRead();
     162           6 : }
     163             : 
     164             : }
     165             : }