Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2021 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 "Optimizer.h"
24 : #include "CoeffsVector.h"
25 : #include "CoeffsMatrix.h"
26 : #include "VesBias.h"
27 : #include "VesTools.h"
28 :
29 : #include "tools/Exception.h"
30 : #include "core/PlumedMain.h"
31 : #include "core/ActionSet.h"
32 : #include "tools/Communicator.h"
33 : #include "tools/File.h"
34 : #include "tools/FileBase.h"
35 :
36 : namespace PLMD {
37 : namespace ves {
38 :
39 79 : Optimizer::Optimizer(const ActionOptions&ao):
40 : Action(ao),
41 : ActionPilot(ao),
42 : ActionWithValue(ao),
43 158 : description_("Undefined"),
44 79 : type_("Undefined"),
45 79 : stepsizes_(0),
46 79 : current_stepsizes(0),
47 79 : fixed_stepsize_(true),
48 79 : iter_counter(0),
49 79 : use_hessian_(false),
50 79 : diagonal_hessian_(true),
51 79 : monitor_instantaneous_gradient_(false),
52 79 : use_mwalkers_mpi_(false),
53 79 : mwalkers_mpi_single_files_(true),
54 79 : dynamic_targetdists_(0),
55 79 : ustride_targetdist_(0),
56 79 : ustride_reweightfactor_(0),
57 79 : coeffssetid_prefix_(""),
58 79 : coeffs_wstride_(100),
59 79 : coeffs_output_fmt_(""),
60 79 : gradient_wstride_(100),
61 79 : gradient_output_fmt_(""),
62 79 : hessian_wstride_(100),
63 79 : hessianOFiles_(0),
64 79 : hessian_output_fmt_(""),
65 79 : targetdist_averages_wstride_(0),
66 79 : targetdist_averages_output_fmt_(""),
67 79 : nbiases_(0),
68 79 : bias_pntrs_(0),
69 79 : ncoeffssets_(0),
70 79 : coeffs_pntrs_(0),
71 79 : gradient_pntrs_(0),
72 79 : hessian_pntrs_(0),
73 79 : targetdist_averages_pntrs_(0),
74 79 : identical_coeffs_shape_(true),
75 79 : bias_output_active_(false),
76 79 : bias_output_stride_(0),
77 79 : fes_output_active_(false),
78 79 : fes_output_stride_(0),
79 79 : fesproj_output_active_(false),
80 79 : fesproj_output_stride_(0),
81 79 : targetdist_output_active_(false),
82 79 : targetdist_output_stride_(0),
83 79 : targetdist_proj_output_active_(false),
84 79 : targetdist_proj_output_stride_(0),
85 316 : isFirstStep(true) {
86 79 : std::vector<std::string> bias_labels(0);
87 158 : parseVector("BIAS",bias_labels);
88 79 : plumed_massert(bias_labels.size()>0,"problem with BIAS keyword");
89 79 : nbiases_ = bias_labels.size();
90 : //
91 79 : std::string error_msg = "";
92 158 : bias_pntrs_ = VesTools::getPointersFromLabels<VesBias*>(bias_labels,plumed.getActionSet(),error_msg);
93 79 : if(error_msg.size()>0) {
94 0 : plumed_merror("Error in keyword BIAS of "+getName()+": "+error_msg);
95 : }
96 :
97 164 : for(unsigned int i=0; i<bias_pntrs_.size(); i++) {
98 85 : bias_pntrs_[i]->linkOptimizer(this);
99 : //
100 85 : std::vector<CoeffsVector*> pntrs_coeffs = bias_pntrs_[i]->getCoeffsPntrs();
101 85 : std::vector<CoeffsVector*> pntrs_gradient = bias_pntrs_[i]->getGradientPntrs();
102 85 : std::vector<CoeffsVector*> pntrs_targetdist_averages = bias_pntrs_[i]->getTargetDistAveragesPntrs();
103 85 : plumed_massert(pntrs_coeffs.size()==pntrs_gradient.size(),"something wrong in the coefficients and gradient passed from VES bias");
104 85 : plumed_massert(pntrs_coeffs.size()==pntrs_targetdist_averages.size(),"something wrong in the coefficients and target distribution averages passed from VES bias");
105 170 : for(unsigned int k=0; k<pntrs_coeffs.size(); k++) {
106 85 : plumed_massert(pntrs_coeffs[k] != NULL,"some coefficient is not linked correctly");
107 85 : plumed_massert(pntrs_gradient[k] != NULL,"some gradient is not linked correctly");
108 85 : plumed_massert(pntrs_targetdist_averages[k] != NULL,"some target distribution average is not linked correctly");
109 : pntrs_coeffs[k]->turnOnIterationCounter();
110 85 : coeffs_pntrs_.push_back(pntrs_coeffs[k]);
111 85 : pntrs_gradient[k]->turnOnIterationCounter();
112 85 : gradient_pntrs_.push_back(pntrs_gradient[k]);
113 85 : pntrs_targetdist_averages[k]->turnOnIterationCounter();
114 85 : targetdist_averages_pntrs_.push_back(pntrs_targetdist_averages[k]);
115 : //
116 85 : auto aux_coeffs_tmp = Tools::make_unique<CoeffsVector>(*pntrs_coeffs[k]);
117 85 : std::string aux_label = pntrs_coeffs[k]->getLabel();
118 85 : if(aux_label.find("coeffs")!=std::string::npos) {
119 170 : aux_label.replace(aux_label.find("coeffs"), std::string("coeffs").length(), "aux_coeffs");
120 : } else {
121 : aux_label += "_aux";
122 : }
123 85 : aux_coeffs_tmp->setLabels(aux_label);
124 85 : aux_coeffs_pntrs_.emplace_back(std::move(aux_coeffs_tmp));
125 85 : AuxCoeffs(i).setValues( Coeffs(i) );
126 85 : }
127 : }
128 79 : ncoeffssets_ = coeffs_pntrs_.size();
129 79 : plumed_massert(aux_coeffs_pntrs_.size()==ncoeffssets_,"problems in linking aux coefficients");
130 79 : plumed_massert(gradient_pntrs_.size()==ncoeffssets_,"problems in linking gradients");
131 79 : setAllCoeffsSetIterationCounters();
132 :
133 :
134 : //
135 79 : identical_coeffs_shape_ = true;
136 85 : for(unsigned int i=1; i<ncoeffssets_; i++) {
137 6 : if(!coeffs_pntrs_[0]->sameShape(*coeffs_pntrs_[i])) {
138 0 : identical_coeffs_shape_ = false;
139 0 : break;
140 : }
141 : }
142 : //
143 79 : if(keywords.exists("STEPSIZE")) {
144 77 : plumed_assert(!keywords.exists("INITIAL_STEPSIZE"));
145 77 : fixed_stepsize_=true;
146 77 : parseMultipleValues("STEPSIZE",stepsizes_);
147 77 : setCurrentStepSizes(stepsizes_);
148 : }
149 79 : if(keywords.exists("INITIAL_STEPSIZE")) {
150 1 : plumed_assert(!keywords.exists("STEPSIZE"));
151 1 : fixed_stepsize_=false;
152 1 : parseMultipleValues("INITIAL_STEPSIZE",stepsizes_);
153 1 : setCurrentStepSizes(stepsizes_);
154 : }
155 : //
156 79 : if(ncoeffssets_==1) {
157 76 : log.printf(" optimizing VES bias %s with label %s: \n",bias_pntrs_[0]->getName().c_str(),bias_pntrs_[0]->getLabel().c_str());
158 76 : log.printf(" KbT: %f\n",bias_pntrs_[0]->getKbT());
159 76 : log.printf(" number of coefficients: %zu\n",coeffs_pntrs_[0]->numberOfCoeffs());
160 76 : if(stepsizes_.size()>0) {
161 75 : if(fixed_stepsize_) {
162 74 : log.printf(" using a constant step size of %f\n",stepsizes_[0]);
163 : } else {
164 1 : log.printf(" using an initial step size of %f\n",stepsizes_[0]);
165 : }
166 : }
167 : } else {
168 3 : log.printf(" optimizing %u coefficient sets from following %u VES biases:\n",ncoeffssets_,nbiases_);
169 12 : for(unsigned int i=0; i<nbiases_; i++) {
170 9 : log.printf(" %s of type %s (KbT: %f) \n",bias_pntrs_[i]->getLabel().c_str(),bias_pntrs_[i]->getName().c_str(),bias_pntrs_[i]->getKbT());
171 : }
172 : size_t tot_ncoeffs = 0;
173 12 : for(unsigned int i=0; i<ncoeffssets_; i++) {
174 9 : log.printf(" coefficient set %u: \n",i);
175 9 : log.printf(" used in bias %s (type %s)\n",coeffs_pntrs_[i]->getPntrToAction()->getLabel().c_str(),coeffs_pntrs_[i]->getPntrToAction()->getName().c_str());
176 9 : log.printf(" number of coefficients: %zu\n",coeffs_pntrs_[i]->numberOfCoeffs());
177 9 : if(stepsizes_.size()>0) {
178 9 : if(fixed_stepsize_) {
179 9 : log.printf(" using a constant step size of %f\n",stepsizes_[i]);
180 : } else {
181 0 : log.printf(" using an initial step size of %f\n",stepsizes_[i]);
182 : }
183 : }
184 9 : tot_ncoeffs += coeffs_pntrs_[i]->numberOfCoeffs();
185 : }
186 3 : log.printf(" total number of coefficients: %zu\n",tot_ncoeffs);
187 3 : if(identical_coeffs_shape_) {
188 3 : log.printf(" the indices shape is identical for all coefficient sets\n");
189 : } else {
190 0 : log.printf(" the indices shape differs between coefficient sets\n");
191 : }
192 : }
193 :
194 : //
195 79 : if(keywords.exists("FULL_HESSIAN")) {
196 0 : bool full_hessian=false;
197 0 : parseFlag("FULL_HESSIAN",full_hessian);
198 0 : diagonal_hessian_ = !full_hessian;
199 : }
200 : //
201 : bool mw_single_files = false;
202 79 : if(keywords.exists("MULTIPLE_WALKERS")) {
203 79 : parseFlag("MULTIPLE_WALKERS",use_mwalkers_mpi_);
204 79 : if(use_mwalkers_mpi_) {
205 : mw_single_files=true;
206 : }
207 : }
208 :
209 79 : int numwalkers=1;
210 79 : int walker_rank=0;
211 79 : if(comm.Get_rank()==0) {
212 70 : numwalkers = multi_sim_comm.Get_size();
213 70 : walker_rank = multi_sim_comm.Get_rank();
214 : }
215 79 : comm.Bcast(numwalkers,0);
216 79 : comm.Bcast(walker_rank,0);
217 79 : if(use_mwalkers_mpi_ && numwalkers==1) {
218 0 : plumed_merror("using the MULTIPLE_WALKERS keyword does not make sense if running the MD code with a single replica");
219 : }
220 79 : if(use_mwalkers_mpi_) {
221 12 : log.printf(" optimization performed using multiple walkers connected via MPI:\n");
222 12 : log.printf(" number of walkers: %d\n",numwalkers);
223 12 : log.printf(" walker number: %d\n",walker_rank);
224 12 : log.printf(" please see and cite ");
225 24 : log << plumed.cite("Raiteri, Laio, Gervasio, Micheletti, and Parrinello, J. Phys. Chem. B 110, 3533 (2006)");
226 12 : log.printf("\n");
227 : }
228 :
229 79 : dynamic_targetdists_.resize(nbiases_,false);
230 79 : if(keywords.exists("TARGETDIST_STRIDE")) {
231 : bool need_stride = false;
232 162 : for(unsigned int i=0; i<nbiases_; i++) {
233 84 : dynamic_targetdists_[i] = bias_pntrs_[i]->dynamicTargetDistribution();
234 84 : if(dynamic_targetdists_[i]) {
235 : need_stride = true;
236 : }
237 : }
238 78 : parse("TARGETDIST_STRIDE",ustride_targetdist_);
239 78 : if(need_stride && ustride_targetdist_==0) {
240 0 : plumed_merror("one of the biases has a dynamic target distribution so you need to give stride for updating it by using the TARGETDIST_STRIDE keyword");
241 : }
242 78 : if(!need_stride && ustride_targetdist_!=0) {
243 0 : plumed_merror("using the TARGETDIST_STRIDE keyword doesn't make sense as there is no dynamic target distribution to update");
244 : }
245 78 : if(ustride_targetdist_>0) {
246 38 : if(nbiases_==1) {
247 38 : log.printf(" the target distribution will be updated very %u coefficient iterations\n",ustride_targetdist_);
248 : } else {
249 0 : log.printf(" the target distribution will be updated very %u coefficient iterations for the following biases\n ",ustride_targetdist_);
250 0 : for(unsigned int i=0; i<nbiases_; i++) {
251 0 : log.printf("%s ",bias_pntrs_[i]->getLabel().c_str());
252 : }
253 0 : log.printf("\n");
254 : }
255 38 : log.printf(" See and cite ");
256 76 : log << plumed.cite("Valsson and Parrinello, J. Chem. Theory Comput. 11, 1996-2002 (2015)");
257 38 : log.printf("\n");
258 : }
259 : }
260 :
261 79 : if(keywords.exists("REWEIGHT_FACTOR_STRIDE")) {
262 : bool reweightfactor_calculated = false;
263 0 : for(unsigned int i=0; i<nbiases_; i++) {
264 0 : reweightfactor_calculated = bias_pntrs_[i]->isReweightFactorCalculated();
265 : }
266 0 : parse("REWEIGHT_FACTOR_STRIDE",ustride_reweightfactor_);
267 0 : if(ustride_reweightfactor_==0 && reweightfactor_calculated) {
268 0 : plumed_merror("the calculation of the reweight factor is enabled, You need to use the REWEIGHT_FACTOR_STRIDE keyword to specfiy how often it should be updated.");
269 : }
270 0 : if(ustride_reweightfactor_>0) {
271 0 : if(!reweightfactor_calculated) {
272 0 : plumed_merror("In order to use the REWEIGHT_FACTOR_STRIDE keyword you need to enable the calculation of the reweight factor in the VES bias by using the CALC_REWEIGHT_FACTOR flag.");
273 : }
274 0 : log.printf(" the reweight factor c(t) will be updated very %u coefficient iterations\n",ustride_reweightfactor_);
275 : }
276 : }
277 :
278 79 : if(keywords.exists("MONITOR_INSTANTANEOUS_GRADIENT")) {
279 158 : parseFlag("MONITOR_INSTANTANEOUS_GRADIENT",monitor_instantaneous_gradient_);
280 : }
281 :
282 79 : if(keywords.exists("MONITOR_AVERAGE_GRADIENT")) {
283 76 : bool monitor_aver_gradient = false;
284 76 : parseFlag("MONITOR_AVERAGE_GRADIENT",monitor_aver_gradient);
285 76 : if(monitor_aver_gradient) {
286 2 : unsigned int averaging_exp_decay=0;
287 4 : parse("MONITOR_AVERAGES_GRADIENT_EXP_DECAY",averaging_exp_decay);
288 : aver_gradient_pntrs_.clear();
289 4 : for(unsigned int i=0; i<ncoeffssets_; i++) {
290 2 : auto aver_gradient_tmp = Tools::make_unique<CoeffsVector>(*gradient_pntrs_[i]);
291 2 : aver_gradient_tmp->clear();
292 : std::string aver_grad_label = aver_gradient_tmp->getLabel();
293 2 : if(aver_grad_label.find("gradient")!=std::string::npos) {
294 4 : aver_grad_label.replace(aver_grad_label.find("gradient"), std::string("gradient").length(), "aver_gradient");
295 : } else {
296 : aver_grad_label += "_aver";
297 : }
298 2 : aver_gradient_tmp->setLabels(aver_grad_label);
299 2 : if(averaging_exp_decay>0) {
300 : aver_gradient_tmp->setupExponentiallyDecayingAveraging(averaging_exp_decay);
301 : }
302 2 : aver_gradient_pntrs_.emplace_back(std::move(aver_gradient_tmp));
303 2 : }
304 : }
305 : }
306 :
307 :
308 79 : if(ncoeffssets_>1) {
309 : coeffssetid_prefix_="c-";
310 3 : if(keywords.exists("COEFFS_SET_ID_PREFIX")) {
311 6 : parse("COEFFS_SET_ID_PREFIX",coeffssetid_prefix_);
312 : }
313 : } else {
314 : coeffssetid_prefix_="";
315 76 : if(keywords.exists("COEFFS_SET_ID_PREFIX")) {
316 152 : parse("COEFFS_SET_ID_PREFIX",coeffssetid_prefix_);
317 : }
318 76 : if(coeffssetid_prefix_.size()>0) {
319 0 : plumed_merror("COEFFS_SET_ID_PREFIX should only be given if optimizing multiple coefficient sets");
320 : }
321 : }
322 :
323 79 : if(keywords.exists("INITIAL_COEFFS")) {
324 : std::vector<std::string> initial_coeffs_fnames;
325 158 : parseFilenames("INITIAL_COEFFS",initial_coeffs_fnames);
326 79 : if(initial_coeffs_fnames.size()>0) {
327 1 : readCoeffsFromFiles(initial_coeffs_fnames,false);
328 1 : comm.Barrier();
329 1 : if(comm.Get_rank()==0 && use_mwalkers_mpi_) {
330 0 : multi_sim_comm.Barrier();
331 : }
332 1 : setAllCoeffsSetIterationCounters();
333 : }
334 79 : }
335 : //
336 :
337 : std::vector<std::string> coeffs_fnames;
338 79 : if(keywords.exists("COEFFS_FILE")) {
339 158 : parseFilenames("COEFFS_FILE",coeffs_fnames,"coeffs.data");
340 79 : bool start_opt_afresh=false;
341 79 : if(keywords.exists("START_OPTIMIZATION_AFRESH")) {
342 76 : parseFlag("START_OPTIMIZATION_AFRESH",start_opt_afresh);
343 76 : if(start_opt_afresh && !getRestart()) {
344 0 : plumed_merror("the START_OPTIMIZATION_AFRESH keyword should only be used when a restart has been triggered by the RESTART keyword or the MD code");
345 : }
346 : }
347 79 : if(getRestart()) {
348 28 : for(unsigned int i=0; i<coeffs_fnames.size(); i++) {
349 15 : IFile ifile;
350 15 : ifile.link(*this);
351 15 : if(use_mwalkers_mpi_) {
352 8 : ifile.enforceSuffix("");
353 : }
354 15 : bool file_exist = ifile.FileExist(coeffs_fnames[i]);
355 15 : if(!file_exist) {
356 0 : std::string fname = FileBase::appendSuffix(coeffs_fnames[i],ifile.getSuffix());
357 0 : plumed_merror("Cannot find coefficient file " + fname + " when trying to restart an optimzation. If you don't want to restart the optimzation please remove the RESTART keyword or use the RESTART=NO within the "+getName()+" action to locally disable the restart.");
358 : }
359 15 : }
360 13 : readCoeffsFromFiles(coeffs_fnames,true);
361 13 : comm.Barrier();
362 13 : if(comm.Get_rank()==0 && use_mwalkers_mpi_) {
363 4 : multi_sim_comm.Barrier();
364 : }
365 13 : unsigned int iter_opt_tmp = coeffs_pntrs_[0]->getIterationCounter();
366 15 : for(unsigned int i=1; i<ncoeffssets_; i++) {
367 2 : plumed_massert(coeffs_pntrs_[i]->getIterationCounter()==iter_opt_tmp,"the iteraton counter should be the same for all files when restarting from previous coefficient files\n");
368 : }
369 13 : if(start_opt_afresh) {
370 : setIterationCounter(0);
371 1 : log.printf(" Optimization started afresh at iteration %u\n",getIterationCounter());
372 2 : for(unsigned int i=0; i<ncoeffssets_; i++) {
373 1 : AuxCoeffs(i).setValues( Coeffs(i) );
374 : }
375 : } else {
376 : setIterationCounter(coeffs_pntrs_[0]->getIterationCounter());
377 12 : log.printf(" Optimization restarted at iteration %u\n",getIterationCounter());
378 : }
379 13 : setAllCoeffsSetIterationCounters();
380 : }
381 :
382 79 : std::string coeffs_wstride_tmpstr="";
383 158 : parse("COEFFS_OUTPUT",coeffs_wstride_tmpstr);
384 79 : if(coeffs_wstride_tmpstr!="OFF" && coeffs_wstride_tmpstr.size()>0) {
385 79 : Tools::convert(coeffs_wstride_tmpstr,coeffs_wstride_);
386 : }
387 79 : if(coeffs_wstride_tmpstr=="OFF") {
388 : coeffs_fnames.clear();
389 : }
390 79 : setupOFiles(coeffs_fnames,coeffsOFiles_,mw_single_files);
391 158 : parse("COEFFS_FMT",coeffs_output_fmt_);
392 79 : if(coeffs_output_fmt_.size()>0) {
393 162 : for(unsigned int i=0; i<ncoeffssets_; i++) {
394 84 : coeffs_pntrs_[i]->setOutputFmt(coeffs_output_fmt_);
395 : }
396 : }
397 79 : if(!getRestart()) {
398 136 : for(unsigned int i=0; i<coeffsOFiles_.size(); i++) {
399 70 : coeffs_pntrs_[i]->writeToFile(*coeffsOFiles_[i],aux_coeffs_pntrs_[i].get(),false);
400 : }
401 : }
402 79 : if(coeffs_fnames.size()>0) {
403 79 : if(ncoeffssets_==1) {
404 152 : log.printf(" Coefficients will be written out to file %s every %u iterations\n",coeffsOFiles_[0]->getPath().c_str(),coeffs_wstride_);
405 : } else {
406 3 : log.printf(" Coefficients will be written out to the following files every %u iterations:\n",coeffs_wstride_);
407 12 : for(unsigned int i=0; i<coeffs_fnames.size(); i++) {
408 18 : log.printf(" coefficient set %u: %s\n",i,coeffsOFiles_[i]->getPath().c_str());
409 : }
410 : }
411 : } else {
412 0 : log.printf(" Output of coefficients to file has been disabled\n");
413 : }
414 : }
415 :
416 : std::vector<std::string> gradient_fnames;
417 79 : if(keywords.exists("GRADIENT_FILE")) {
418 79 : parseFilenames("GRADIENT_FILE",gradient_fnames);
419 158 : parse("GRADIENT_OUTPUT",gradient_wstride_);
420 :
421 79 : if(coeffs_fnames.size()>0) {
422 151 : for(unsigned int i=0; i<gradient_fnames.size(); i++) {
423 72 : plumed_massert(gradient_fnames[i]!=coeffs_fnames[i],"COEFFS_FILE and GRADIENT_FILE cannot be the same");
424 : }
425 : }
426 79 : setupOFiles(gradient_fnames,gradientOFiles_,mw_single_files);
427 158 : parse("GRADIENT_FMT",gradient_output_fmt_);
428 79 : if(gradient_output_fmt_.size()>0) {
429 138 : for(unsigned int i=0; i<ncoeffssets_; i++) {
430 72 : gradient_pntrs_[i]->setOutputFmt(gradient_output_fmt_);
431 : }
432 : }
433 :
434 79 : if(gradient_fnames.size()>0) {
435 66 : if(ncoeffssets_==1) {
436 126 : log.printf(" Gradient will be written out to file %s every %u iterations\n",gradientOFiles_[0]->getPath().c_str(),gradient_wstride_);
437 : } else {
438 3 : log.printf(" Gradient will be written out to the following files every %u iterations:\n",gradient_wstride_);
439 12 : for(unsigned int i=0; i<gradient_fnames.size(); i++) {
440 18 : log.printf(" coefficient set %u: %s\n",i,gradientOFiles_[i]->getPath().c_str());
441 : }
442 : }
443 : }
444 : }
445 :
446 : std::vector<std::string> hessian_fnames;
447 79 : if(keywords.exists("HESSIAN_FILE")) {
448 76 : parseFilenames("HESSIAN_FILE",hessian_fnames);
449 152 : parse("HESSIAN_OUTPUT",hessian_wstride_);
450 :
451 76 : if(coeffs_fnames.size()>0) {
452 142 : for(unsigned int i=0; i<hessian_fnames.size(); i++) {
453 66 : plumed_massert(hessian_fnames[i]!=coeffs_fnames[i],"COEFFS_FILE and HESSIAN_FILE cannot be the same");
454 : }
455 : }
456 76 : if(gradient_fnames.size()>0) {
457 131 : for(unsigned int i=0; i<hessian_fnames.size(); i++) {
458 66 : plumed_massert(hessian_fnames[i]!=gradient_fnames[i],"GRADIENT_FILE and HESSIAN_FILE cannot be the same");
459 : }
460 : }
461 76 : setupOFiles(hessian_fnames,hessianOFiles_,mw_single_files);
462 152 : parse("HESSIAN_FMT",hessian_output_fmt_);
463 :
464 76 : if(hessian_fnames.size()>0) {
465 62 : if(ncoeffssets_==1) {
466 120 : log.printf(" Hessian will be written out to file %s every %u iterations\n",hessianOFiles_[0]->getPath().c_str(),hessian_wstride_);
467 : } else {
468 2 : log.printf(" Hessian will be written out to the following files every %u iterations:\n",hessian_wstride_);
469 8 : for(unsigned int i=0; i<hessian_fnames.size(); i++) {
470 12 : log.printf(" coefficient set %u: %s\n",i,hessianOFiles_[i]->getPath().c_str());
471 : }
472 : }
473 : }
474 : }
475 :
476 :
477 : //
478 79 : if(keywords.exists("MASK_FILE")) {
479 : std::vector<std::string> mask_fnames_in;
480 156 : parseVector("MASK_FILE",mask_fnames_in);
481 78 : if(mask_fnames_in.size()==1 && ncoeffssets_>1) {
482 0 : if(identical_coeffs_shape_) {
483 0 : mask_fnames_in.resize(ncoeffssets_,mask_fnames_in[0]);
484 : } else {
485 0 : plumed_merror("the coefficients indices shape differs between biases so you need to give a separate file for each coefficient set\n");
486 : }
487 : }
488 78 : if(mask_fnames_in.size()>0 && mask_fnames_in.size()!=ncoeffssets_) {
489 0 : plumed_merror("Error in MASK_FILE keyword: either give one value for all biases or a separate value for each coefficient set");
490 : }
491 :
492 78 : coeffs_mask_pntrs_.resize(ncoeffssets_);
493 162 : for(unsigned int i=0; i<ncoeffssets_; i++) {
494 168 : coeffs_mask_pntrs_[i] = Tools::make_unique<CoeffsVector>(*coeffs_pntrs_[i]);
495 168 : coeffs_mask_pntrs_[i]->setLabels("mask");
496 84 : coeffs_mask_pntrs_[i]->setValues(1.0);
497 168 : coeffs_mask_pntrs_[i]->setOutputFmt("%f");
498 : }
499 :
500 78 : if(mask_fnames_in.size()>0) {
501 1 : if(ncoeffssets_==1) {
502 1 : size_t nread = coeffs_mask_pntrs_[0]->readFromFile(mask_fnames_in[0],true,true);
503 1 : log.printf(" read %zu values from mask file %s\n",nread,mask_fnames_in[0].c_str());
504 1 : size_t ndeactived = coeffs_mask_pntrs_[0]->countValues(0.0);
505 1 : log.printf(" deactived optimization of %zu coefficients\n",ndeactived);
506 : } else {
507 0 : for(unsigned int i=0; i<ncoeffssets_; i++) {
508 0 : size_t nread = coeffs_mask_pntrs_[i]->readFromFile(mask_fnames_in[i],true,true);
509 0 : log.printf(" mask for coefficient set %u:\n",i);
510 0 : log.printf(" read %zu values from file %s\n",nread,mask_fnames_in[i].c_str());
511 0 : size_t ndeactived = coeffs_mask_pntrs_[0]->countValues(0.0);
512 0 : log.printf(" deactived optimization of %zu coefficients\n",ndeactived);
513 : }
514 : }
515 : }
516 :
517 : std::vector<std::string> mask_fnames_out;
518 78 : parseFilenames("OUTPUT_MASK_FILE",mask_fnames_out);
519 :
520 79 : for(unsigned int i=0; i<mask_fnames_out.size(); i++) {
521 1 : if(mask_fnames_in.size()>0) {
522 1 : plumed_massert(mask_fnames_out[i]!=mask_fnames_in[i],"MASK_FILE and OUTPUT_MASK_FILE cannot be the same");
523 : }
524 1 : OFile maskOFile;
525 1 : maskOFile.link(*this);
526 1 : maskOFile.enforceBackup();
527 1 : if(use_mwalkers_mpi_ && mwalkers_mpi_single_files_) {
528 0 : unsigned int r=0;
529 0 : if(comm.Get_rank()==0) {
530 0 : r=multi_sim_comm.Get_rank();
531 : }
532 0 : comm.Bcast(r,0);
533 0 : if(r>0) {
534 : mask_fnames_out[i]="/dev/null";
535 : }
536 0 : maskOFile.enforceSuffix("");
537 : }
538 1 : maskOFile.open(mask_fnames_out[i]);
539 1 : coeffs_mask_pntrs_[i]->writeToFile(maskOFile,true);
540 1 : maskOFile.close();
541 1 : }
542 78 : }
543 :
544 79 : if(getRestart() && ustride_targetdist_>0) {
545 16 : for(unsigned int i=0; i<nbiases_; i++) {
546 8 : if(dynamic_targetdists_[i]) {
547 8 : bias_pntrs_[i]->restartTargetDistributions();
548 : }
549 : }
550 : }
551 :
552 :
553 : std::vector<std::string> targetdist_averages_fnames;
554 79 : if(keywords.exists("TARGETDIST_AVERAGES_FILE")) {
555 158 : parseFilenames("TARGETDIST_AVERAGES_FILE",targetdist_averages_fnames,"targetdist-averages.data");
556 158 : parse("TARGETDIST_AVERAGES_OUTPUT",targetdist_averages_wstride_);
557 :
558 79 : if(coeffs_fnames.size()>0) {
559 164 : for(unsigned int i=0; i<targetdist_averages_fnames.size(); i++) {
560 85 : plumed_massert(targetdist_averages_fnames[i]!=coeffs_fnames[i],"COEFFS_FILE and TARGETDIST_AVERAGES_FILE cannot be the same");
561 : }
562 : }
563 79 : if(gradient_fnames.size()>0) {
564 138 : for(unsigned int i=0; i<targetdist_averages_fnames.size(); i++) {
565 72 : plumed_massert(targetdist_averages_fnames[i]!=gradient_fnames[i],"GRADIENT_FILE and TARGETDIST_AVERAGES_FILE cannot be the same");
566 : }
567 : }
568 79 : if(hessian_fnames.size()>0) {
569 128 : for(unsigned int i=0; i<targetdist_averages_fnames.size(); i++) {
570 66 : plumed_massert(targetdist_averages_fnames[i]!=hessian_fnames[i],"HESSIAN_FILE and TARGETDIST_AVERAGES_FILE cannot be the same");
571 : }
572 : }
573 79 : setupOFiles(targetdist_averages_fnames,targetdist_averagesOFiles_,mw_single_files);
574 158 : parse("TARGETDIST_AVERAGES_FMT",targetdist_averages_output_fmt_);
575 79 : if(targetdist_averages_output_fmt_.size()>0) {
576 156 : for(unsigned int i=0; i<ncoeffssets_; i++) {
577 81 : targetdist_averages_pntrs_[i]->setOutputFmt(targetdist_averages_output_fmt_);
578 : }
579 : }
580 :
581 164 : for(unsigned int i=0; i<targetdist_averagesOFiles_.size(); i++) {
582 85 : targetdist_averages_pntrs_[i]->writeToFile(*targetdist_averagesOFiles_[i]);
583 : }
584 :
585 79 : if(targetdist_averages_wstride_==0) {
586 : targetdist_averagesOFiles_.clear();
587 : }
588 :
589 79 : if(targetdist_averages_fnames.size()>0 && targetdist_averages_wstride_ > 0) {
590 34 : if(ncoeffssets_==1) {
591 68 : log.printf(" Target distribution averages will be written out to file %s every %u iterations\n",targetdist_averagesOFiles_[0]->getPath().c_str(),targetdist_averages_wstride_);
592 : } else {
593 0 : log.printf(" Target distribution averages will be written out to the following files every %u iterations:\n",targetdist_averages_wstride_);
594 0 : for(unsigned int i=0; i<targetdist_averages_fnames.size(); i++) {
595 0 : log.printf(" coefficient set %u: %s\n",i,targetdist_averagesOFiles_[i]->getPath().c_str());
596 : }
597 : }
598 : }
599 : }
600 :
601 :
602 79 : if(keywords.exists("BIAS_OUTPUT")) {
603 79 : parse("BIAS_OUTPUT",bias_output_stride_);
604 79 : if(bias_output_stride_>0) {
605 74 : bias_output_active_=true;
606 154 : for(unsigned int i=0; i<nbiases_; i++) {
607 80 : bias_pntrs_[i]->enableBiasFileOutput();
608 80 : bias_pntrs_[i]->setupBiasFileOutput();
609 80 : bias_pntrs_[i]->writeBiasToFile();
610 : }
611 : } else {
612 5 : bias_output_active_=false;
613 5 : bias_output_stride_=1000;
614 : }
615 : }
616 :
617 79 : if(keywords.exists("FES_OUTPUT")) {
618 79 : parse("FES_OUTPUT",fes_output_stride_);
619 79 : if(fes_output_stride_>0) {
620 74 : fes_output_active_=true;
621 154 : for(unsigned int i=0; i<nbiases_; i++) {
622 80 : bias_pntrs_[i]->enableFesFileOutput();
623 80 : bias_pntrs_[i]->setupFesFileOutput();
624 80 : bias_pntrs_[i]->writeFesToFile();
625 : }
626 : } else {
627 5 : fes_output_active_=false;
628 5 : fes_output_stride_=1000;
629 : }
630 : }
631 :
632 79 : if(keywords.exists("FES_PROJ_OUTPUT")) {
633 79 : parse("FES_PROJ_OUTPUT",fesproj_output_stride_);
634 79 : if(fesproj_output_stride_>0) {
635 16 : fesproj_output_active_=true;
636 32 : for(unsigned int i=0; i<nbiases_; i++) {
637 16 : bias_pntrs_[i]->enableFesProjFileOutput();
638 16 : bias_pntrs_[i]->setupFesProjFileOutput();
639 16 : bias_pntrs_[i]->writeFesProjToFile();
640 : }
641 : } else {
642 63 : fesproj_output_active_=false;
643 63 : fesproj_output_stride_=1000;
644 : }
645 : }
646 :
647 164 : for(unsigned int i=0; i<nbiases_; i++) {
648 85 : if(!dynamic_targetdists_[i] && bias_pntrs_[i]->isStaticTargetDistFileOutputActive()) {
649 4 : bias_pntrs_[i]->setupTargetDistFileOutput();
650 4 : bias_pntrs_[i]->writeTargetDistToFile();
651 4 : bias_pntrs_[i]->setupTargetDistProjFileOutput();
652 4 : bias_pntrs_[i]->writeTargetDistProjToFile();
653 : }
654 : }
655 :
656 79 : if(keywords.exists("TARGETDIST_OUTPUT")) {
657 78 : parse("TARGETDIST_OUTPUT",targetdist_output_stride_);
658 78 : if(targetdist_output_stride_>0) {
659 37 : if(ustride_targetdist_==0) {
660 0 : plumed_merror("it doesn't make sense to use the TARGETDIST_OUTPUT keyword if you don't have a target distribution that needs to be updated");
661 : }
662 37 : if(targetdist_output_stride_%ustride_targetdist_!=0) {
663 0 : plumed_merror("the value given in TARGETDIST_OUTPUT doesn't make sense, it should be multiple of TARGETDIST_STRIDE");
664 : }
665 37 : if(targetdist_output_stride_%coeffs_wstride_!=0) {
666 0 : plumed_merror("the value given in TARGETDIST_OUTPUT doesn't make sense, it should be multiple of COEFFS_OUTPUT");
667 : }
668 :
669 37 : targetdist_output_active_=true;
670 74 : for(unsigned int i=0; i<nbiases_; i++) {
671 37 : if(dynamic_targetdists_[i]) {
672 37 : bias_pntrs_[i]->enableDynamicTargetDistFileOutput();
673 37 : bias_pntrs_[i]->setupTargetDistFileOutput();
674 37 : bias_pntrs_[i]->writeTargetDistToFile();
675 : }
676 : }
677 : } else {
678 41 : targetdist_output_active_=false;
679 41 : targetdist_output_stride_=1000;
680 : }
681 : }
682 :
683 79 : if(keywords.exists("TARGETDIST_PROJ_OUTPUT")) {
684 78 : parse("TARGETDIST_PROJ_OUTPUT",targetdist_proj_output_stride_);
685 78 : if(targetdist_proj_output_stride_>0) {
686 3 : if(ustride_targetdist_==0) {
687 0 : plumed_merror("it doesn't make sense to use the TARGETDIST_PROJ_OUTPUT keyword if you don't have a target distribution that needs to be updated");
688 : }
689 3 : if(targetdist_proj_output_stride_%ustride_targetdist_!=0) {
690 0 : plumed_merror("the value given in TARGETDIST_PROJ_OUTPUT doesn't make sense, it should be multiple of TARGETDIST_STRIDE");
691 : }
692 :
693 3 : targetdist_proj_output_active_=true;
694 6 : for(unsigned int i=0; i<nbiases_; i++) {
695 3 : if(dynamic_targetdists_[i]) {
696 3 : bias_pntrs_[i]->enableDynamicTargetDistFileOutput();
697 3 : bias_pntrs_[i]->setupTargetDistProjFileOutput();
698 3 : bias_pntrs_[i]->writeTargetDistProjToFile();
699 : }
700 : }
701 : } else {
702 75 : targetdist_proj_output_active_=false;
703 75 : targetdist_proj_output_stride_=1000;
704 : }
705 : }
706 :
707 79 : if(ncoeffssets_==1) {
708 76 : log.printf(" Output Components:\n");
709 76 : log.printf(" ");
710 76 : if(monitor_instantaneous_gradient_) {
711 6 : addComponent("gradrms");
712 3 : componentIsNotPeriodic("gradrms");
713 3 : log.printf(" ");
714 6 : addComponent("gradmax");
715 6 : componentIsNotPeriodic("gradmax");
716 : }
717 76 : if(aver_gradient_pntrs_.size()>0) {
718 2 : log.printf(" ");
719 4 : addComponent("avergradrms");
720 2 : componentIsNotPeriodic("avergradrms");
721 2 : log.printf(" ");
722 4 : addComponent("avergradmax");
723 4 : componentIsNotPeriodic("avergradmax");
724 : }
725 76 : if(!fixed_stepsize_) {
726 1 : log.printf(" ");
727 2 : addComponent("stepsize");
728 1 : componentIsNotPeriodic("stepsize");
729 2 : getPntrToComponent("stepsize")->set( getCurrentStepSize(0) );
730 : }
731 : } else {
732 12 : for(unsigned int i=0; i<ncoeffssets_; i++) {
733 9 : log.printf(" Output Components for coefficient set %u:\n",i);
734 9 : std::string is="";
735 9 : Tools::convert(i,is);
736 18 : is = "-" + coeffssetid_prefix_ + is;
737 9 : log.printf(" ");
738 9 : if(monitor_instantaneous_gradient_) {
739 6 : addComponent("gradrms"+is);
740 3 : componentIsNotPeriodic("gradrms"+is);
741 3 : log.printf(" ");
742 6 : addComponent("gradmax"+is);
743 6 : componentIsNotPeriodic("gradmax"+is);
744 : }
745 9 : if(aver_gradient_pntrs_.size()>0) {
746 0 : log.printf(" ");
747 0 : addComponent("avergradrms"+is);
748 0 : componentIsNotPeriodic("avergradrms"+is);
749 0 : log.printf(" ");
750 0 : addComponent("avergradmax"+is);
751 0 : componentIsNotPeriodic("avergradmax"+is);
752 : }
753 9 : if(!fixed_stepsize_) {
754 0 : log.printf(" ");
755 0 : addComponent("stepsize"+is);
756 0 : componentIsNotPeriodic("stepsize"+is);
757 0 : getPntrToComponent("stepsize"+is)->set( getCurrentStepSize(i) );
758 : }
759 : }
760 : }
761 :
762 158 : }
763 :
764 :
765 79 : Optimizer::~Optimizer() {
766 : //
767 164 : for(unsigned int i=0; i<ncoeffssets_; i++) {
768 85 : if(coeffsOFiles_.size()>0 && getIterationCounter()%coeffs_wstride_!=0) {
769 2 : coeffs_pntrs_[i]->writeToFile(*coeffsOFiles_[i],aux_coeffs_pntrs_[i].get(),false);
770 : }
771 85 : if(targetdist_averagesOFiles_.size()>0 && iter_counter%targetdist_averages_wstride_!=0) {
772 1 : targetdist_averages_pntrs_[i]->writeToFile(*targetdist_averagesOFiles_[i]);
773 : }
774 : }
775 : //
776 79 : if(!isTargetDistOutputActive()) {
777 90 : for(unsigned int i=0; i<nbiases_; i++) {
778 48 : if(dynamic_targetdists_[i]) {
779 1 : bias_pntrs_[i]->enableDynamicTargetDistFileOutput();
780 1 : bias_pntrs_[i]->setupTargetDistFileOutput();
781 1 : bias_pntrs_[i]->writeTargetDistToFile();
782 : }
783 : }
784 : }
785 : //
786 79 : if(isBiasOutputActive() && getIterationCounter()%getBiasOutputStride()!=0) {
787 1 : writeBiasOutputFiles();
788 : }
789 79 : if(isFesOutputActive() && getIterationCounter()%getFesOutputStride()!=0) {
790 1 : writeFesOutputFiles();
791 : }
792 79 : if(isFesProjOutputActive() && getIterationCounter()%getFesProjOutputStride()!=0) {
793 1 : writeFesProjOutputFiles();
794 : }
795 79 : if(isTargetDistOutputActive() && getIterationCounter()%getTargetDistOutputStride()!=0) {
796 2 : writeTargetDistOutputFiles();
797 : }
798 79 : if(isTargetDistProjOutputActive() && getIterationCounter()%getTargetDistProjOutputStride()!=0) {
799 1 : writeTargetDistProjOutputFiles();
800 : }
801 395 : }
802 :
803 :
804 87 : void Optimizer::registerKeywords( Keywords& keys ) {
805 87 : Action::registerKeywords(keys);
806 87 : ActionPilot::registerKeywords(keys);
807 87 : ActionWithValue::registerKeywords(keys);
808 : //
809 174 : keys.remove("NUMERICAL_DERIVATIVES");
810 : // Default always active keywords
811 87 : keys.add("compulsory","BIAS","the label of the VES bias to be optimized");
812 87 : keys.add("compulsory","STRIDE","the frequency of updating the coefficients given in the number of MD steps.");
813 87 : keys.add("compulsory","COEFFS_FILE","coeffs.data","the name of output file for the coefficients");
814 87 : keys.add("compulsory","COEFFS_OUTPUT","100","how often the coefficients should be written to file. This parameter is given as the number of iterations.");
815 87 : keys.add("optional","COEFFS_FMT","specify format for coefficient file(s) (useful for decrease the number of digits in regtests)");
816 87 : keys.add("optional","COEFFS_SET_ID_PREFIX","suffix to add to the filename given in FILE to identify the bias, should only be given if a single filename is given in FILE when optimizing multiple biases.");
817 : //
818 87 : keys.add("optional","INITIAL_COEFFS","the name(s) of file(s) with the initial coefficients");
819 : // Hidden keywords to output the gradient to a file.
820 87 : keys.add("hidden","GRADIENT_FILE","the name of output file for the gradient");
821 87 : keys.add("hidden","GRADIENT_OUTPUT","how often the gradient should be written to file. This parameter is given as the number of bias iterations. It is by default 100 if GRADIENT_FILE is specficed");
822 87 : keys.add("hidden","GRADIENT_FMT","specify format for gradient file(s) (useful for decrease the number of digits in regtests)");
823 : // Either use a fixed stepsize (useFixedStepSizeKeywords) or changing stepsize (useDynamicsStepSizeKeywords)
824 87 : keys.reserve("compulsory","STEPSIZE","the step size used for the optimization");
825 87 : keys.reserve("compulsory","INITIAL_STEPSIZE","the initial step size used for the optimization");
826 : // Keywords related to the Hessian, actived with the useHessianKeywords function
827 174 : keys.reserveFlag("FULL_HESSIAN",false,"if the full Hessian matrix should be used for the optimization, otherwise only the diagonal part of the Hessian is used");
828 87 : keys.reserve("hidden","HESSIAN_FILE","the name of output file for the Hessian");
829 87 : keys.reserve("hidden","HESSIAN_OUTPUT","how often the Hessian should be written to file. This parameter is given as the number of bias iterations. It is by default 100 if HESSIAN_FILE is specficed");
830 87 : keys.reserve("hidden","HESSIAN_FMT","specify format for hessian file(s) (useful for decrease the number of digits in regtests)");
831 : // Keywords related to the multiple walkers, actived with the useMultipleWalkersKeywords function
832 174 : keys.reserveFlag("MULTIPLE_WALKERS",false,"if optimization is to be performed using multiple walkers connected via MPI");
833 : // Keywords related to the mask file, actived with the useMaskKeywords function
834 87 : keys.reserve("optional","MASK_FILE","read in a mask file which allows one to employ different step sizes for different coefficients and/or deactivate the optimization of certain coefficients (by putting values of 0.0). One can write out the resulting mask by using the OUTPUT_MASK_FILE keyword.");
835 87 : keys.reserve("optional","OUTPUT_MASK_FILE","Name of the file to write out the mask resulting from using the MASK_FILE keyword. Can also be used to generate a template mask file.");
836 : //
837 174 : keys.reserveFlag("START_OPTIMIZATION_AFRESH",false,"if the iterations should be started afresh when a restart has been triggered by the RESTART keyword or the MD code.");
838 : //
839 87 : keys.addFlag("MONITOR_INSTANTANEOUS_GRADIENT",false,"if quantities related to the instantaneous gradient should be outputted.");
840 : //
841 174 : keys.reserveFlag("MONITOR_AVERAGE_GRADIENT",false,"if the averaged gradient should be monitored and quantities related to it should be outputted.");
842 87 : keys.reserve("optional","MONITOR_AVERAGES_GRADIENT_EXP_DECAY","use an exponentially decaying averaging with a given time constant when monitoring the averaged gradient");
843 : //
844 87 : keys.reserve("optional","TARGETDIST_STRIDE","stride for updating a target distribution that is iteratively updated during the optimization. Note that the value is given in terms of coefficient iterations.");
845 87 : keys.reserve("optional","TARGETDIST_OUTPUT","how often the dynamic target distribution(s) should be written out to file. Note that the value is given in terms of coefficient iterations.");
846 87 : keys.reserve("optional","TARGETDIST_PROJ_OUTPUT","how often the projections of the dynamic target distribution(s) should be written out to file. Note that the value is given in terms of coefficient iterations.");
847 : //
848 87 : keys.add("optional","TARGETDIST_AVERAGES_FILE","the name of output file for the target distribution averages. By default it is targetdist-averages.data.");
849 87 : keys.add("optional","TARGETDIST_AVERAGES_OUTPUT","how often the target distribution averages should be written out to file. Note that the value is given in terms of coefficient iterations. If no value is given are the averages only written at the beginning of the optimization");
850 87 : keys.add("hidden","TARGETDIST_AVERAGES_FMT","specify format for target distribution averages file(s) (useful for decrease the number of digits in regtests)");
851 : //
852 87 : keys.add("optional","BIAS_OUTPUT","how often the bias(es) should be written out to file. Note that the value is given in terms of coefficient iterations.");
853 87 : keys.add("optional","FES_OUTPUT","how often the FES(s) should be written out to file. Note that the value is given in terms of coefficient iterations.");
854 87 : keys.add("optional","FES_PROJ_OUTPUT","how often the projections of the FES(s) should be written out to file. Note that the value is given in terms of coefficient iterations.");
855 : //
856 87 : keys.reserve("optional","REWEIGHT_FACTOR_STRIDE","stride for updating the reweighting factor c(t). Note that the value is given in terms of coefficient iterations.");
857 : //
858 87 : keys.use("RESTART");
859 : //
860 87 : keys.use("UPDATE_FROM");
861 87 : keys.use("UPDATE_UNTIL");
862 : // Components that are always active
863 174 : keys.addOutputComponent("gradrms","MONITOR_INSTANTANEOUS_GRADIENT","scalar","the root mean square value of the coefficient gradient. For multiple biases this component is labeled using the number of the bias as gradrms-#.");
864 174 : keys.addOutputComponent("gradmax","MONITOR_INSTANTANEOUS_GRADIENT","scalar","the largest absolute value of the coefficient gradient. For multiple biases this component is labeled using the number of the bias as gradmax-#.");
865 : // keys.addOutputComponent("gradmaxidx","default","the index of the maximum absolute value of the gradient");
866 174 : keys.setValueDescription("scalar","a scalar");
867 87 : }
868 :
869 :
870 80 : void Optimizer::useHessianKeywords(Keywords& keys) {
871 : // keys.use("FULL_HESSIAN");
872 80 : keys.use("HESSIAN_FILE");
873 80 : keys.use("HESSIAN_OUTPUT");
874 80 : keys.use("HESSIAN_FMT");
875 80 : }
876 :
877 :
878 87 : void Optimizer::useMultipleWalkersKeywords(Keywords& keys) {
879 87 : keys.use("MULTIPLE_WALKERS");
880 87 : }
881 :
882 :
883 81 : void Optimizer::useFixedStepSizeKeywords(Keywords& keys) {
884 81 : keys.use("STEPSIZE");
885 81 : }
886 :
887 :
888 3 : void Optimizer::useDynamicStepSizeKeywords(Keywords& keys) {
889 3 : keys.use("INITIAL_STEPSIZE");
890 6 : keys.addOutputComponent("stepsize","default","scalar","the current value of step size used to update the coefficients. For multiple biases this component is labeled using the number of the bias as stepsize-#.");
891 3 : }
892 :
893 :
894 84 : void Optimizer::useMaskKeywords(Keywords& keys) {
895 84 : keys.use("MASK_FILE");
896 84 : keys.use("OUTPUT_MASK_FILE");
897 84 : }
898 :
899 :
900 80 : void Optimizer::useRestartKeywords(Keywords& keys) {
901 80 : keys.use("START_OPTIMIZATION_AFRESH");
902 80 : }
903 :
904 :
905 80 : void Optimizer::useMonitorAverageGradientKeywords(Keywords& keys) {
906 80 : keys.use("MONITOR_AVERAGE_GRADIENT");
907 80 : keys.use("MONITOR_AVERAGES_GRADIENT_EXP_DECAY");
908 160 : keys.addOutputComponent("avergradrms","MONITOR_AVERAGE_GRADIENT","scalar","the root mean square value of the averaged coefficient gradient. For multiple biases this component is labeled using the number of the bias as gradrms-#.");
909 160 : keys.addOutputComponent("avergradmax","MONITOR_AVERAGE_GRADIENT","scalar","the largest absolute value of the averaged coefficient gradient. For multiple biases this component is labeled using the number of the bias as gradmax-#.");
910 80 : }
911 :
912 :
913 84 : void Optimizer::useDynamicTargetDistributionKeywords(Keywords& keys) {
914 84 : keys.use("TARGETDIST_STRIDE");
915 84 : keys.use("TARGETDIST_OUTPUT");
916 84 : keys.use("TARGETDIST_PROJ_OUTPUT");
917 84 : }
918 :
919 :
920 0 : void Optimizer::useReweightFactorKeywords(Keywords& keys) {
921 0 : keys.use("REWEIGHT_FACTOR_STRIDE");
922 0 : }
923 :
924 :
925 76 : void Optimizer::turnOnHessian() {
926 76 : plumed_massert(hessian_pntrs_.size()==0,"turnOnHessian() should only be run during initialization");
927 76 : use_hessian_=true;
928 76 : hessian_pntrs_.clear();
929 158 : for(unsigned int i=0; i<nbiases_; i++) {
930 82 : std::vector<CoeffsMatrix*> pntrs_hessian = enableHessian(bias_pntrs_[i],diagonal_hessian_);
931 164 : for(unsigned int k=0; k<pntrs_hessian.size(); k++) {
932 82 : pntrs_hessian[k]->turnOnIterationCounter();
933 82 : pntrs_hessian[k]->setIterationCounterAndTime(getIterationCounter(),getTime());
934 82 : hessian_pntrs_.push_back(pntrs_hessian[k]);
935 : }
936 : }
937 76 : plumed_massert(hessian_pntrs_.size()==ncoeffssets_,"problems in linking Hessians");
938 76 : if(diagonal_hessian_) {
939 76 : log.printf(" Optimization performed using diagonal Hessian matrix\n");
940 : } else {
941 0 : log.printf(" Optimization performed using full Hessian matrix\n");
942 : }
943 : //
944 76 : if(hessian_output_fmt_.size()>0) {
945 128 : for(unsigned int i=0; i<ncoeffssets_; i++) {
946 66 : hessian_pntrs_[i]->setOutputFmt(hessian_output_fmt_);
947 : }
948 : }
949 :
950 76 : }
951 :
952 :
953 0 : void Optimizer::turnOffHessian() {
954 0 : use_hessian_=false;
955 0 : for(unsigned int i=0; i<nbiases_; i++) {
956 0 : bias_pntrs_[i]->disableHessian();
957 : }
958 : hessian_pntrs_.clear();
959 0 : hessianOFiles_.clear();
960 0 : }
961 :
962 :
963 82 : std::vector<CoeffsMatrix*> Optimizer::enableHessian(VesBias* bias_pntr_in, const bool diagonal_hessian) {
964 82 : plumed_massert(use_hessian_,"the Hessian should not be used");
965 82 : bias_pntr_in->enableHessian(diagonal_hessian);
966 : std::vector<CoeffsMatrix*> hessian_pntrs_out = bias_pntr_in->getHessianPntrs();
967 164 : for(unsigned int k=0; k<hessian_pntrs_out.size(); k++) {
968 82 : plumed_massert(hessian_pntrs_out[k] != NULL,"Hessian is needed but not linked correctly");
969 : }
970 82 : return hessian_pntrs_out;
971 : }
972 :
973 :
974 : // CoeffsMatrix* Optimizer::switchToDiagonalHessian(VesBias* bias_pntr_in) {
975 : // plumed_massert(use_hessian_,"it does not make sense to switch to diagonal Hessian if it Hessian is not used");
976 : // diagonal_hessian_=true;
977 : // bias_pntr_in->enableHessian(diagonal_hessian_);
978 : // CoeffsMatrix* hessian_pntr_out = bias_pntr_in->getHessianPntr();
979 : // plumed_massert(hessian_pntr_out != NULL,"Hessian is needed but not linked correctly");
980 : // //
981 : // log.printf(" %s (with label %s): switching to a diagonal Hessian for VES bias %s (with label %s) at time %f\n",getName().c_str(),getLabel().c_str(),bias_pntr_in->getName().c_str(),bias_pntr_in->getLabel().c_str(),getTime());
982 : // return hessian_pntr_out;
983 : // }
984 :
985 :
986 : // CoeffsMatrix* Optimizer::switchToFullHessian(VesBias* bias_pntr_in) {
987 : // plumed_massert(use_hessian_,"it does not make sense to switch to diagonal Hessian if it Hessian is not used");
988 : // diagonal_hessian_=false;
989 : // bias_pntr_in->enableHessian(diagonal_hessian_);
990 : // CoeffsMatrix* hessian_pntr_out = bias_pntr_in->getHessianPntr();
991 : // plumed_massert(hessian_pntr_out != NULL,"Hessian is needed but not linked correctly");
992 : // //
993 : // log.printf(" %s (with label %s): switching to a diagonal Hessian for VES bias %s (with label %s) at time %f\n",getName().c_str(),getLabel().c_str(),bias_pntr_in->getName().c_str(),bias_pntr_in->getLabel().c_str(),getTime());
994 : // return hessian_pntr_out;
995 : // }
996 :
997 :
998 22879 : void Optimizer::update() {
999 22879 : if(onStep() && !isFirstStep) {
1000 45560 : for(unsigned int i=0; i<nbiases_; i++) {
1001 22810 : bias_pntrs_[i]->updateGradientAndHessian(use_mwalkers_mpi_);
1002 : }
1003 45560 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1004 22810 : if(gradient_pntrs_[i]->isActive()) {
1005 22715 : coeffsUpdate(i);
1006 : } else {
1007 190 : std::string msg = "iteration " + getIterationCounterStr(+1) +
1008 285 : " for " + bias_pntrs_[i]->getLabel() +
1009 95 : " - the coefficients are not updated as CV values are outside the bias intervals";
1010 95 : warning(msg);
1011 : }
1012 :
1013 : // +1 as this is done before increaseIterationCounter() is used
1014 22810 : unsigned int curr_iter = getIterationCounter()+1;
1015 22810 : double curr_time = getTime();
1016 22810 : coeffs_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1017 : aux_coeffs_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1018 22810 : gradient_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1019 22810 : targetdist_averages_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1020 22810 : if(use_hessian_) {
1021 22780 : hessian_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1022 : }
1023 22810 : if(aver_gradient_pntrs_.size()>0) {
1024 : aver_gradient_pntrs_[i]->setIterationCounterAndTime(curr_iter,curr_time);
1025 20 : aver_gradient_pntrs_[i]->addToAverage(*gradient_pntrs_[i]);
1026 : }
1027 : }
1028 : increaseIterationCounter();
1029 :
1030 22750 : if(ustride_targetdist_>0 && getIterationCounter()%ustride_targetdist_==0) {
1031 708 : for(unsigned int i=0; i<nbiases_; i++) {
1032 354 : if(dynamic_targetdists_[i]) {
1033 354 : bias_pntrs_[i]->updateTargetDistributions();
1034 : }
1035 : }
1036 : }
1037 22750 : if(ustride_reweightfactor_>0 && getIterationCounter()%ustride_reweightfactor_==0) {
1038 0 : for(unsigned int i=0; i<nbiases_; i++) {
1039 0 : bias_pntrs_[i]->updateReweightFactor();
1040 : }
1041 : }
1042 :
1043 22750 : updateOutputComponents();
1044 45560 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1045 22810 : writeOutputFiles(i);
1046 : }
1047 : //
1048 22750 : if(isBiasOutputActive() && getIterationCounter()%getBiasOutputStride()==0) {
1049 74 : writeBiasOutputFiles();
1050 : }
1051 22750 : if(isFesOutputActive() && getIterationCounter()%getFesOutputStride()==0) {
1052 74 : writeFesOutputFiles();
1053 : }
1054 22750 : if(isFesProjOutputActive() && getIterationCounter()%getFesProjOutputStride()==0) {
1055 16 : writeFesProjOutputFiles();
1056 : }
1057 22750 : if(isTargetDistOutputActive() && getIterationCounter()%getTargetDistOutputStride()==0) {
1058 36 : writeTargetDistOutputFiles();
1059 : }
1060 22750 : if(isTargetDistProjOutputActive() && getIterationCounter()%getTargetDistProjOutputStride()==0) {
1061 3 : writeTargetDistProjOutputFiles();
1062 : }
1063 : } else {
1064 129 : isFirstStep=false;
1065 : }
1066 22879 : }
1067 :
1068 :
1069 22750 : void Optimizer::updateOutputComponents() {
1070 22750 : if(ncoeffssets_==1) {
1071 22720 : if(!fixed_stepsize_) {
1072 20 : getPntrToComponent("stepsize")->set( getCurrentStepSize(0) );
1073 : }
1074 22720 : if(monitor_instantaneous_gradient_) {
1075 30 : getPntrToComponent("gradrms")->set( gradient_pntrs_[0]->getRMS() );
1076 30 : size_t gradient_maxabs_idx=0;
1077 60 : getPntrToComponent("gradmax")->set( gradient_pntrs_[0]->getMaxAbsValue(gradient_maxabs_idx) );
1078 : }
1079 22720 : if(aver_gradient_pntrs_.size()>0) {
1080 20 : getPntrToComponent("avergradrms")->set( aver_gradient_pntrs_[0]->getRMS() );
1081 20 : size_t avergradient_maxabs_idx=0;
1082 40 : getPntrToComponent("avergradmax")->set( aver_gradient_pntrs_[0]->getMaxAbsValue(avergradient_maxabs_idx) );
1083 : }
1084 : } else {
1085 120 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1086 90 : std::string is="";
1087 90 : Tools::convert(i,is);
1088 180 : is = "-" + coeffssetid_prefix_ + is;
1089 90 : if(!fixed_stepsize_) {
1090 0 : getPntrToComponent("stepsize"+is)->set( getCurrentStepSize(i) );
1091 : }
1092 90 : if(monitor_instantaneous_gradient_) {
1093 30 : getPntrToComponent("gradrms"+is)->set( gradient_pntrs_[i]->getRMS() );
1094 30 : size_t gradient_maxabs_idx=0;
1095 60 : getPntrToComponent("gradmax"+is)->set( gradient_pntrs_[i]->getMaxAbsValue(gradient_maxabs_idx) );
1096 : }
1097 90 : if(aver_gradient_pntrs_.size()>0) {
1098 0 : getPntrToComponent("avergradrms"+is)->set( aver_gradient_pntrs_[0]->getRMS() );
1099 0 : size_t avergradient_maxabs_idx=0;
1100 0 : getPntrToComponent("avergradmax"+is)->set( aver_gradient_pntrs_[0]->getMaxAbsValue(avergradient_maxabs_idx) );
1101 : }
1102 : }
1103 : }
1104 22750 : }
1105 :
1106 :
1107 1 : void Optimizer::turnOffCoeffsOutputFiles() {
1108 1 : coeffsOFiles_.clear();
1109 1 : }
1110 :
1111 :
1112 22810 : void Optimizer::writeOutputFiles(const unsigned int coeffs_id) {
1113 22810 : if(coeffsOFiles_.size()>0 && iter_counter%coeffs_wstride_==0) {
1114 789 : coeffs_pntrs_[coeffs_id]->writeToFile(*coeffsOFiles_[coeffs_id],aux_coeffs_pntrs_[coeffs_id].get(),false);
1115 : }
1116 22810 : if(gradientOFiles_.size()>0 && iter_counter%gradient_wstride_==0) {
1117 720 : if(aver_gradient_pntrs_.size()==0) {
1118 700 : gradient_pntrs_[coeffs_id]->writeToFile(*gradientOFiles_[coeffs_id],false);
1119 : } else {
1120 20 : gradient_pntrs_[coeffs_id]->writeToFile(*gradientOFiles_[coeffs_id],aver_gradient_pntrs_[coeffs_id].get(),false);
1121 : }
1122 : }
1123 22810 : if(hessianOFiles_.size()>0 && iter_counter%hessian_wstride_==0) {
1124 660 : hessian_pntrs_[coeffs_id]->writeToFile(*hessianOFiles_[coeffs_id]);
1125 : }
1126 22810 : if(targetdist_averagesOFiles_.size()>0 && iter_counter%targetdist_averages_wstride_==0) {
1127 333 : targetdist_averages_pntrs_[coeffs_id]->writeToFile(*targetdist_averagesOFiles_[coeffs_id]);
1128 : }
1129 22810 : }
1130 :
1131 :
1132 388 : void Optimizer::setupOFiles(std::vector<std::string>& fnames, std::vector<std::unique_ptr<OFile>>& OFiles, const bool multi_sim_single_files) {
1133 388 : plumed_assert(ncoeffssets_>0);
1134 388 : OFiles.resize(fnames.size());
1135 702 : for(unsigned int i=0; i<fnames.size(); i++) {
1136 314 : OFiles[i] = Tools::make_unique<OFile>();
1137 314 : OFiles[i]->link(*this);
1138 314 : if(multi_sim_single_files) {
1139 48 : unsigned int r=0;
1140 48 : if(comm.Get_rank()==0) {
1141 48 : r=multi_sim_comm.Get_rank();
1142 : }
1143 48 : comm.Bcast(r,0);
1144 48 : if(r>0) {
1145 : fnames[i]="/dev/null";
1146 : }
1147 96 : OFiles[i]->enforceSuffix("");
1148 : }
1149 314 : OFiles[i]->open(fnames[i]);
1150 314 : OFiles[i]->setHeavyFlush();
1151 : }
1152 388 : }
1153 :
1154 :
1155 14 : void Optimizer::readCoeffsFromFiles(const std::vector<std::string>& fnames, const bool read_aux_coeffs) {
1156 14 : plumed_assert(ncoeffssets_>0);
1157 14 : plumed_assert(fnames.size()==ncoeffssets_);
1158 14 : if(ncoeffssets_==1) {
1159 13 : log.printf(" Read in coefficients from file ");
1160 : } else {
1161 1 : log.printf(" Read in coefficients from files:\n");
1162 : }
1163 30 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1164 16 : IFile ifile;
1165 16 : ifile.link(*this);
1166 16 : if(use_mwalkers_mpi_ && mwalkers_mpi_single_files_) {
1167 8 : ifile.enforceSuffix("");
1168 : }
1169 16 : ifile.open(fnames[i]);
1170 32 : if(!ifile.FieldExist(coeffs_pntrs_[i]->getDataLabel())) {
1171 0 : std::string error_msg = "Problem with reading coefficients from file " + ifile.getPath() + ": no field with name " + coeffs_pntrs_[i]->getDataLabel() + "\n";
1172 0 : plumed_merror(error_msg);
1173 : }
1174 16 : size_t ncoeffs_read = coeffs_pntrs_[i]->readFromFile(ifile,false,false);
1175 16 : if(ncoeffssets_==1) {
1176 26 : log.printf("%s (read %zu of %zu values)\n", ifile.getPath().c_str(),ncoeffs_read,coeffs_pntrs_[i]->numberOfCoeffs());
1177 : } else {
1178 6 : log.printf(" coefficient set %u: %s (read %zu of %zu values)\n",i,ifile.getPath().c_str(),ncoeffs_read,coeffs_pntrs_[i]->numberOfCoeffs());
1179 : }
1180 16 : ifile.close();
1181 16 : if(read_aux_coeffs) {
1182 15 : ifile.open(fnames[i]);
1183 30 : if(!ifile.FieldExist(aux_coeffs_pntrs_[i]->getDataLabel())) {
1184 0 : std::string error_msg = "Problem with reading coefficients from file " + ifile.getPath() + ": no field with name " + aux_coeffs_pntrs_[i]->getDataLabel() + "\n";
1185 0 : plumed_merror(error_msg);
1186 : }
1187 15 : aux_coeffs_pntrs_[i]->readFromFile(ifile,false,false);
1188 15 : ifile.close();
1189 : } else {
1190 1 : AuxCoeffs(i).setValues( Coeffs(i) );
1191 : }
1192 16 : }
1193 14 : }
1194 :
1195 :
1196 85 : void Optimizer::addCoeffsSetIDsToFilenames(std::vector<std::string>& fnames, std::string& coeffssetid_prefix) {
1197 85 : if(ncoeffssets_==1) {
1198 : return;
1199 : }
1200 : //
1201 9 : if(fnames.size()==1) {
1202 0 : fnames.resize(ncoeffssets_,fnames[0]);
1203 : }
1204 9 : plumed_assert(fnames.size()==ncoeffssets_);
1205 : //
1206 36 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1207 27 : std::string is="";
1208 27 : Tools::convert(i,is);
1209 54 : fnames[i] = FileBase::appendSuffix(fnames[i],"."+coeffssetid_prefix_+is);
1210 : }
1211 : }
1212 :
1213 :
1214 93 : void Optimizer::setAllCoeffsSetIterationCounters() {
1215 194 : for(unsigned int i=0; i<ncoeffssets_; i++) {
1216 101 : coeffs_pntrs_[i]->setIterationCounterAndTime(getIterationCounter(),getTime());
1217 101 : aux_coeffs_pntrs_[i]->setIterationCounterAndTime(getIterationCounter(),getTime());
1218 101 : gradient_pntrs_[i]->setIterationCounterAndTime(getIterationCounter(),getTime());
1219 101 : targetdist_averages_pntrs_[i]->setIterationCounterAndTime(getIterationCounter(),getTime());
1220 101 : if(use_hessian_) {
1221 0 : hessian_pntrs_[i]->setIterationCounterAndTime(getIterationCounter(),getTime());
1222 : }
1223 : }
1224 93 : }
1225 :
1226 :
1227 95 : std::string Optimizer::getIterationCounterStr(const int offset) const {
1228 : std::string s;
1229 95 : Tools::convert(getIterationCounter()+offset,s);
1230 95 : return s;
1231 : }
1232 :
1233 :
1234 75 : void Optimizer::writeBiasOutputFiles() const {
1235 156 : for(unsigned int i=0; i<nbiases_; i++) {
1236 81 : bias_pntrs_[i]->writeBiasToFile();
1237 : }
1238 75 : }
1239 :
1240 :
1241 75 : void Optimizer::writeFesOutputFiles() const {
1242 156 : for(unsigned int i=0; i<nbiases_; i++) {
1243 81 : bias_pntrs_[i]->writeFesToFile();
1244 : }
1245 75 : }
1246 :
1247 :
1248 17 : void Optimizer::writeFesProjOutputFiles() const {
1249 34 : for(unsigned int i=0; i<nbiases_; i++) {
1250 17 : bias_pntrs_[i]->writeFesProjToFile();
1251 : }
1252 17 : }
1253 :
1254 :
1255 38 : void Optimizer::writeTargetDistOutputFiles() const {
1256 76 : for(unsigned int i=0; i<nbiases_; i++) {
1257 38 : if(dynamic_targetdists_[i]) {
1258 38 : bias_pntrs_[i]->writeTargetDistToFile();
1259 : }
1260 : }
1261 38 : }
1262 :
1263 :
1264 4 : void Optimizer::writeTargetDistProjOutputFiles() const {
1265 8 : for(unsigned int i=0; i<nbiases_; i++) {
1266 4 : if(dynamic_targetdists_[i]) {
1267 4 : bias_pntrs_[i]->writeTargetDistProjToFile();
1268 : }
1269 : }
1270 4 : }
1271 :
1272 :
1273 :
1274 : }
1275 : }
|
