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       51206 : ActionOptions::ActionOptions(PlumedMain&p,const std::vector<std::string>&l):
      44       51206 :   plumed(p),
      45       51206 :   line(l),
      46       51206 :   keys(emptyKeys)
      47             : {
      48       51206 : }
      49             : 
      50       51204 : ActionOptions::ActionOptions(const ActionOptions&ao,const Keywords&keys):
      51       51204 :   plumed(ao.plumed),
      52       51204 :   line(ao.line),
      53       51204 :   keys(keys)
      54             : {
      55       51204 : }
      56             : 
      57       74188 : void Action::registerKeywords( Keywords& keys ) {
      58       74188 :   plumed_assert( keys.size()==0 );
      59      148376 :   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" );
      60      148376 :   keys.reserve("optional","UPDATE_FROM","Only update this action from this time");
      61      148376 :   keys.reserve("optional","UPDATE_UNTIL","Only update this action until this time");
      62      148376 :   keys.reserve("optional","RESTART","allows per-action setting of restart (YES/NO/AUTO)");
      63       74188 : }
      64             : 
      65       51204 : Action::Action(const ActionOptions&ao):
      66       51204 :   name(ao.line[0]),
      67       51204 :   line(ao.line),
      68       51204 :   update_from(std::numeric_limits<double>::max()),
      69       51204 :   update_until(std::numeric_limits<double>::max()),
      70       51204 :   timestep(0),
      71       51204 :   active(false),
      72       51204 :   restart(ao.plumed.getRestart()),
      73       51204 :   doCheckPoint(ao.plumed.getCPT()),
      74       51204 :   never_activate(name=="CONSTANT"),
      75       51204 :   plumed(ao.plumed),
      76       51204 :   log(plumed.getLog()),
      77       51204 :   comm(plumed.comm),
      78       51204 :   multi_sim_comm(plumed.multi_sim_comm),
      79      102408 :   keywords(ao.keys)
      80             : {
      81             :   // Retrieve the timestep and save it
      82       51204 :   resetStoredTimestep();
      83             : 
      84             :   line.erase(line.begin());
      85      102408 :   if( !keywords.exists("NO_ACTION_LOG") ) {
      86       35784 :     log.printf("Action %s\n",name.c_str());
      87       35784 :     if(ao.fullPath.length()>0) log<<"  from library: "<<ao.fullPath<<"\n";
      88             :   }
      89             : 
      90       51204 :   if(comm.Get_rank()==0) {
      91       37295 :     replica_index=multi_sim_comm.Get_rank();
      92             :   }
      93       51204 :   comm.Bcast(replica_index,0);
      94             : 
      95      153398 :   if ( keywords.exists("LABEL") ) { parse("LABEL",label); }
      96       51204 :   if(label.length()==0) {
      97        3128 :     std::string s; Tools::convert(plumed.getActionSet().size()-plumed.getActionSet().select<ActionForInterface*>().size(),s);
      98        6256 :     label="@"+s;
      99       48082 :   } else if ( label.find(".")!=std::string::npos ) warning("using full stop in an action label should be avaoided as . has a special meaning in PLUMED action labels");
     100       51204 :   if( plumed.getActionSet().selectWithLabel<Action*>(label) ) error("label " + label + " has been already used");
     101      102408 :   if( !keywords.exists("NO_ACTION_LOG") ) log.printf("  with label %s\n",label.c_str());
     102      104647 :   if ( keywords.exists("UPDATE_FROM") ) parse("UPDATE_FROM",update_from);
     103      102408 :   if( !keywords.exists("NO_ACTION_LOG") && update_from!=std::numeric_limits<double>::max()) log.printf("  only update from time %f\n",update_from);
     104      104647 :   if ( keywords.exists("UPDATE_UNTIL") ) parse("UPDATE_UNTIL",update_until);
     105      102408 :   if( !keywords.exists("NO_ACTION_LOG") && update_until!=std::numeric_limits<double>::max()) log.printf("  only update until time %f\n",update_until);
     106      102408 :   if ( keywords.exists("RESTART") ) {
     107        2167 :     std::string srestart="AUTO";
     108        2166 :     parse("RESTART",srestart);
     109        2166 :     if( plumed.parseOnlyMode() ) restart=false;
     110        2166 :     else if(srestart=="YES") restart=true;
     111        2005 :     else if(srestart=="NO")  restart=false;
     112        1983 :     else if(srestart=="AUTO") {
     113             :       // do nothing, this is the default
     114           2 :     } else error("RESTART should be either YES, NO, or AUTO");
     115             :   }
     116       51204 : }
     117             : 
     118      100890 : void Action::resetStoredTimestep() {
     119      100890 :   ActionWithValue* ts = plumed.getActionSet().selectWithLabel<ActionWithValue*>("timestep");
     120      100890 :   if( ts ) timestep = (ts->copyOutput(0))->get();
     121      100890 : }
     122             : 
     123      102406 : Action::~Action() {
     124       51203 :   if(files.size()!=0) {
     125           0 :     std::cerr<<"WARNING: some files open in action "+getLabel()+" where not properly closed. This could lead to data loss!!\n";
     126             :   }
     127      102406 : }
     128             : 
     129          82 : FILE* Action::fopen(const char *path, const char *mode) {
     130             :   bool write(false);
     131         164 :   for(const char*p=mode; *p; p++) if(*p=='w' || *p=='a' || *p=='+') write=true;
     132             :   FILE* fp;
     133          82 :   if(write && comm.Get_rank()!=0) fp=plumed.fopen("/dev/null",mode);
     134          82 :   else      fp=plumed.fopen(path,mode);
     135          81 :   files.insert(fp);
     136          81 :   return fp;
     137             : }
     138             : 
     139          99 : int Action::fclose(FILE*fp) {
     140             :   files.erase(fp);
     141          99 :   return plumed.fclose(fp);
     142             : }
     143             : 
     144       61094 : void Action::fflush() {
     145       61094 :   for(const auto & p : files) {
     146           0 :     std::fflush(p);
     147             :   }
     148       61094 : }
     149             : 
     150           0 : std::string Action::getKeyword(const std::string& key) {
     151             :   // Check keyword has been registered
     152           0 :   plumed_massert(keywords.exists(key), "keyword " + key + " has not been registered");
     153             : 
     154             :   std::string outkey;
     155           0 :   if( Tools::getKey(line,key,outkey ) ) return key + outkey;
     156             : 
     157           0 :   if( keywords.style(key,"compulsory") ) {
     158           0 :     if( keywords.getDefaultValue(key,outkey) ) {
     159           0 :       if( outkey.length()==0 ) error("keyword " + key + " has weird default value");
     160           0 :       return key + "=" +  outkey;
     161             :     } else {
     162           0 :       error("keyword " + key + " is compulsory for this action");
     163             :     }
     164             :   }
     165           0 :   return "";
     166             : }
     167             : 
     168       95780 : void Action::parseFlag(const std::string&key,bool & t) {
     169             :   // Check keyword has been registered
     170       95780 :   plumed_massert(keywords.exists(key), "keyword " + key + " has not been registered");
     171             :   // Check keyword is a flag
     172      191560 :   if(!keywords.style(key,"nohtml")) {
     173      191560 :     plumed_massert( keywords.style(key,"vessel") || keywords.style(key,"flag") || keywords.style(key,"hidden"), "keyword " + key + " is not a flag");
     174             :   }
     175             : 
     176             :   // Read in the flag otherwise get the default value from the keywords object
     177       95780 :   if(!Tools::parseFlag(line,key,t)) {
     178      168008 :     if( keywords.style(key,"nohtml") || keywords.style(key,"vessel") ) {
     179           0 :       t=false;
     180       84004 :     } else if ( !keywords.getLogicalDefault(key,t) ) {
     181           0 :       log.printf("ERROR in action %s with label %s : flag %s has no default",name.c_str(),label.c_str(),key.c_str() );
     182           0 :       plumed_error();
     183             :     }
     184             :   }
     185       95780 : }
     186             : 
     187     1263912 : void Action::addDependency(Action*action) {
     188             :   bool found=false;
     189    40069541 :   for(const auto & d : after ) {
     190    39068588 :     if( action==d ) { found=true; break; }
     191             :   }
     192     1263912 :   if( !found ) after.push_back(action);
     193     1263912 : }
     194             : 
     195      571368 : bool Action::checkForDependency( Action* action ) {
     196     1124742 :   for(const auto & d : after) {
     197      563928 :     if( action==d ) { return true; }
     198      557559 :     if( d->checkForDependency(action) ) { return true; }
     199             :   }
     200             :   return false;
     201             : }
     202             : 
     203     5195048 : void Action::activate() {
     204             : // This is set to true if actions are only need to be computed in setup (during checkRead)
     205     5195048 :   if( never_activate ) return;
     206             : // preparation step is called only the first time an Action is activated.
     207             : // since it could change its dependences (e.g. in an ActionAtomistic which is
     208             : // accessing to a virtual atom), this is done just before dependencies are
     209             : // activated
     210     5080682 :   if(!active) {
     211     2624180 :     this->unlockRequests();
     212     2624180 :     prepare();
     213     2624180 :     this->lockRequests();
     214             :   } else return;
     215     6476010 :   for(const auto & p : after) p->activate();
     216     2624180 :   active=true;
     217             : }
     218             : 
     219        2228 : void Action::setOption(const std::string &s) {
     220             : // This overloads the action and activate some options
     221        2228 :   options.insert(s);
     222        4370 :   for(const auto & p : after) p->setOption(s);
     223        2228 : }
     224             : 
     225           0 : void Action::clearOptions() {
     226             : // This overloads the action and activate some options
     227             :   options.clear();
     228           0 : }
     229             : 
     230             : 
     231      199672 : void Action::clearDependencies() {
     232             :   after.clear();
     233      199672 : }
     234             : 
     235       67297 : void Action::checkRead() {
     236       67297 :   if(!line.empty()) {
     237           1 :     std::string msg="cannot understand the following words from the input line : ";
     238           2 :     for(unsigned i=0; i<line.size(); i++) {
     239           1 :       if(i>0) msg = msg + ", ";
     240           2 :       msg = msg + line[i];
     241             :     }
     242           1 :     error(msg);
     243             :   }
     244       67296 :   setupConstantValues(false);
     245       67296 : }
     246             : 
     247       95247 : void Action::setupConstantValues( const bool& have_atoms ) {
     248       95247 :   if( have_atoms ) {
     249             :     // This ensures that we switch off actions that only depend on constant when passed from the
     250             :     // MD code on the first step
     251       27951 :     ActionAtomistic* at = castToActionAtomistic();
     252       27951 :     ActionWithValue* av = castToActionWithValue();
     253       27951 :     if( at && av ) {
     254       15930 :       never_activate=av->getNumberOfComponents()>0;
     255       15977 :       for(unsigned i=0; i<av->getNumberOfComponents(); ++i) {
     256       15930 :         if( !av->copyOutput(i)->isConstant() ) { never_activate=false; break; }
     257             :       }
     258             :     }
     259             :   }
     260       95247 :   ActionWithArguments* aa = castToActionWithArguments();
     261      114271 :   if( aa && aa->getNumberOfArguments()>0 && getName()!="BIASVALUE" ) never_activate = aa->calculateConstantValues( have_atoms );
     262       95247 : }
     263             : 
     264     5402881 : long long int Action::getStep()const {
     265     5402881 :   return plumed.getStep();
     266             : }
     267             : 
     268     3671868 : double Action::getTime()const {
     269     3671868 :   return timestep*getStep();
     270             : }
     271             : 
     272      218036 : double Action::getTimeStep()const {
     273      218036 :   return timestep;
     274             : }
     275             : 
     276         635 : double Action::getkBT() {
     277         635 :   double temp=-1.0;
     278        1891 :   if( keywords.exists("TEMP") ) parse("TEMP",temp);
     279        1579 :   if(temp>=0.0 && keywords.style("TEMP","optional") ) return getKBoltzmann()*temp;
     280         250 :   ActionForInterface* kb=plumed.getActionSet().selectWithLabel<ActionForInterface*>("kBT");
     281         250 :   double kbt=0; if(kb) kbt=(kb->copyOutput(0))->get();
     282         366 :   if( temp>=0 && keywords.style("TEMP","compulsory") ) {
     283          58 :     double kB=getKBoltzmann();
     284          58 :     if( kbt>0 && std::abs(kbt-kB*temp)>1e-4) {
     285           0 :       std::string strt1, strt2; Tools::convert( temp, strt1 ); Tools::convert( kbt/kB, strt2 );
     286           0 :       warning("using TEMP=" + strt1 + " while MD engine uses " + strt2 + "\n");
     287             :     }
     288          58 :     kbt = kB*temp;
     289          58 :     plumed_massert(kbt>0,"your MD engine does not pass the temperature to plumed, you must specify it using TEMP");
     290             :     return kbt;
     291             :   }
     292             :   return kbt;
     293             : }
     294             : 
     295           0 : void Action::exit(int c) {
     296           0 :   plumed.exit(c);
     297           0 : }
     298             : 
     299           0 : void Action::calculateNumericalDerivatives( ActionWithValue* a ) {
     300           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");
     301             : }
     302             : 
     303     1004472 : void Action::prepare() {
     304     1004472 :   return;
     305             : }
     306             : 
     307          28 : [[noreturn]] void Action::error( const std::string & msg ) const {
     308          56 :   if( !keywords.exists("NO_ACTION_LOG") ) log.printf("ERROR in input to action %s with label %s : %s \n \n", name.c_str(), label.c_str(), msg.c_str() );
     309          84 :   plumed_merror("ERROR in input to action " + name + " with label " + label + " : " + msg );
     310             : }
     311             : 
     312        4287 : void Action::warning( const std::string & msg ) {
     313        4287 :   log.printf("WARNING for action %s with label %s : %s \n", name.c_str(), label.c_str(), msg.c_str() );
     314        4287 : }
     315             : 
     316           0 : void Action::calculateFromPDB( const PDB& pdb ) {
     317           0 :   activate();
     318           0 :   for(const auto & p : after) {
     319           0 :     ActionWithValue*av=castToActionWithValue();
     320           0 :     if(av) { av->clearInputForces(); av->clearDerivatives(); }
     321           0 :     p->readAtomsFromPDB( pdb );
     322           0 :     p->calculate();
     323             :   }
     324           0 :   readAtomsFromPDB( pdb );
     325           0 :   calculate();
     326           0 : }
     327             : 
     328       28109 : bool Action::getExchangeStep()const {
     329       28109 :   return plumed.getExchangeStep();
     330             : }
     331             : 
     332          40 : std::string Action::cite(const std::string&s) {
     333          40 :   return plumed.cite(s);
     334             : }
     335             : 
     336             : /// Check if action should be updated.
     337     2524620 : bool Action::checkUpdate()const {
     338     2524620 :   double t=getTime();
     339     2524620 :   if(t<update_until && (update_from==std::numeric_limits<double>::max() || t>=update_from)) return true;
     340         510 :   else return false;
     341             : }
     342             : 
     343        1095 : bool Action::getCPT() const {
     344        1095 :   return plumed.getCPT();
     345             : }
     346             : 
     347    16246561 : const Units& Action::getUnits() const {
     348    16246561 :   return plumed.getUnits();
     349             : }
     350             : 
     351        2142 : bool Action::usingNaturalUnits() const {
     352        2142 :   return plumed.usingNaturalUnits();
     353             : }
     354             : 
     355         678 : double Action::getKBoltzmann() const {
     356         678 :   if( usingNaturalUnits() ) return 1.0;
     357         678 :   else return kBoltzmann/getUnits().getEnergy();
     358             : }
     359             : 
     360          88 : std::string Action::writeInGraph() const {
     361          44 :   std::string nam=getName();
     362          44 :   std::size_t u=nam.find_last_of("_"); std::string sub=nam.substr(u+1);
     363         118 :   if( sub=="SCALAR" || sub=="VECTOR" || sub=="GRID" ) return nam.substr(0,u);
     364             :   return nam;
     365             : }
     366             : 
     367             : }
     368             :