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 : #ifndef __PLUMED_ves_VesBias_h
23 : #define __PLUMED_ves_VesBias_h
24 :
25 : #include "CoeffsVector.h"
26 : #include "CoeffsMatrix.h"
27 :
28 : #include "core/ActionPilot.h"
29 : #include "core/ActionWithValue.h"
30 : #include "core/ActionWithArguments.h"
31 : #include "bias/Bias.h"
32 :
33 : #include <vector>
34 : #include <string>
35 : #include <cmath>
36 :
37 :
38 : #define PLUMED_VES_VESBIAS_INIT(ao) Action(ao),VesBias(ao)
39 :
40 : namespace PLMD {
41 :
42 : class Value;
43 :
44 : namespace ves {
45 :
46 : class CoeffsVector;
47 : class CoeffsMatrix;
48 : class BasisFunctions;
49 : class Optimizer;
50 : class TargetDistribution;
51 : class FermiSwitchingFunction;
52 :
53 : /**
54 : \ingroup INHERIT
55 : Abstract base class for implementing biases the extents the normal Bias.h class
56 : to include functions related to the variational approach.
57 : */
58 :
59 : class VesBias:
60 : public bias::Bias
61 : {
62 : private:
63 : unsigned int ncoeffssets_;
64 : std::vector<CoeffsVector*> coeffs_pntrs_;
65 : std::vector<CoeffsVector*> targetdist_averages_pntrs_;
66 : std::vector<CoeffsVector*> gradient_pntrs_;
67 : std::vector<CoeffsMatrix*> hessian_pntrs_;
68 : std::vector<std::vector<double> > sampled_averages;
69 : std::vector<std::vector<double> > sampled_cross_averages;
70 : bool use_multiple_coeffssets_;
71 : //
72 : std::vector<std::string> coeffs_fnames;
73 : //
74 : size_t ncoeffs_total_;
75 : //
76 : Optimizer* optimizer_pntr_;
77 : bool optimize_coeffs_;
78 : //
79 : bool compute_hessian_;
80 : bool diagonal_hessian_;
81 : //
82 : std::vector<unsigned int> aver_counters;
83 : //
84 : double kbt_;
85 : //
86 : std::vector<TargetDistribution*> targetdist_pntrs_;
87 : bool dynamic_targetdist_;
88 : //
89 : std::vector<unsigned int> grid_bins_;
90 : std::vector<double> grid_min_;
91 : std::vector<double> grid_max_;
92 : //
93 : std::string bias_filename_;
94 : std::string fes_filename_;
95 : std::string targetdist_filename_;
96 : std::string targetdist_averages_filename_;
97 : std::string coeffs_id_prefix_;
98 : //
99 : std::string bias_file_fmt_;
100 : std::string fes_file_fmt_;
101 : std::string targetdist_file_fmt_;
102 : std::string targetdist_restart_file_fmt_;
103 : //
104 : bool filenames_have_iteration_number_;
105 : //
106 : bool bias_fileoutput_active_;
107 : bool fes_fileoutput_active_;
108 : bool fesproj_fileoutput_active_;
109 : bool dynamic_targetdist_fileoutput_active_;
110 : bool static_targetdist_fileoutput_active_;
111 : //
112 :
113 : bool bias_cutoff_active_;
114 : double bias_cutoff_value_;
115 : double bias_current_max_value;
116 : FermiSwitchingFunction* bias_cutoff_swfunc_pntr_;
117 : //
118 : std::vector< std::vector<std::string> > projection_args_;
119 : //
120 : bool calc_reweightfactor_;
121 : private:
122 : void initializeCoeffs(CoeffsVector*);
123 : std::vector<double> computeCovarianceFromAverages(const unsigned int) const;
124 : void multiSimSumAverages(const unsigned int, const double walker_weight=1.0);
125 : protected:
126 : //
127 : void checkThatTemperatureIsGiven();
128 : //
129 : void addCoeffsSet(const std::vector<std::string>&,const std::vector<unsigned int>&);
130 : void addCoeffsSet(std::vector<Value*>&,std::vector<BasisFunctions*>&);
131 : void addCoeffsSet(CoeffsVector*);
132 : //
133 : std::string getCoeffsSetLabelString(const std::string&, const unsigned int coeffs_id = 0) const;
134 : void clearCoeffsPntrsVector() {coeffs_pntrs_.clear();}
135 : void addToSampledAverages(const std::vector<double>&, const unsigned int c_id = 0);
136 : void setTargetDistAverages(const std::vector<double>&, const unsigned int coeffs_id = 0);
137 : void setTargetDistAverages(const CoeffsVector&, const unsigned int coeffs_id= 0);
138 : void setTargetDistAveragesToZero(const unsigned int coeffs_id= 0);
139 : //
140 : bool readCoeffsFromFiles();
141 : //
142 : template<class T>
143 : bool parseMultipleValues(const std::string&, std::vector<T>&, unsigned int);
144 : template<class T>
145 : bool parseMultipleValues(const std::string&, std::vector<T>&, unsigned int, const T&);
146 : //
147 : public:
148 : static void registerKeywords(Keywords&);
149 : explicit VesBias(const ActionOptions&ao);
150 : ~VesBias();
151 : //
152 : void apply();
153 : //
154 : static void useInitialCoeffsKeywords(Keywords&);
155 : static void useTargetDistributionKeywords(Keywords&);
156 : static void useMultipleTargetDistributionKeywords(Keywords&);
157 : static void useGridBinKeywords(Keywords&);
158 : static void useGridLimitsKeywords(Keywords&);
159 : static void useBiasCutoffKeywords(Keywords&);
160 : static void useProjectionArgKeywords(Keywords&);
161 : static void useReweightFactorKeywords(Keywords&);
162 : //
163 251 : std::vector<CoeffsVector*> getCoeffsPntrs() const {return coeffs_pntrs_;}
164 77 : std::vector<CoeffsVector*> getTargetDistAveragesPntrs() const {return targetdist_averages_pntrs_;}
165 77 : std::vector<CoeffsVector*> getGradientPntrs()const {return gradient_pntrs_;}
166 77 : std::vector<CoeffsMatrix*> getHessianPntrs() const {return hessian_pntrs_;}
167 123 : std::vector<TargetDistribution*> getTargetDistributionPntrs() const {return targetdist_pntrs_;}
168 : //
169 82 : CoeffsVector* getCoeffsPntr(const unsigned int coeffs_id = 0) const {return coeffs_pntrs_[coeffs_id];}
170 : CoeffsVector* getTargetDistAveragesPntr(const unsigned int coeffs_id = 0) const {return targetdist_averages_pntrs_[coeffs_id];}
171 : CoeffsVector* getGradientPntr(const unsigned int coeffs_id = 0)const {return gradient_pntrs_[coeffs_id];}
172 : CoeffsMatrix* getHessianPntr(const unsigned int coeffs_id = 0) const {return hessian_pntrs_[coeffs_id];}
173 : //
174 82 : unsigned int getNumberOfTargetDistributionPntrs() const {return targetdist_pntrs_.size();}
175 : //
176 3095 : size_t numberOfCoeffs(const unsigned int coeffs_id = 0) const {return coeffs_pntrs_[coeffs_id]->numberOfCoeffs();}
177 : size_t totalNumberOfCoeffs() const {return ncoeffs_total_;}
178 3 : unsigned int numberOfCoeffsSets() const {return ncoeffssets_;}
179 77 : double getKbT() const {return kbt_;}
180 : double getBeta() const;
181 : //
182 : CoeffsVector& Coeffs(const unsigned int coeffs_id = 0) const {return *coeffs_pntrs_[coeffs_id];}
183 2033 : CoeffsVector& TargetDistAverages(const unsigned int coeffs_id = 0) const {return *targetdist_averages_pntrs_[coeffs_id];}
184 1720 : CoeffsVector& Gradient(const unsigned int coeffs_id = 0) const {return *gradient_pntrs_[coeffs_id];}
185 2490 : CoeffsMatrix& Hessian(const unsigned int coeffs_id = 0) const {return *hessian_pntrs_[coeffs_id];}
186 : //
187 : size_t getCoeffsIndex(const std::vector<unsigned int>& indices, const unsigned int coeffs_id = 0) const;
188 : std::vector<unsigned int> getCoeffsIndices(const size_t index, const unsigned int coeffs_id = 0) const;
189 : size_t getHessianIndex(const size_t index1, const size_t index2, const unsigned int coeffs_id = 0) const;
190 : //
191 : bool computeHessian() const {return compute_hessian_;}
192 : bool diagonalHessian() const {return diagonal_hessian_;}
193 : //
194 1468 : bool optimizeCoeffs() const {return optimize_coeffs_;}
195 1330 : Optimizer* getOptimizerPntr() const {return optimizer_pntr_;}
196 : bool useMultipleWalkers() const;
197 : //
198 : unsigned int getIterationCounter() const;
199 : //
200 : void updateGradientAndHessian(const bool);
201 : void clearGradientAndHessian() {};
202 : //
203 0 : virtual void updateTargetDistributions() {};
204 0 : virtual void restartTargetDistributions() {};
205 : //
206 : void linkOptimizer(Optimizer*);
207 : void enableHessian(const bool diagonal_hessian=true);
208 : void disableHessian();
209 : //
210 : void enableMultipleCoeffsSets() {use_multiple_coeffssets_=true;}
211 : //
212 37 : void enableDynamicTargetDistribution() {dynamic_targetdist_=true;}
213 : void disableDynamicTargetDistribution() {dynamic_targetdist_=false;}
214 245 : bool dynamicTargetDistribution() const {return dynamic_targetdist_;}
215 : //
216 82 : std::vector<unsigned int> getGridBins() const {return grid_bins_;}
217 : void setGridBins(const std::vector<unsigned int>&);
218 : void setGridBins(const unsigned int);
219 : std::vector<double> getGridMax() const {return grid_max_;}
220 : void setGridMax(const std::vector<double>&);
221 : std::vector<double> getGridMin() const {return grid_min_;}
222 : void setGridMin(const std::vector<double>&);
223 : //
224 539 : bool filenamesIncludeIterationNumber() const {return filenames_have_iteration_number_;}
225 3 : void enableIterationNumberInFilenames() {filenames_have_iteration_number_=true;}
226 : //
227 : std::string getIterationFilenameSuffix() const;
228 : std::string getCoeffsSetFilenameSuffix(const unsigned int coeffs_id) const;
229 : std::string getCurrentOutputFilename(const std::string&, const std::string& suffix="") const;
230 : std::string getBiasOutputFilename() const {return bias_filename_;}
231 : std::string getCurrentBiasOutputFilename(const std::string& suffix="") const;
232 : std::string getFesOutputFilename() const {return fes_filename_;}
233 : std::string getCurrentFesOutputFilename(const std::string& suffix="") const;
234 : std::string getTargetDistOutputFilename() const {return targetdist_filename_;}
235 : std::string getCurrentTargetDistOutputFilename(const std::string& suffix="") const;
236 : //
237 79 : void enableBiasFileOutput() {bias_fileoutput_active_=true;}
238 : void disableBiasFileOutput() {bias_fileoutput_active_=false;}
239 : bool isBiasFileOutputActive() const {return bias_fileoutput_active_;}
240 : std::string getBiasFileFmt() const {return bias_file_fmt_;}
241 : //
242 79 : void enableFesFileOutput() {fes_fileoutput_active_=true;}
243 : void disableFesFileOutput() {fes_fileoutput_active_=false;}
244 : bool isFesFileOutputActive() const {return fes_fileoutput_active_;}
245 : std::string getFesFileFmt() const {return fes_file_fmt_;}
246 : //
247 17 : void enableFesProjFileOutput() {fesproj_fileoutput_active_=true;}
248 : void disableFesFileProjOutput() {fesproj_fileoutput_active_=false;}
249 : bool isFesProjFileOutputActive() const {return fesproj_fileoutput_active_;}
250 : //
251 36 : void enableDynamicTargetDistFileOutput() {dynamic_targetdist_fileoutput_active_=true;}
252 : void disableDynamicTargetDistFileOutput() {dynamic_targetdist_fileoutput_active_=false;}
253 : bool isDynamicTargetDistFileOutputActive() const {return dynamic_targetdist_fileoutput_active_;}
254 : std::string getTargetDistFileFmt() const {return targetdist_file_fmt_;}
255 : std::string getTargetDistRestartFileFmt() const {return targetdist_restart_file_fmt_;}
256 : //
257 : void enableStaticTargetDistFileOutput() {static_targetdist_fileoutput_active_=true;}
258 42 : void disableStaticTargetDistFileOutput() {static_targetdist_fileoutput_active_=false;}
259 47 : bool isStaticTargetDistFileOutputActive() const {return static_targetdist_fileoutput_active_;}
260 : //
261 : std::vector< std::vector<std::string> > getProjectionArguments() const {return projection_args_;}
262 112 : std::vector<std::string> getProjectionArgument(unsigned int i) const {return projection_args_[i];}
263 122 : unsigned int getNumberOfProjectionArguments() const {return projection_args_.size();}
264 : //
265 : void setupBiasCutoff(const double, const double);
266 1456547 : bool biasCutoffActive() const {return bias_cutoff_active_;}
267 45 : double getBiasCutoffValue() const {return bias_cutoff_value_;}
268 467 : void setCurrentBiasMaxValue(const double max_value) {bias_current_max_value=max_value;}
269 : double getCurrentBiasMaxValue() const {return bias_current_max_value;}
270 : double getBiasCutoffSwitchingFunction(const double, double&) const;
271 : double getBiasCutoffSwitchingFunction(const double) const;
272 : void applyBiasCutoff(double&, std::vector<double>&) const;
273 : void applyBiasCutoff(double&, std::vector<double>&, std::vector<double>&) const;
274 : //
275 : OFile* getOFile(const std::string& filename, const bool multi_sim_single_file=false, const bool enforce_backup=true);
276 : //
277 0 : virtual void setupBiasFileOutput() {};
278 0 : virtual void writeBiasToFile() {};
279 0 : virtual void resetBiasFileOutput() {};
280 : //
281 0 : virtual void setupFesFileOutput() {};
282 0 : virtual void writeFesToFile() {};
283 0 : virtual void resetFesFileOutput() {};
284 : //
285 0 : virtual void setupFesProjFileOutput() {};
286 0 : virtual void writeFesProjToFile() {};
287 0 : virtual void resetFesProjFileOutput() {};
288 : //
289 39 : virtual void setupTargetDistFileOutput() {};
290 0 : virtual void writeTargetDistToFile() {};
291 0 : virtual void resetTargetDistFileOutput() {};
292 : //
293 9 : virtual void setupTargetDistProjFileOutput() {};
294 0 : virtual void writeTargetDistProjToFile() {};
295 0 : virtual void resetTargetDistProjFileOutput() {};
296 : //
297 : void updateReweightFactor();
298 : virtual double calculateReweightFactor() const;
299 0 : bool isReweightFactorCalculated() const {return calc_reweightfactor_;}
300 : };
301 :
302 :
303 : inline
304 : size_t VesBias::getCoeffsIndex(const std::vector<unsigned int>& indices, const unsigned int coeffs_id) const {return coeffs_pntrs_[coeffs_id]->getIndex(indices);}
305 :
306 : inline
307 : std::vector<unsigned int> VesBias::getCoeffsIndices(const size_t index, const unsigned int coeffs_id) const {return coeffs_pntrs_[coeffs_id]->getIndices(index);}
308 :
309 : inline
310 64954 : size_t VesBias::getHessianIndex(const size_t index1, const size_t index2, const unsigned int coeffs_id) const {return hessian_pntrs_[coeffs_id]->getMatrixIndex(index1,index2);}
311 :
312 :
313 : inline
314 1204 : double VesBias::getBeta() const {
315 1204 : plumed_massert(kbt_!=0.0,"you are requesting beta=1/(kB*T) when kB*T has not been defined. You need to give the temperature using the TEMP keyword as the MD engine does not pass it to PLUMED.");
316 1204 : return 1.0/kbt_;
317 : }
318 :
319 :
320 : inline
321 : std::string VesBias::getCurrentBiasOutputFilename(const std::string& suffix) const {
322 170 : return getCurrentOutputFilename(bias_filename_,suffix);
323 : }
324 :
325 :
326 : inline
327 : std::string VesBias::getCurrentFesOutputFilename(const std::string& suffix) const {
328 197 : return getCurrentOutputFilename(fes_filename_,suffix);
329 : }
330 :
331 :
332 : inline
333 : double VesBias::getBiasCutoffSwitchingFunction(const double bias) const {
334 : double dummy=0.0;
335 : return getBiasCutoffSwitchingFunction(bias,dummy);
336 : }
337 :
338 :
339 : inline
340 600 : void VesBias::applyBiasCutoff(double& bias, std::vector<double>& forces) const {
341 600 : std::vector<double> dummy(0);
342 600 : applyBiasCutoff(bias,forces,dummy);
343 600 : }
344 :
345 :
346 : inline
347 663 : void VesBias::applyBiasCutoff(double& bias, std::vector<double>& forces, std::vector<double>& coeffsderivs_values) const {
348 663 : double deriv_factor_sf=0.0;
349 663 : double value_sf = getBiasCutoffSwitchingFunction(bias,deriv_factor_sf);
350 663 : bias *= value_sf;
351 3315 : for(unsigned int i=0; i<forces.size(); i++) {
352 663 : forces[i] *= deriv_factor_sf;
353 : }
354 : //
355 3531 : for(unsigned int i=0; i<coeffsderivs_values.size(); i++) {
356 735 : coeffsderivs_values[i] *= deriv_factor_sf;
357 : }
358 663 : }
359 :
360 :
361 : inline
362 770 : std::vector<double> VesBias::computeCovarianceFromAverages(const unsigned int c_id) const {
363 : size_t ncoeffs = numberOfCoeffs(c_id);
364 1540 : std::vector<double> covariance(sampled_cross_averages[c_id].size(),0.0);
365 : // diagonal part
366 61430 : for(size_t i=0; i<ncoeffs; i++) {
367 : size_t midx = getHessianIndex(i,i,c_id);
368 90990 : covariance[midx] = sampled_cross_averages[c_id][midx] - sampled_averages[c_id][i]*sampled_averages[c_id][i];
369 : }
370 770 : if(!diagonal_hessian_) {
371 0 : for(size_t i=0; i<ncoeffs; i++) {
372 0 : for(size_t j=(i+1); j<ncoeffs; j++) {
373 : size_t midx = getHessianIndex(i,j,c_id);
374 0 : covariance[midx] = sampled_cross_averages[c_id][midx] - sampled_averages[c_id][i]*sampled_averages[c_id][j];
375 : }
376 : }
377 : }
378 770 : return covariance;
379 : }
380 :
381 :
382 : template<class T>
383 164 : bool VesBias::parseMultipleValues(const std::string& keyword, std::vector<T>& values, unsigned int nvalues) {
384 164 : plumed_assert(nvalues>0);
385 164 : plumed_assert(values.size()==0);
386 : bool identical_values=false;
387 : //
388 164 : parseVector(keyword,values);
389 164 : if(values.size()==1 && nvalues>1) {
390 0 : values.resize(nvalues,values[0]);
391 : identical_values=true;
392 : }
393 164 : if(values.size()>0 && values.size()!=nvalues) {
394 0 : std::string s1; Tools::convert(nvalues,s1);
395 0 : plumed_merror("Error in " + keyword + " keyword: either give 1 common parameter value or " + s1 + " separate parameter values");
396 : }
397 164 : return identical_values;
398 : }
399 :
400 : template<class T>
401 82 : bool VesBias::parseMultipleValues(const std::string& keyword, std::vector<T>& values, unsigned int nvalues, const T& default_value) {
402 82 : bool identical_values = parseMultipleValues(keyword,values,nvalues);
403 82 : if(values.size()==0) {
404 0 : values.resize(nvalues,default_value);
405 : identical_values=true;
406 : }
407 82 : return identical_values;
408 : }
409 :
410 :
411 : }
412 : }
413 :
414 : #endif
|
