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 "core/ActionWithArguments.h"
      23             : #include "core/ActionWithValue.h"
      24             : #include "core/ActionRegister.h"
      25             : 
      26             : //+PLUMEDOC COLVAR SUM
      27             : /*
      28             : Calculate the sum of the arguments
      29             : 
      30             : This action takes a single vector, a single matrix or a single grid in input. The output is a scalar
      31             : that contains the sum of all the elements in the input vector/matrix/grid. This action is
      32             : very useful if you want to do calculations like the one illustrated in this example:
      33             : 
      34             : ```plumed
      35             : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
      36             : b: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.1}
      37             : s: SUM ARG=b PERIODIC=NO
      38             : PRINT ARG=s FILE=colvar
      39             : ```
      40             : 
      41             : This example calculates and outputs the number of distances that are less than 0.1 nm.
      42             : 
      43             : You can do something similar by summing the elements of a matrix as shown below:
      44             : 
      45             : ```plumed
      46             : c: CONTACT_MATRIX GROUP=1-100 SWITCH={RATIONAL R_0=0.1}
      47             : s: SUM ARG=c PERIODIC=NO
      48             : PRINT ARG=s FILE=colvar
      49             : ```
      50             : 
      51             : If you want to sum all the elements in a grid you can. However, if the input is a function on a grid
      52             : it is more likely that you will want to compute the integral using [INTEGRATE_GRID](INTEGRATE_GRID.md).
      53             : 
      54             : */
      55             : //+ENDPLUMEDOC
      56             : 
      57             : //+PLUMEDOC COLVAR MEAN
      58             : /*
      59             : Calculate the arithmetic mean of the elements in a vector
      60             : 
      61             : This action takes a single vector or a single matrix in input.  The output is a scalar
      62             : that contains the mean of all the elements in the input vector/matrix. This action is
      63             : useful if you want do to calculations like the one illustrated in this example:
      64             : 
      65             : ```plumed
      66             : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
      67             : m: MEAN ARG=d PERIODIC=NO
      68             : PRINT ARG=m FILE=colvar
      69             : ```
      70             : 
      71             : The output from the MEAN action here is the average over the four distances that are evaluated by the
      72             : DISTANCE action.  Notice that you can also do the calculation above using:
      73             : 
      74             : ```plumed
      75             : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
      76             : s: SUM ARG=d PERIODIC=NO
      77             : m: CUSTOM ARG=s FUNC=x/4 PERIODIC=NO
      78             : PRINT ARG=m FILE=colvar
      79             : ```
      80             : 
      81             : */
      82             : //+ENDPLUMEDOC
      83             : 
      84             : namespace PLMD {
      85             : namespace function {
      86             : 
      87             : class Sum :
      88             :   public ActionWithValue,
      89             :   public ActionWithArguments {
      90             : private:
      91             :   bool ismatrix;
      92             : public:
      93             :   static void registerKeywords(Keywords& keys);
      94             :   explicit Sum(const ActionOptions&);
      95             :   unsigned getNumberOfDerivatives() override ;
      96             :   void calculate() override ;
      97             :   void apply() override ;
      98             : };
      99             : 
     100             : PLUMED_REGISTER_ACTION(Sum,"SUM")
     101             : PLUMED_REGISTER_ACTION(Sum,"MEAN")
     102             : 
     103        1495 : void Sum::registerKeywords(Keywords& keys) {
     104        1495 :   Action::registerKeywords(keys);
     105        1495 :   ActionWithValue::registerKeywords(keys);
     106        1495 :   ActionWithArguments::registerKeywords(keys);
     107        2990 :   keys.addInputKeyword("compulsory","ARG","scalar/vector/matrix/grid","the vector/matrix/grid whose elements should be added together (the sum of a scalar is the scalar itself)");
     108        1495 :   keys.add("compulsory","PERIODIC","if the output of your function is periodic then you should specify the periodicity of the function.  If the output is not periodic you must state this using PERIODIC=NO");
     109        5980 :   keys.setValueDescription("scalar","the " + keys.getDisplayName() + " of the elements in the input value");
     110        1495 :   keys.add("hidden","MASKED_INPUT_ALLOWED","turns on that you are allowed to use masked inputs ");
     111        2990 :   keys.reset_style("NUMERICAL_DERIVATIVES","hidden");
     112        1495 : }
     113             : 
     114         871 : Sum::Sum(const ActionOptions&ao):
     115             :   Action(ao),
     116             :   ActionWithValue(ao),
     117             :   ActionWithArguments(ao),
     118         871 :   ismatrix(false) {
     119             : 
     120         871 :   if( getNumberOfArguments()!=1 ) {
     121           0 :     error("should only be one argument in input for this action");
     122             :   }
     123         875 :   ismatrix = getPntrToArgument(0)->getRank()==2 && !getPntrToArgument(0)->hasDerivatives();
     124             : 
     125         871 :   addValue();
     126         871 :   std::vector<std::string> period;
     127        1742 :   parseVector("PERIODIC",period);
     128         871 :   if( period.size()==1 ) {
     129         870 :     if( period[0]!="NO") {
     130           0 :       error("input to PERIODIC keyword does not make sense");
     131             :     }
     132         870 :     setNotPeriodic();
     133           1 :   } else if( period.size()==2 ) {
     134           1 :     setPeriodic( period[0], period[1] );
     135             :   } else {
     136           0 :     error("input to PERIODIC keyword does not make sense");
     137             :   }
     138         871 : }
     139             : 
     140        1056 : unsigned Sum::getNumberOfDerivatives() {
     141        1056 :   return 0;
     142             : }
     143             : 
     144       50056 : void Sum::calculate() {
     145             :   double myval = 0;
     146             :   const Value* myarg = getPntrToArgument(0);
     147       50056 :   if( ismatrix ) {
     148         371 :     unsigned ncols = myarg->getNumberOfColumns();
     149       13936 :     for(unsigned i=0; i<myarg->getShape()[0]; ++i) {
     150      521346 :       for(unsigned j=0; j<myarg->getRowLength(i); ++j) {
     151      507781 :         myval += myarg->get(i*ncols + j, false);
     152             :       }
     153             :     }
     154             :   } else {
     155       49685 :     unsigned nargs = myarg->getNumberOfStoredValues();
     156     3201945 :     for(unsigned i=0; i<nargs; ++i) {
     157     3152260 :       myval += myarg->get(i);
     158             :     }
     159             :   }
     160       50056 :   if( getName()=="MEAN" ) {
     161        3817 :     getPntrToComponent(0)->set( myval / myarg->getNumberOfValues() );
     162             :   } else {
     163       46239 :     getPntrToComponent(0)->set( myval );
     164             :   }
     165       50056 : }
     166             : 
     167       50036 : void Sum::apply() {
     168       50036 :   if( !getPntrToComponent(0)->forcesWereAdded() ) {
     169             :     return ;
     170             :   }
     171             : 
     172             :   Value* myarg = getPntrToArgument(0);
     173       40549 :   double force = getPntrToComponent(0)->getForce(0);
     174       40549 :   if( getName()=="MEAN" ) {
     175        1861 :     force = force / myarg->getNumberOfValues();
     176             :   }
     177       40549 :   if( ismatrix ) {
     178         140 :     unsigned ncols = myarg->getNumberOfColumns();
     179        4068 :     for(unsigned i=0; i<myarg->getShape()[0]; ++i) {
     180      231324 :       for(unsigned j=0; j<myarg->getRowLength(i); ++j) {
     181      227396 :         myarg->addForce(i*ncols + j, force, false);
     182             :       }
     183             :     }
     184             :   } else {
     185       40409 :     unsigned nargs = myarg->getNumberOfStoredValues();
     186     2488594 :     for(unsigned i=0; i<nargs; ++i) {
     187     2448185 :       myarg->addForce( i, force );
     188             :     }
     189             :   }
     190             : }
     191             : 
     192             : }
     193             : 
     194             : 
     195             : }