Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2018 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_Grid_h
23 : #define __PLUMED_tools_Grid_h
24 :
25 : #include <vector>
26 : #include <string>
27 : #include <map>
28 : #include <cmath>
29 :
30 : namespace PLMD {
31 :
32 :
33 : // simple function to enable various weighting
34 :
35 2 : class WeightBase {
36 : public:
37 : virtual double projectInnerLoop(double &input, double &v)=0;
38 : virtual double projectOuterLoop(double &v)=0;
39 2 : virtual ~WeightBase() {}
40 : };
41 :
42 4 : class BiasWeight:public WeightBase {
43 : public:
44 : double beta,invbeta;
45 2 : explicit BiasWeight(double v) {beta=v; invbeta=1./beta;}
46 13120 : double projectInnerLoop(double &input, double &v) {return input+exp(beta*v);}
47 164 : double projectOuterLoop(double &v) {return -invbeta*std::log(v);}
48 : };
49 :
50 : class ProbWeight:public WeightBase {
51 : public:
52 : double beta,invbeta;
53 : explicit ProbWeight(double v) {beta=v; invbeta=1./beta;}
54 : double projectInnerLoop(double &input, double &v) {return input+v;}
55 : double projectOuterLoop(double &v) {return -invbeta*std::log(v);}
56 : };
57 :
58 :
59 :
60 :
61 :
62 :
63 : class Value;
64 : class IFile;
65 : class OFile;
66 : class KernelFunctions;
67 : class Communicator;
68 :
69 : /// \ingroup TOOLBOX
70 8 : class Grid
71 : {
72 : public:
73 : // we use a size_t here
74 : // should be 8 bytes on all 64-bit machines
75 : // and more portable than "unsigned long long"
76 : typedef size_t index_t;
77 : // to restore old implementation (unsigned) use the following instead:
78 : // typedef unsigned index_t;
79 : private:
80 : double contour_location;
81 : std::vector<double> grid_;
82 : std::vector< std::vector<double> > der_;
83 : protected:
84 : std::string funcname;
85 : std::vector<std::string> argnames;
86 : std::vector<std::string> str_min_, str_max_;
87 : std::vector<double> min_,max_,dx_;
88 : std::vector<unsigned> nbin_;
89 : std::vector<bool> pbc_;
90 : index_t maxsize_;
91 : unsigned dimension_;
92 : bool dospline_, usederiv_;
93 : std::string fmt_; // format for output
94 : /// get "neighbors" for spline
95 : std::vector<index_t> getSplineNeighbors(const std::vector<unsigned> & indices)const;
96 :
97 :
98 : public:
99 : /// clear grid
100 : virtual void clear();
101 : /// this constructor here is Value-aware
102 : Grid(const std::string& funcl, const std::vector<Value*> & args, const std::vector<std::string> & gmin,
103 : const std::vector<std::string> & gmax, const std::vector<unsigned> & nbin, bool dospline,
104 : bool usederiv, bool doclear=true);
105 : /// this constructor here is not Value-aware
106 : Grid(const std::string& funcl, const std::vector<std::string> &names, const std::vector<std::string> & gmin,
107 : const std::vector<std::string> & gmax, const std::vector<unsigned> & nbin, bool dospline,
108 : bool usederiv, bool doclear, const std::vector<bool> &isperiodic, const std::vector<std::string> &pmin,
109 : const std::vector<std::string> &pmax );
110 : /// this is the real initializator
111 : void Init(const std::string & funcl, const std::vector<std::string> &names, const std::vector<std::string> & gmin,
112 : const std::vector<std::string> & gmax, const std::vector<unsigned> & nbin, bool dospline, bool usederiv,
113 : bool doclear, const std::vector<bool> &isperiodic, const std::vector<std::string> &pmin, const std::vector<std::string> &pmax);
114 : /// get lower boundary
115 : std::vector<std::string> getMin() const;
116 : /// get upper boundary
117 : std::vector<std::string> getMax() const;
118 : /// get bin size
119 : std::vector<double> getDx() const;
120 : /// get bin volume
121 : double getBinVolume() const;
122 : /// get number of bins
123 : std::vector<unsigned> getNbin() const;
124 : /// get if periodic
125 : std::vector<bool> getIsPeriodic() const;
126 : /// get grid dimension
127 : unsigned getDimension() const;
128 : /// get argument names of this grid
129 : std::vector<std::string> getArgNames() const;
130 :
131 : /// methods to handle grid indices
132 : std::vector<unsigned> getIndices(index_t index) const;
133 : std::vector<unsigned> getIndices(const std::vector<double> & x) const;
134 : index_t getIndex(const std::vector<unsigned> & indices) const;
135 : index_t getIndex(const std::vector<double> & x) const;
136 : std::vector<double> getPoint(index_t index) const;
137 : std::vector<double> getPoint(const std::vector<unsigned> & indices) const;
138 : std::vector<double> getPoint(const std::vector<double> & x) const;
139 : /// faster versions relying on preallocated vectors
140 : void getPoint(index_t index,std::vector<double> & point) const;
141 : void getPoint(const std::vector<unsigned> & indices,std::vector<double> & point) const;
142 : void getPoint(const std::vector<double> & x,std::vector<double> & point) const;
143 :
144 : /// get neighbors
145 : std::vector<index_t> getNeighbors(index_t index,const std::vector<unsigned> & neigh) const;
146 : std::vector<index_t> getNeighbors(const std::vector<unsigned> & indices,const std::vector<unsigned> & neigh) const;
147 : std::vector<index_t> getNeighbors(const std::vector<double> & x,const std::vector<unsigned> & neigh) const;
148 :
149 : /// write header for grid file
150 : void writeHeader(OFile& file);
151 :
152 : /// read grid from file
153 : static Grid* create(const std::string&,const std::vector<Value*>&,IFile&,bool,bool,bool);
154 : /// read grid from file and check boundaries are what is expected from input
155 : static Grid* create(const std::string&,const std::vector<Value*>&, IFile&,
156 : const std::vector<std::string>&,const std::vector<std::string>&,
157 : const std::vector<unsigned>&,bool,bool,bool);
158 : /// get grid size
159 : virtual index_t getSize() const;
160 : /// get grid value
161 : virtual double getValue(index_t index) const;
162 : virtual double getValue(const std::vector<unsigned> & indices) const;
163 : virtual double getValue(const std::vector<double> & x) const;
164 : /// get minimum value
165 : virtual double getMinValue() const;
166 : /// get maximum value
167 : virtual double getMaxValue() const;
168 : /// get grid value and derivatives
169 : virtual double getValueAndDerivatives(index_t index, std::vector<double>& der) const ;
170 : virtual double getValueAndDerivatives(const std::vector<unsigned> & indices, std::vector<double>& der) const;
171 : virtual double getValueAndDerivatives(const std::vector<double> & x, std::vector<double>& der) const;
172 : /// Get the difference from the contour
173 : double getDifferenceFromContour(const std::vector<double> & x, std::vector<double>& der) const ;
174 : /// Find a set of points on a contour in the function
175 : void findSetOfPointsOnContour(const double& target, const std::vector<bool>& nosearch, unsigned& npoints, std::vector<std::vector<double> >& points );
176 :
177 : /// set grid value
178 : virtual void setValue(index_t index, double value);
179 : virtual void setValue(const std::vector<unsigned> & indices, double value);
180 : /// set grid value and derivatives
181 : virtual void setValueAndDerivatives(index_t index, double value, std::vector<double>& der);
182 : virtual void setValueAndDerivatives(const std::vector<unsigned> & indices, double value, std::vector<double>& der);
183 : /// add to grid value
184 : virtual void addValue(index_t index, double value);
185 : virtual void addValue(const std::vector<unsigned> & indices, double value);
186 : /// add to grid value and derivatives
187 : virtual void addValueAndDerivatives(index_t index, double value, std::vector<double>& der);
188 : virtual void addValueAndDerivatives(const std::vector<unsigned> & indices, double value, std::vector<double>& der);
189 : /// Scale all grid values and derivatives by a constant factor
190 : virtual void scaleAllValuesAndDerivatives( const double& scalef );
191 : /// Takes the scalef times the logarithm of all grid values and derivatives
192 : virtual void logAllValuesAndDerivatives( const double& scalef );
193 : /// Set the minimum value of the grid to zero and translates accordingly
194 : virtual void setMinToZero();
195 : /// apply function: takes pointer to function that accepts a double and apply
196 : virtual void applyFunctionAllValuesAndDerivatives( double (*func)(double val), double (*funcder)(double valder) );
197 : /// add a kernel function to the grid
198 : void addKernel( const KernelFunctions& kernel );
199 :
200 : /// dump grid on file
201 : virtual void writeToFile(OFile&);
202 : /// dump grid to gaussian cube file
203 : void writeCubeFile(OFile&, const double& lunit);
204 :
205 95 : virtual ~Grid() {}
206 :
207 : /// project a high dimensional grid onto a low dimensional one: this should be changed at some time
208 : /// to enable many types of weighting
209 : Grid project( const std::vector<std::string> & proj, WeightBase *ptr2obj );
210 : void projectOnLowDimension(double &val, std::vector<int> &varHigh, WeightBase* ptr2obj );
211 : /// set output format
212 8 : void setOutputFmt(const std::string & ss) {fmt_=ss;}
213 : /// Integrate the function calculated on the grid
214 : double integrate( std::vector<unsigned>& npoints );
215 : ///
216 : void mpiSumValuesAndDerivatives( Communicator& comm );
217 : };
218 :
219 :
220 : class SparseGrid : public Grid
221 : {
222 :
223 : std::map<index_t,double> map_;
224 : typedef std::map<index_t,double>::const_iterator iterator;
225 : std::map< index_t,std::vector<double> > der_;
226 : typedef std::map<index_t,std::vector<double> >::const_iterator iterator_der;
227 :
228 : protected:
229 : void clear();
230 :
231 : public:
232 10 : SparseGrid(const std::string& funcl, const std::vector<Value*> & args, const std::vector<std::string> & gmin,
233 : const std::vector<std::string> & gmax,
234 : const std::vector<unsigned> & nbin, bool dospline, bool usederiv):
235 10 : Grid(funcl,args,gmin,gmax,nbin,dospline,usederiv,false) {}
236 :
237 : index_t getSize() const;
238 : index_t getMaxSize() const;
239 :
240 : /// this is to access to Grid:: version of these methods (allowing overloading of virtual methods)
241 : using Grid::getValue;
242 : using Grid::getValueAndDerivatives;
243 : using Grid::setValue;
244 : using Grid::setValueAndDerivatives;
245 : using Grid::addValue;
246 : using Grid::addValueAndDerivatives;
247 :
248 : /// get grid value
249 : double getValue(index_t index) const;
250 : /// get grid value and derivatives
251 : double getValueAndDerivatives(index_t index, std::vector<double>& der) const;
252 :
253 : /// set grid value
254 : void setValue(index_t index, double value);
255 : /// set grid value and derivatives
256 : void setValueAndDerivatives(index_t index, double value, std::vector<double>& der);
257 : /// add to grid value
258 : void addValue(index_t index, double value);
259 : /// add to grid value and derivatives
260 : void addValueAndDerivatives(index_t index, double value, std::vector<double>& der);
261 :
262 : /// dump grid on file
263 : void writeToFile(OFile&);
264 :
265 10 : virtual ~SparseGrid() {}
266 : };
267 : }
268 :
269 : #endif
|
