Line data    Source code

       1             : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       2             :    Copyright (c) 2011-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 "Action.h"
      23             : #include "ActionAtomistic.h"
      24             : #include "ActionWithValue.h"
      25             : #include "ActionWithArguments.h"
      26             : #include "ActionWithVirtualAtom.h"
      27             : #include "ActionForInterface.h"
      28             : #include "DomainDecomposition.h"
      29             : #include "PbcAction.h"
      30             : #include "ActionToPutData.h"
      31             : #include "ActionToGetData.h"
      32             : #include "PlumedMain.h"
      33             : #include "tools/Log.h"
      34             : #include "tools/Exception.h"
      35             : #include "tools/Communicator.h"
      36             : #include "ActionSet.h"
      37             : #include <iostream>
      38             : 
      39             : namespace PLMD {
      40             : 
      41             : Keywords ActionOptions::emptyKeys;
      42             : 
      43       51099 : ActionOptions::ActionOptions(PlumedMain&p,const std::vector<std::string>&l):
      44       51099 :   plumed(p),
      45       51099 :   line(l),
      46       51099 :   keys(emptyKeys) {
      47       51099 : }
      48             : 
      49       51097 : ActionOptions::ActionOptions(const ActionOptions&ao,const Keywords&keys):
      50       51097 :   plumed(ao.plumed),
      51       51097 :   line(ao.line),
      52       51097 :   keys(keys) {
      53       51097 : }
      54             : 
      55       74737 : void Action::registerKeywords( Keywords& keys ) {
      56       74737 :   plumed_assert( keys.size()==0 );
      57      149474 :   keys.add( "hidden", "LABEL", "a label for the action so that its output can be referenced in the input to other actions.  Actions with scalar output are referenced using their label only.  Actions with vector output must have a separate label for every component.  Individual components are then referred to using label.component" );
      58      149474 :   keys.reserve("optional","UPDATE_FROM","Only update this action from this time");
      59      149474 :   keys.reserve("optional","UPDATE_UNTIL","Only update this action until this time");
      60      149474 :   keys.reserve("optional","RESTART","allows per-action setting of restart (YES/NO/AUTO)");
      61       74737 : }
      62             : 
      63       51097 : Action::Action(const ActionOptions&ao):
      64       51097 :   name(ao.line[0]),
      65       51097 :   line(ao.line),
      66       51097 :   update_from(std::numeric_limits<double>::max()),
      67       51097 :   update_until(std::numeric_limits<double>::max()),
      68       51097 :   timestep(0),
      69       51097 :   active(false),
      70       51097 :   restart(ao.plumed.getRestart()),
      71       51097 :   doCheckPoint(ao.plumed.getCPT()),
      72       51097 :   never_activate(name=="CONSTANT"),
      73       51097 :   plumed(ao.plumed),
      74       51097 :   log(plumed.getLog()),
      75       51097 :   comm(plumed.comm),
      76       51097 :   multi_sim_comm(plumed.multi_sim_comm),
      77      102194 :   keywords(ao.keys) {
      78             :   // Retrieve the timestep and save it
      79       51097 :   resetStoredTimestep();
      80             : 
      81             :   line.erase(line.begin());
      82      102194 :   if( !keywords.exists("NO_ACTION_LOG") ) {
      83       35417 :     log.printf("Action %s\n",name.c_str());
      84       35417 :     if(ao.fullPath.length()>0) {
      85          41 :       log<<"  from library: "<<ao.fullPath<<"\n";
      86             :     }
      87             :   }
      88             : 
      89       51097 :   if(comm.Get_rank()==0) {
      90       37148 :     replica_index=multi_sim_comm.Get_rank();
      91             :   }
      92       51097 :   comm.Bcast(replica_index,0);
      93             : 
      94      102194 :   if ( keywords.exists("LABEL") ) {
      95      101728 :     parse("LABEL",label);
      96             :   }
      97       51097 :   if(label.length()==0) {
      98             :     std::string s;
      99        3244 :     Tools::convert(plumed.getActionSet().size()-plumed.getActionSet().select<ActionForInterface*>().size(),s);
     100        6488 :     label="@"+s;
     101       47853 :   } else if ( label.find(".")!=std::string::npos ) {
     102          12 :     warning("using full stop in an action label should be avaoided as . has a special meaning in PLUMED action labels");
     103             :   }
     104       51097 :   if( plumed.getActionSet().selectWithLabel<Action*>(label) ) {
     105           0 :     error("label " + label + " has been already used");
     106             :   }
     107      102194 :   if( !keywords.exists("NO_ACTION_LOG") ) {
     108       35417 :     log.printf("  with label %s\n",label.c_str());
     109             :   }
     110      102194 :   if ( keywords.exists("UPDATE_FROM") ) {
     111        4592 :     parse("UPDATE_FROM",update_from);
     112             :   }
     113      102194 :   if( !keywords.exists("NO_ACTION_LOG") && update_from!=std::numeric_limits<double>::max()) {
     114           5 :     log.printf("  only update from time %f\n",update_from);
     115             :   }
     116      102194 :   if ( keywords.exists("UPDATE_UNTIL") ) {
     117        4592 :     parse("UPDATE_UNTIL",update_until);
     118             :   }
     119      102194 :   if( !keywords.exists("NO_ACTION_LOG") && update_until!=std::numeric_limits<double>::max()) {
     120           5 :     log.printf("  only update until time %f\n",update_until);
     121             :   }
     122      102194 :   if ( keywords.exists("RESTART") ) {
     123        2261 :     std::string srestart="AUTO";
     124        2260 :     parse("RESTART",srestart);
     125        2260 :     if( plumed.parseOnlyMode() ) {
     126           0 :       restart=false;
     127        2260 :     } else if(srestart=="YES") {
     128         172 :       restart=true;
     129        2088 :     } else if(srestart=="NO") {
     130          22 :       restart=false;
     131        2066 :     } else if(srestart=="AUTO") {
     132             :       // do nothing, this is the default
     133             :     } else {
     134           2 :       error("RESTART should be either YES, NO, or AUTO");
     135             :     }
     136             :   }
     137       51097 : }
     138             : 
     139      100670 : void Action::resetStoredTimestep() {
     140      100670 :   ActionWithValue* ts = plumed.getActionSet().selectWithLabel<ActionWithValue*>("timestep");
     141      100670 :   if( ts ) {
     142       94367 :     timestep = (ts->copyOutput(0))->get();
     143             :   }
     144      100670 : }
     145             : 
     146      102192 : Action::~Action() {
     147       51096 :   if(files.size()!=0) {
     148           0 :     std::cerr<<"WARNING: some files open in action "+getLabel()+" where not properly closed. This could lead to data loss!!\n";
     149             :   }
     150      102192 : }
     151             : 
     152          82 : FILE* Action::fopen(const char *path, const char *mode) {
     153             :   bool write(false);
     154         164 :   for(const char*p=mode; *p; p++)
     155          82 :     if(*p=='w' || *p=='a' || *p=='+') {
     156             :       write=true;
     157             :     }
     158             :   FILE* fp;
     159          82 :   if(write && comm.Get_rank()!=0) {
     160           0 :     fp=plumed.fopen("/dev/null",mode);
     161             :   } else {
     162          82 :     fp=plumed.fopen(path,mode);
     163             :   }
     164          81 :   files.insert(fp);
     165          81 :   return fp;
     166             : }
     167             : 
     168          99 : int Action::fclose(FILE*fp) {
     169             :   files.erase(fp);
     170          99 :   return plumed.fclose(fp);
     171             : }
     172             : 
     173       60971 : void Action::fflush() {
     174       60971 :   for(const auto & p : files) {
     175           0 :     std::fflush(p);
     176             :   }
     177       60971 : }
     178             : 
     179           0 : std::string Action::getKeyword(const std::string& key) {
     180             :   // Check keyword has been registered
     181           0 :   plumed_massert(keywords.exists(key), "keyword " + key + " has not been registered");
     182             : 
     183             :   std::string outkey;
     184           0 :   if( Tools::getKey(line,key,outkey ) ) {
     185           0 :     return key + outkey;
     186             :   }
     187             : 
     188           0 :   if( keywords.style(key,"compulsory") ) {
     189           0 :     if( keywords.getDefaultValue(key,outkey) ) {
     190           0 :       if( outkey.length()==0 ) {
     191           0 :         error("keyword " + key + " has weird default value");
     192             :       }
     193           0 :       return key + "=" +  outkey;
     194             :     } else {
     195           0 :       error("keyword " + key + " is compulsory for this action");
     196             :     }
     197             :   }
     198           0 :   return "";
     199             : }
     200             : 
     201       95772 : void Action::parseFlag(const std::string&key,bool & t) {
     202             :   // Check keyword has been registered
     203       95772 :   plumed_massert(keywords.exists(key), "keyword " + key + " has not been registered");
     204             :   // Check keyword is a flag
     205      191544 :   if(!keywords.style(key,"nohtml")) {
     206      191544 :     plumed_massert( keywords.style(key,"vessel") || keywords.style(key,"flag") || keywords.style(key,"hidden"), "keyword " + key + " is not a flag");
     207             :   }
     208             : 
     209             :   // Read in the flag otherwise get the default value from the keywords object
     210       95772 :   if(!Tools::parseFlag(line,key,t)) {
     211      167406 :     if( keywords.style(key,"nohtml") || keywords.style(key,"vessel") ) {
     212           0 :       t=false;
     213       83703 :     } else if ( !keywords.getLogicalDefault(key,t) ) {
     214           0 :       log.printf("ERROR in action %s with label %s : flag %s has no default",name.c_str(),label.c_str(),key.c_str() );
     215           0 :       plumed_error();
     216             :     }
     217             :   }
     218       95772 : }
     219             : 
     220     1263000 : void Action::addDependency(Action*action) {
     221             :   bool found=false;
     222    40065387 :   for(const auto & d : after ) {
     223    39065479 :     if( action==d ) {
     224             :       found=true;
     225             :       break;
     226             :     }
     227             :   }
     228     1263000 :   if( !found ) {
     229      999908 :     after.push_back(action);
     230             :   }
     231     1263000 : }
     232             : 
     233      571383 : bool Action::checkForDependency( Action* action ) {
     234     1124757 :   for(const auto & d : after) {
     235      563943 :     if( action==d ) {
     236             :       return true;
     237             :     }
     238      557562 :     if( d->checkForDependency(action) ) {
     239             :       return true;
     240             :     }
     241             :   }
     242             :   return false;
     243             : }
     244             : 
     245     5247031 : void Action::activate() {
     246             : // This is set to true if actions are only need to be computed in setup (during checkRead)
     247     5247031 :   if( never_activate ) {
     248             :     return;
     249             :   }
     250             : // preparation step is called only the first time an Action is activated.
     251             : // since it could change its dependences (e.g. in an ActionAtomistic which is
     252             : // accessing to a virtual atom), this is done just before dependencies are
     253             : // activated
     254     5132663 :   if(!active) {
     255     2649967 :     this->unlockRequests();
     256     2649967 :     prepare();
     257     2649967 :     this->lockRequests();
     258             :   } else {
     259             :     return;
     260             :   }
     261     6545200 :   for(const auto & p : after) {
     262     3895233 :     p->activate();
     263             :   }
     264     2649967 :   active=true;
     265             : }
     266             : 
     267        2228 : void Action::setOption(const std::string &s) {
     268             : // This overloads the action and activate some options
     269        2228 :   options.insert(s);
     270        4370 :   for(const auto & p : after) {
     271        2142 :     p->setOption(s);
     272             :   }
     273        2228 : }
     274             : 
     275           0 : void Action::clearOptions() {
     276             : // This overloads the action and activate some options
     277             :   options.clear();
     278           0 : }
     279             : 
     280             : 
     281      199566 : void Action::clearDependencies() {
     282             :   after.clear();
     283      199566 : }
     284             : 
     285       66960 : void Action::checkRead() {
     286       66960 :   if(!line.empty()) {
     287           1 :     std::string msg="cannot understand the following words from the input line : ";
     288           2 :     for(unsigned i=0; i<line.size(); i++) {
     289           1 :       if(i>0) {
     290           0 :         msg = msg + ", ";
     291             :       }
     292           2 :       msg = msg + line[i];
     293             :     }
     294           1 :     error(msg);
     295             :   }
     296       66959 :   setupConstantValues(false);
     297       66959 : }
     298             : 
     299       95294 : void Action::setupConstantValues( const bool& have_atoms ) {
     300       95294 :   if( have_atoms ) {
     301             :     // This ensures that we switch off actions that only depend on constant when passed from the
     302             :     // MD code on the first step
     303       28335 :     ActionAtomistic* at = castToActionAtomistic();
     304       28335 :     ActionWithValue* av = castToActionWithValue();
     305       28335 :     if( at && av ) {
     306       16003 :       never_activate=av->getNumberOfComponents()>0;
     307       16050 :       for(unsigned i=0; i<av->getNumberOfComponents(); ++i) {
     308       16003 :         if( !av->copyOutput(i)->isConstant() ) {
     309       15956 :           never_activate=false;
     310       15956 :           break;
     311             :         }
     312             :       }
     313             :     }
     314             :   }
     315       95294 :   ActionWithArguments* aa = castToActionWithArguments();
     316      114527 :   if( aa && aa->getNumberOfArguments()>0 && getName()!="BIASVALUE" ) {
     317       18592 :     never_activate = aa->calculateConstantValues( have_atoms );
     318             :   }
     319       95294 : }
     320             : 
     321     5469625 : long long int Action::getStep()const {
     322     5469625 :   return plumed.getStep();
     323             : }
     324             : 
     325     3721740 : double Action::getTime()const {
     326     3721740 :   return timestep*getStep();
     327             : }
     328             : 
     329      219887 : double Action::getTimeStep()const {
     330      219887 :   return timestep;
     331             : }
     332             : 
     333         682 : double Action::getkBT() {
     334         682 :   double temp=-1.0;
     335        1364 :   if( keywords.exists("TEMP") ) {
     336        1336 :     parse("TEMP",temp);
     337             :   }
     338        1698 :   if(temp>=0.0 && keywords.style("TEMP","optional") ) {
     339         421 :     return getKBoltzmann()*temp;
     340             :   }
     341         261 :   ActionForInterface* kb=plumed.getActionSet().selectWithLabel<ActionForInterface*>("kBT");
     342             :   double kbt=0;
     343         261 :   if(kb) {
     344           3 :     kbt=(kb->copyOutput(0))->get();
     345             :   }
     346         377 :   if( temp>=0 && keywords.style("TEMP","compulsory") ) {
     347          58 :     double kB=getKBoltzmann();
     348          58 :     if( kbt>0 && std::abs(kbt-kB*temp)>1e-4) {
     349             :       std::string strt1, strt2;
     350           0 :       Tools::convert( temp, strt1 );
     351           0 :       Tools::convert( kbt/kB, strt2 );
     352           0 :       warning("using TEMP=" + strt1 + " while MD engine uses " + strt2 + "\n");
     353             :     }
     354          58 :     kbt = kB*temp;
     355          58 :     plumed_massert(kbt>0,"your MD engine does not pass the temperature to plumed, you must specify it using TEMP");
     356             :     return kbt;
     357             :   }
     358             :   return kbt;
     359             : }
     360             : 
     361           0 : void Action::exit(int c) {
     362           0 :   plumed.exit(c);
     363           0 : }
     364             : 
     365           0 : void Action::calculateNumericalDerivatives( ActionWithValue* a ) {
     366           0 :   plumed_merror("if you get here it means that you are trying to use numerical derivatives for a class that does not implement them");
     367             : }
     368             : 
     369     1025834 : void Action::prepare() {
     370     1025834 :   return;
     371             : }
     372             : 
     373          28 : [[noreturn]] void Action::error( const std::string & msg ) const {
     374          28 :   log.printf("ERROR in input to action %s with label %s : %s \n \n", name.c_str(), label.c_str(), msg.c_str() );
     375          84 :   plumed_merror("ERROR in input to action " + name + " with label " + label + " : " + msg );
     376             : }
     377             : 
     378        4303 : void Action::warning( const std::string & msg ) {
     379        4303 :   log.printf("WARNING for action %s with label %s : %s \n", name.c_str(), label.c_str(), msg.c_str() );
     380        4303 : }
     381             : 
     382           0 : void Action::calculateFromPDB( const PDB& pdb ) {
     383           0 :   activate();
     384           0 :   for(const auto & p : after) {
     385           0 :     ActionWithValue*av=castToActionWithValue();
     386           0 :     if(av) {
     387           0 :       av->clearInputForces();
     388           0 :       av->clearDerivatives();
     389             :     }
     390           0 :     p->readAtomsFromPDB( pdb );
     391           0 :     p->calculate();
     392             :   }
     393           0 :   readAtomsFromPDB( pdb );
     394           0 :   calculate();
     395           0 : }
     396             : 
     397       30198 : bool Action::getExchangeStep()const {
     398       30198 :   return plumed.getExchangeStep();
     399             : }
     400             : 
     401          41 : std::string Action::cite(const std::string&s) {
     402          41 :   return plumed.cite(s);
     403             : }
     404             : 
     405             : /// Check if action should be updated.
     406     2550407 : bool Action::checkUpdate()const {
     407     2550407 :   double t=getTime();
     408     2550407 :   if(t<update_until && (update_from==std::numeric_limits<double>::max() || t>=update_from)) {
     409             :     return true;
     410             :   } else {
     411         510 :     return false;
     412             :   }
     413             : }
     414             : 
     415        1340 : bool Action::getCPT() const {
     416        1340 :   return plumed.getCPT();
     417             : }
     418             : 
     419    15615296 : const Units& Action::getUnits() const {
     420    15615296 :   return plumed.getUnits();
     421             : }
     422             : 
     423        2197 : bool Action::usingNaturalUnits() const {
     424        2197 :   return plumed.usingNaturalUnits();
     425             : }
     426             : 
     427         714 : double Action::getKBoltzmann() const {
     428         714 :   if( usingNaturalUnits() ) {
     429             :     return 1.0;
     430             :   } else {
     431         714 :     return kBoltzmann/getUnits().getEnergy();
     432             :   }
     433             : }
     434             : 
     435          44 : std::string Action::writeInGraph() const {
     436          44 :   std::string nam=getName();
     437          44 :   std::size_t u=nam.find_last_of("_");
     438          44 :   std::string sub=nam.substr(u+1);
     439         118 :   if( sub=="SCALAR" || sub=="VECTOR" || sub=="GRID" ) {
     440           8 :     return nam.substr(0,u);
     441             :   }
     442             :   return nam;
     443             : }
     444             : 
     445             : }
     446             :