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