Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 : Copyright (c) 2018-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 "core/ActionShortcut.h" 23 : #include "core/ActionRegister.h" 24 : 25 : namespace PLMD { 26 : namespace wham { 27 : 28 : //+PLUMEDOC REWEIGHTING WHAM_WEIGHTS 29 : /* 30 : Calculate and output weights for configurations using the weighted histogram analysis method. 31 : 32 : This shortcut action allows you to calculate and output weights computed using the weighted histogram 33 : analysis technique. For more detail on how this technique works see \ref REWEIGHT_WHAM 34 : 35 : \par Examples 36 : 37 : The following input can be used to analyze the output from a series of umbrella sampling calculations. 38 : The trajectory from each of the simulations run with the different biases should be concatenated into a 39 : single trajectory before running the following analysis script on the concatenated trajectory using PLUMED 40 : driver. The umbrella sampling simulations that will be analyzed using the script below applied a harmonic 41 : restraint that restrained the torsional angle involving atoms 5, 7, 9 and 15 to particular values. The script 42 : below calculates the reweighting weights for each of the trajectories and then applies the binless WHAM algorithm 43 : to determine a weight for each configuration in the concatenated trajectory. 44 : 45 : \plumedfile 46 : #SETTINGS NREPLICAS=4 47 : phi: TORSION ATOMS=5,7,9,15 48 : rp: RESTRAINT ARG=phi KAPPA=50.0 ... 49 : AT=@replicas:{ 50 : -3.00000000000000000000 51 : -1.45161290322580645168 52 : .09677419354838709664 53 : 1.64516129032258064496 54 : } 55 : ... 56 : 57 : WHAM_WEIGHTS BIAS=rp.bias TEMP=300 FILE=wham-weights 58 : \endplumedfile 59 : 60 : The script above must be run with multiple replicas using the following command: 61 : 62 : \verbatim 63 : mpirun -np 4 plumed driver --mf_xtc alltraj.xtc --multi 4 64 : \endverbatim 65 : 66 : */ 67 : //+ENDPLUMEDOC 68 : 69 : class WhamWeights : public ActionShortcut { 70 : public: 71 : static void registerKeywords( Keywords& keys ); 72 : explicit WhamWeights( const ActionOptions& ); 73 : }; 74 : 75 : PLUMED_REGISTER_ACTION(WhamWeights,"WHAM_WEIGHTS") 76 : 77 8 : void WhamWeights::registerKeywords( Keywords& keys ) { 78 8 : ActionShortcut::registerKeywords( keys ); keys.remove("LABEL"); 79 16 : keys.add("compulsory","BIAS","*.bias","the value of the biases to use when performing WHAM"); 80 16 : keys.add("optional","TEMP","the temperature at which the simulation was run"); 81 16 : keys.add("compulsory","STRIDE","1","the frequency with which the bias should be stored to perform WHAM"); 82 16 : keys.add("compulsory","FILE","the file on which to output the WHAM weights"); 83 16 : keys.add("optional","FMT","the format to use for the real numbers in the output file"); 84 8 : keys.setValueDescription("the weights that were calculated using WHAM"); 85 16 : keys.needsAction("GATHER_REPLICAS"); keys.needsAction("CONCATENATE"); 86 24 : keys.needsAction("COLLECT"); keys.needsAction("WHAM"); keys.needsAction("DUMPVECTOR"); 87 8 : } 88 : 89 6 : WhamWeights::WhamWeights( const ActionOptions& ao ) : 90 : Action(ao), 91 6 : ActionShortcut(ao) 92 : { 93 : // Input for collection of weights for WHAM 94 12 : std::string bias; parse("BIAS",bias); 95 6 : std::string stride; parse("STRIDE",stride); 96 : // Input for GATHER_REPLICAS 97 12 : readInputLine( getShortcutLabel() + "_gather: GATHER_REPLICAS ARG=" + bias ); 98 : // Put all the replicas in a single vector 99 12 : readInputLine( getShortcutLabel() + "_gatherv: CONCATENATE ARG=" + getShortcutLabel() + "_gather.*"); 100 : // Input for COLLECT_FRAMES 101 12 : readInputLine( getShortcutLabel() + "_collect: COLLECT TYPE=vector ARG=" + getShortcutLabel() + "_gatherv STRIDE=" + stride); 102 : // Input for WHAM 103 18 : std::string temp, tempstr=""; parse("TEMP",temp); if( temp.length()>0 ) tempstr="TEMP=" + temp; 104 12 : readInputLine( getShortcutLabel() + ": WHAM ARG=" + getShortcutLabel() + "_collect " + tempstr ); 105 : // Input for PRINT (will just output at end of calc 106 12 : std::string filename, fmt; parse("FILE",filename); parse("FMT",fmt); 107 12 : readInputLine( "DUMPVECTOR STRIDE=0 ARG=" + getShortcutLabel() + " FILE=" + filename + " FMT=" + fmt ); 108 6 : } 109 : 110 : } 111 : }