Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2018 The VES code team
3 : (see the PEOPLE-VES file at the root of this folder for a list of names)
4 :
5 : See http://www.ves-code.org for more information.
6 :
7 : This file is part of VES code module.
8 :
9 : The VES code module 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 : The VES code module 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 the VES code module. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 :
23 : #include "BasisFunctions.h"
24 :
25 : #include "core/ActionRegister.h"
26 :
27 :
28 : namespace PLMD {
29 : namespace ves {
30 :
31 : //+PLUMEDOC VES_BASISF BF_CHEBYSHEV
32 : /*
33 : Chebyshev polynomial basis functions.
34 :
35 : Use as basis functions [Chebyshev polynomials](https://en.wikipedia.org/wiki/Chebyshev_polynomials)
36 : of the first kind \f$T_{n}(x)\f$ defined on a bounded interval.
37 : You need to provide the interval \f$[a,b]\f$
38 : on which the basis functions are to be used, and the order of the
39 : expansion \f$N\f$ (i.e. the highest order polynomial used).
40 : The total number of basis functions is \f$N+1\f$ as the constant \f$T_{0}(x)=1\f$
41 : is also included.
42 : These basis functions should not be used for periodic CVs.
43 :
44 : Intrinsically the Chebyshev polynomials are defined on the interval \f$[-1,1]\f$.
45 : A variable \f$t\f$ in the interval \f$[a,b]\f$ is transformed to a variable \f$x\f$
46 : in the intrinsic interval \f$[-1,1]\f$ by using the transform function
47 : \f[
48 : x(t) = \frac{t-(a+b)/2}
49 : {(b-a)/2}
50 : \f]
51 :
52 : The Chebyshev polynomials are given by the recurrence relation
53 : \f{align}{
54 : T_{0}(x) &= 1 \\
55 : T_{1}(x) &= x \\
56 : T_{n+1}(x) &= 2 \, x \, T_{n}(x) - T_{n-1}(x)
57 : \f}
58 :
59 : The first 6 polynomials are shown below
60 : \image html ves_basisf-chebyshev.png
61 :
62 : The Chebyshev polynomial are orthogonal over the interval \f$[-1,1]\f$
63 : with respect to the weight \f$\frac{1}{\sqrt{1-x^2}}\f$
64 : \f[
65 : \int_{-1}^{1} dx \, T_{n}(x)\, T_{m}(x) \, \frac{1}{\sqrt{1-x^2}} =
66 : \begin{cases}
67 : 0 & n \neq m \\
68 : \pi & n = m = 0 \\
69 : \pi/2 & n = m \neq 0
70 : \end{cases}
71 : \f]
72 :
73 : For further mathematical properties of the Chebyshev polynomials see for example
74 : the [Wikipedia page](https://en.wikipedia.org/wiki/Chebyshev_polynomials).
75 :
76 : \par Examples
77 :
78 : Here we employ a Chebyshev expansion of order 20 over the interval 0.0 to 10.0.
79 : This results in a total number of 21 basis functions.
80 : The label used to identify the basis function action can then be
81 : referenced later on in the input file.
82 : \plumedfile
83 : bfC: BF_CHEBYSHEV MINIMUM=0.0 MAXIMUM=10.0 ORDER=20
84 : \endplumedfile
85 :
86 : */
87 : //+ENDPLUMEDOC
88 :
89 :
90 4 : class BF_Chebyshev : public BasisFunctions {
91 : virtual void setupUniformIntegrals();
92 : public:
93 : static void registerKeywords(Keywords&);
94 : explicit BF_Chebyshev(const ActionOptions&);
95 : void getAllValues(const double, double&, bool&, std::vector<double>&, std::vector<double>&) const;
96 : };
97 :
98 :
99 6456 : PLUMED_REGISTER_ACTION(BF_Chebyshev,"BF_CHEBYSHEV")
100 :
101 :
102 5 : void BF_Chebyshev::registerKeywords(Keywords& keys) {
103 5 : BasisFunctions::registerKeywords(keys);
104 5 : }
105 :
106 4 : BF_Chebyshev::BF_Chebyshev(const ActionOptions&ao):
107 4 : PLUMED_VES_BASISFUNCTIONS_INIT(ao)
108 : {
109 4 : setNumberOfBasisFunctions(getOrder()+1);
110 12 : setIntrinsicInterval("-1.0","+1.0");
111 : setNonPeriodic();
112 : setIntervalBounded();
113 8 : setType("chebyshev-1st-kind");
114 8 : setDescription("Chebyshev polynomials of the first kind");
115 8 : setLabelPrefix("T");
116 4 : setupBF();
117 4 : checkRead();
118 4 : }
119 :
120 :
121 8484 : void BF_Chebyshev::getAllValues(const double arg, double& argT, bool& inside_range, std::vector<double>& values, std::vector<double>& derivs) const {
122 : // plumed_assert(values.size()==numberOfBasisFunctions());
123 : // plumed_assert(derivs.size()==numberOfBasisFunctions());
124 8484 : inside_range=true;
125 8484 : argT=translateArgument(arg, inside_range);
126 8484 : std::vector<double> derivsT(derivs.size());
127 : //
128 8484 : values[0]=1.0;
129 8484 : derivsT[0]=0.0;
130 8484 : derivs[0]=0.0;
131 8484 : values[1]=argT;
132 8484 : derivsT[1]=1.0;
133 8484 : derivs[1]=intervalDerivf();
134 285796 : for(unsigned int i=1; i < getOrder(); i++) {
135 554624 : values[i+1] = 2.0*argT*values[i]-values[i-1];
136 554624 : derivsT[i+1] = 2.0*values[i]+2.0*argT*derivsT[i]-derivsT[i-1];
137 277312 : derivs[i+1] = intervalDerivf()*derivsT[i+1];
138 : }
139 12430 : if(!inside_range) {for(unsigned int i=0; i<derivs.size(); i++) {derivs[i]=0.0;}}
140 8484 : }
141 :
142 :
143 3 : void BF_Chebyshev::setupUniformIntegrals() {
144 109 : for(unsigned int i=0; i<numberOfBasisFunctions(); i++) {
145 53 : double io = i;
146 : double value = 0.0;
147 53 : if(i % 2 == 0) {
148 28 : value = -2.0/( pow(io,2.0)-1.0)*0.5;
149 : }
150 : setUniformIntegral(i,value);
151 : }
152 3 : }
153 :
154 :
155 : }
156 4839 : }
|