1 /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2    Copyright (c) 2011-2021 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 #ifndef __PLUMED_tools_SwitchingFunction_h
23 #define __PLUMED_tools_SwitchingFunction_h
24 
25 #include <string>
26 #include <vector>
27 #include "lepton/Lepton.h"
28 
29 namespace PLMD {
30 
31 class Log;
32 class Keywords;
33 
34 /// \ingroup TOOLBOX
35 /// Small class to compure switching functions.
36 /// Switching functions are created using set() and
37 /// then can be used with function calculate() or calculateSqr().
38 /// Since this is typically computed on a distance vector,
39 /// the second all (calculateSqr()) allows to skip the calculation
40 /// of a square root in some case, thus potentially increasing
41 /// performances.
42 class SwitchingFunction {
43 /// This is to check that switching function has been initialized
44   bool init=false;
45 /// Type of function
46   enum {rational,exponential,gaussian,smap,cubic,tanh,cosinus,matheval,leptontype,nativeq} type=rational;
47 /// Inverse of scaling length.
48 /// We store the inverse to avoid a division
49   double invr0=0.0;
50 /// Minimum distance (before this, function is one)
51   double d0=0.0;
52 /// Maximum distance (after this, function is zero)
53   double dmax=0.0;
54 /// Exponents for rational function
55   int nn=6;
56   int mm=0;
57 /// Parameters for smap function
58   int a=0;
59   int b=0;
60   double c=0.0;
61   double d=0.0;
62 // nativeq
63   double lambda=0.0;
64   double beta=0.0;
65   double ref=0.0;
66 /// Square of invr0, useful in calculateSqr()
67   double invr0_2=0.0;
68 /// Square of dmax, useful in calculateSqr()
69   double dmax_2=0.0;
70 /// Parameters for stretching the function to zero at d_max
71   double stretch=1.0;
72   double shift=0.0;
73 /// Low-level tool to compute rational functions.
74 /// It is separated since it is called both by calculate() and calculateSqr()
75   double do_rational(double rdist,double&dfunc,int nn,int mm)const;
76 /// Function for lepton;
77   std::string lepton_func;
78 /// Lepton expression.
79 /// \warning Since lepton::CompiledExpression is mutable, a vector is necessary for multithreading!
80   std::vector<lepton::CompiledExpression> expression;
81 /// Lepton expression for derivative
82 /// \warning Since lepton::CompiledExpression is mutable, a vector is necessary for multithreading!
83   std::vector<lepton::CompiledExpression> expression_deriv;
84   std::vector<double*> lepton_ref;
85   std::vector<double*> lepton_ref_deriv;
86 /// Set to true for fast rational functions (depending on x**2 only)
87   bool fastrational=false;
88 /// Set to true if lepton only uses x2
89   bool leptonx2=false;
90 public:
91   static void registerKeywords( Keywords& keys );
92 /// Set a "rational" switching function.
93 /// Notice that a d_max is set automatically to a value such that
94 /// f(d_max)=0.00001.
95   void set(int nn,int mm,double r_0,double d_0);
96 /// Set an arbitrary switching function.
97 /// Parse the string in definition and possibly returns errors
98 /// in the errormsg string
99   void set(const std::string& definition, std::string& errormsg);
100 /// Returns a string with a description of the switching function
101   std::string description() const ;
102 /// Compute the switching function.
103 /// Returns s(x). df will be set to the value of the derivative
104 /// of the switching function _divided_by_x
105   double calculate(double x,double&df)const;
106 /// Compute the switching function.
107 /// Returns \f$ s(\sqrt{x})\f$ .
108 /// df will be set to the \f$ \frac{1}{\sqrt{x}}\frac{ds}{d\sqrt{x}}= 2 \frac{ds}{dx}\f$
109 /// (same as calculate()).
110 /// The advantage is that in some case the expensive square root can be avoided
111 /// (namely for rational functions, if nn and mm are even and d0 is zero)
112   double calculateSqr(double distance2,double&dfunc)const;
113 /// Returns d0
114   double get_d0() const;
115 /// Returns r0
116   double get_r0() const;
117 /// Return dmax
118   double get_dmax() const;
119 /// Return dmax squared
120   double get_dmax2() const;
121 };
122 
123 }
124 
125 #endif
126 
127