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 "TargetDistribution.h"
24 : #include "VesTools.h"
25 :
26 : #include "core/ActionRegister.h"
27 : #include "core/ActionSet.h"
28 : #include "core/PlumedMain.h"
29 : #include "tools/Grid.h"
30 :
31 : #include <limits>
32 :
33 :
34 : namespace PLMD {
35 :
36 : // class Grid;
37 : class Action;
38 :
39 : namespace ves {
40 :
41 : //+PLUMEDOC VES_TARGETDIST TD_LINEAR_COMBINATION
42 : /*
43 : Target distribution given by linear combination of distributions (static or dynamic).
44 :
45 : Employ a target distribution that is a linear combination of the other
46 : distributions, defined as
47 :
48 : $$
49 : p(\mathbf{s}) = \sum_{i} w_{i} \, p_{i}(\mathbf{s})
50 : $$
51 :
52 : where the weights $w_{i}$ are normalized to 1, $\sum_{i}w_{i}=1$.
53 :
54 : The labels of the distributions $p_{i}(\mathbf{s})$ to be used in the
55 : linear combination are given in the DISTRIBUTIONS keyword.
56 :
57 : The weights $w_{i}$ can be given using
58 : the WEIGHTS keyword. The distributions are weighted equally if no weights are given.
59 :
60 : It is assumed that all the distributions $p_{i}(\mathbf{s})$ are normalized.
61 : If that is not the case for some reason should you
62 : normalize each distribution separately by using the NORMALIZE
63 : keyword when defining them in the input file (i.e. before the
64 : TD_LINEAR_COMBINATION action).
65 : Note that normalizing the overall
66 : linear combination will generally lead to different results than normalizing
67 : each distribution separately.
68 :
69 : The linear combination will be a dynamic target distribution if one or more
70 : of the distributions used is a dynamic distribution, otherwise it will be a
71 : static distribution.
72 :
73 : ## Examples
74 :
75 : Here we employ a linear combination of a uniform and a Gaussian distribution.
76 : No weights are given so the two distributions will be weighted equally.
77 :
78 : ```plumed
79 : td_uni: TD_UNIFORM
80 :
81 : td_gauss: TD_GAUSSIAN CENTER1=-2.0 SIGMA1=0.5
82 :
83 : td_comb: TD_LINEAR_COMBINATION DISTRIBUTIONS=td_uni,td_gauss
84 : ```
85 :
86 : Here we employ a linear combination of a uniform and two Gaussian distribution.
87 : The weights are automatically normalized to 1 such that giving
88 : WEIGHTS=1.0,1.0,2.0 as we do here is equal to giving WEIGHTS=0.25,0.25,0.50.
89 :
90 : ```plumed
91 : td_uni: TD_UNIFORM
92 :
93 : td_gauss1: TD_GAUSSIAN CENTER1=-2.0,-2.0 SIGMA1=0.5,0.3
94 :
95 : td_gauss2: TD_GAUSSIAN CENTER1=+2.0,+2.0 SIGMA1=0.3,0.5
96 :
97 : TD_LINEAR_COMBINATION ...
98 : DISTRIBUTIONS=td_uni,td_gauss1,td_gauss2
99 : WEIGHTS=1.0,1.0,2.0
100 : LABEL=td_comb
101 : ... TD_LINEAR_COMBINATION
102 : ```
103 :
104 : In the above example the two Gaussian kernels are given using two separate
105 : DISTRIBUTION keywords. As the [TD_GAUSSIAN](TD_GAUSSIAN.md) target distribution allows multiple
106 : centers is it also possible to use just one DISTRIBUTION keyword for the two
107 : Gaussian kernels. This is shown in the following example which will give the
108 : exact same result as the one above as the weights have been appropriately
109 : adjusted
110 :
111 : ```plumed
112 : td_uni: TD_UNIFORM
113 :
114 : TD_GAUSSIAN ...
115 : CENTER1=-2.0,-2.0 SIGMA1=0.5,0.3
116 : CENTER2=+2.0,+2.0 SIGMA2=0.3,0.5
117 : WEIGHTS=1.0,2.0
118 : LABEL=td_gauss
119 : ... TD_GAUSSIAN
120 :
121 : TD_LINEAR_COMBINATION ...
122 : DISTRIBUTIONS=td_uni,td_gauss
123 : WEIGHTS=0.25,0.75
124 : LABEL=td_comb
125 : ... TD_LINEAR_COMBINATION
126 : ```
127 :
128 : */
129 : //+ENDPLUMEDOC
130 :
131 : class VesBias;
132 :
133 : class TD_LinearCombination: public TargetDistribution {
134 : private:
135 : std::vector<TargetDistribution*> distribution_pntrs_;
136 : std::vector<Grid*> grid_pntrs_;
137 : std::vector<double> weights_;
138 : unsigned int ndist_;
139 : void setupAdditionalGrids(const std::vector<Value*>&, const std::vector<std::string>&, const std::vector<std::string>&, const std::vector<unsigned int>&) override;
140 : public:
141 : static void registerKeywords(Keywords&);
142 : explicit TD_LinearCombination(const ActionOptions& ao);
143 : void updateGrid() override;
144 : double getValue(const std::vector<double>&) const override;
145 : //
146 : void linkVesBias(VesBias*) override;
147 : void linkAction(Action*) override;
148 : //
149 : void linkBiasGrid(Grid*) override;
150 : void linkBiasWithoutCutoffGrid(Grid*) override;
151 : void linkFesGrid(Grid*) override;
152 : //
153 : };
154 :
155 :
156 : PLUMED_REGISTER_ACTION(TD_LinearCombination,"TD_LINEAR_COMBINATION")
157 :
158 :
159 14 : void TD_LinearCombination::registerKeywords(Keywords& keys) {
160 14 : TargetDistribution::registerKeywords(keys);
161 14 : keys.add("compulsory","DISTRIBUTIONS","The labels of the target distribution actions to be used in the linear combination.");
162 14 : keys.add("optional","WEIGHTS","The weights of target distributions. Have to be as many as the number of target distribution labels given in DISTRIBUTIONS. If no weights are given the distributions are weighted equally. The weights are automatically normalized to 1.");
163 14 : keys.use("WELLTEMPERED_FACTOR");
164 : //keys.use("SHIFT_TO_ZERO");
165 14 : keys.use("NORMALIZE");
166 14 : }
167 :
168 :
169 12 : TD_LinearCombination::TD_LinearCombination(const ActionOptions& ao):
170 : PLUMED_VES_TARGETDISTRIBUTION_INIT(ao),
171 24 : distribution_pntrs_(0),
172 12 : grid_pntrs_(0),
173 12 : weights_(0),
174 24 : ndist_(0) {
175 : std::vector<std::string> targetdist_labels;
176 12 : parseVector("DISTRIBUTIONS",targetdist_labels);
177 :
178 12 : std::string error_msg = "";
179 24 : distribution_pntrs_ = VesTools::getPointersFromLabels<TargetDistribution*>(targetdist_labels,plumed.getActionSet(),error_msg);
180 12 : if(error_msg.size()>0) {
181 0 : plumed_merror("Error in keyword DISTRIBUTIONS of "+getName()+": "+error_msg);
182 : }
183 :
184 40 : for(unsigned int i=0; i<distribution_pntrs_.size(); i++) {
185 28 : if(distribution_pntrs_[i]->isDynamic()) {
186 : setDynamic();
187 : }
188 28 : if(distribution_pntrs_[i]->fesGridNeeded()) {
189 : setFesGridNeeded();
190 : }
191 28 : if(distribution_pntrs_[i]->biasGridNeeded()) {
192 : setBiasGridNeeded();
193 : }
194 : }
195 :
196 12 : ndist_ = distribution_pntrs_.size();
197 12 : grid_pntrs_.assign(ndist_,NULL);
198 12 : if(ndist_==0) {
199 0 : plumed_merror(getName()+ ": no distributions are given.");
200 : }
201 12 : if(ndist_==1) {
202 0 : plumed_merror(getName()+ ": giving only one distribution does not make sense.");
203 : }
204 : //
205 24 : parseVector("WEIGHTS",weights_);
206 12 : if(weights_.size()==0) {
207 4 : weights_.assign(distribution_pntrs_.size(),1.0);
208 : }
209 12 : if(distribution_pntrs_.size()!=weights_.size()) {
210 0 : plumed_merror(getName()+ ": there has to be as many weights given in WEIGHTS as the number of target distribution labels given in DISTRIBUTIONS");
211 : }
212 : //
213 : double sum_weights=0.0;
214 40 : for(unsigned int i=0; i<weights_.size(); i++) {
215 28 : sum_weights+=weights_[i];
216 : }
217 40 : for(unsigned int i=0; i<weights_.size(); i++) {
218 28 : weights_[i]/=sum_weights;
219 : }
220 12 : checkRead();
221 12 : }
222 :
223 :
224 0 : double TD_LinearCombination::getValue(const std::vector<double>& argument) const {
225 0 : plumed_merror("getValue not implemented for TD_LinearCombination");
226 : return 0.0;
227 : }
228 :
229 :
230 12 : void TD_LinearCombination::setupAdditionalGrids(const std::vector<Value*>& arguments, const std::vector<std::string>& min, const std::vector<std::string>& max, const std::vector<unsigned int>& nbins) {
231 40 : for(unsigned int i=0; i<ndist_; i++) {
232 28 : distribution_pntrs_[i]->setupGrids(arguments,min,max,nbins);
233 28 : if(distribution_pntrs_[i]->getDimension()!=this->getDimension()) {
234 0 : plumed_merror(getName() + ": all target distribution must have the same dimension");
235 : }
236 28 : grid_pntrs_[i]=distribution_pntrs_[i]->getTargetDistGridPntr();
237 : }
238 12 : }
239 :
240 :
241 22 : void TD_LinearCombination::updateGrid() {
242 70 : for(unsigned int i=0; i<ndist_; i++) {
243 48 : distribution_pntrs_[i]->updateTargetDist();
244 : }
245 162033 : for(Grid::index_t l=0; l<targetDistGrid().getSize(); l++) {
246 : double value = 0.0;
247 506837 : for(unsigned int i=0; i<ndist_; i++) {
248 344826 : value += weights_[i]*grid_pntrs_[i]->getValue(l);
249 : }
250 162011 : if(value==0.0) {
251 : value=std::numeric_limits<double>::denorm_min();
252 : }
253 162011 : targetDistGrid().setValue(l,value);
254 162011 : logTargetDistGrid().setValue(l,-std::log(value));
255 : }
256 22 : logTargetDistGrid().setMinToZero();
257 22 : }
258 :
259 :
260 1 : void TD_LinearCombination::linkVesBias(VesBias* vesbias_pntr_in) {
261 1 : TargetDistribution::linkVesBias(vesbias_pntr_in);
262 3 : for(unsigned int i=0; i<ndist_; i++) {
263 2 : distribution_pntrs_[i]->linkVesBias(vesbias_pntr_in);
264 : }
265 1 : }
266 :
267 :
268 0 : void TD_LinearCombination::linkAction(Action* action_pntr_in) {
269 0 : TargetDistribution::linkAction(action_pntr_in);
270 0 : for(unsigned int i=0; i<ndist_; i++) {
271 0 : distribution_pntrs_[i]->linkAction(action_pntr_in);
272 : }
273 0 : }
274 :
275 :
276 0 : void TD_LinearCombination::linkBiasGrid(Grid* bias_grid_pntr_in) {
277 0 : TargetDistribution::linkBiasGrid(bias_grid_pntr_in);
278 0 : for(unsigned int i=0; i<ndist_; i++) {
279 0 : distribution_pntrs_[i]->linkBiasGrid(bias_grid_pntr_in);
280 : }
281 0 : }
282 :
283 :
284 0 : void TD_LinearCombination::linkBiasWithoutCutoffGrid(Grid* bias_withoutcutoff_grid_pntr_in) {
285 0 : TargetDistribution::linkBiasWithoutCutoffGrid(bias_withoutcutoff_grid_pntr_in);
286 0 : for(unsigned int i=0; i<ndist_; i++) {
287 0 : distribution_pntrs_[i]->linkBiasWithoutCutoffGrid(bias_withoutcutoff_grid_pntr_in);
288 : }
289 0 : }
290 :
291 :
292 1 : void TD_LinearCombination::linkFesGrid(Grid* fes_grid_pntr_in) {
293 1 : TargetDistribution::linkFesGrid(fes_grid_pntr_in);
294 3 : for(unsigned int i=0; i<ndist_; i++) {
295 2 : distribution_pntrs_[i]->linkFesGrid(fes_grid_pntr_in);
296 : }
297 1 : }
298 :
299 :
300 : }
301 : }
|