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