Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2018 The VES code team
3 : (see the PEOPLE-VES file at the root of this folder for a list of names)
4 :
5 : See http://www.ves-code.org for more information.
6 :
7 : This file is part of VES code module.
8 :
9 : The VES code module is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : The VES code module is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with the VES code module. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 :
23 : #include "TargetDistribution.h"
24 : #include "VesTools.h"
25 :
26 :
27 : #include "core/ActionRegister.h"
28 : #include "core/ActionSet.h"
29 : #include "core/PlumedMain.h"
30 : #include "tools/Grid.h"
31 :
32 : #include "GridIntegrationWeights.h"
33 :
34 :
35 : namespace PLMD {
36 :
37 : // class Grid;
38 : class Action;
39 :
40 : namespace ves {
41 :
42 : //+PLUMEDOC VES_TARGETDIST TD_PRODUCT_COMBINATION
43 : /*
44 : Target distribution given by product combination of distributions (static or dynamic).
45 :
46 : Employ a target distribution that is a product combination of the other
47 : distributions, defined as
48 : \f[
49 : p(\mathbf{s}) =
50 : \frac{\prod_{i} p_{i}(\mathbf{s})}
51 : {\int d \mathbf{s} \prod_{i} p_{i}(\mathbf{s})}
52 : \f]
53 : where the distributions \f$p_{i}(\mathbf{s})\f$ are in full dimensional space
54 : of the arguments used.
55 :
56 : Note the difference between this target distribution and the one defined in
57 : \ref TD_PRODUCT_DISTRIBUTION. Here we have a non-separable distribution given
58 : as a product of distribution \f$p_{i}(\mathbf{s})\f$ which are in full dimensional
59 : space of the arguments used.
60 :
61 : The labels of the distributions \f$p_{i}(\mathbf{s})\f$ to be used in the
62 : product combination are given in the DISTRIBUTIONS keyword.
63 :
64 : The target distribution resulting from the product combination will be
65 : automatically normalized. Therefore, the product combination needs to
66 : be a proper distribution that is non-negative and normalizable. The
67 : code will perform checks to make sure that this is indeed the case.
68 :
69 : The product 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 : \par Examples
74 :
75 : In the following example the overall interval on which the
76 : target distribution is defined is from 0.23 to 0.8.
77 : We employ a product combination of a well-tempered
78 : distribution and a uniform distribution that decays to
79 : zero at 0.6. This results in a target distribution that
80 : is well-tempered from 0.23 to 0.6 and then decays to zero.
81 : In other words, we cut off the tail of the well-tempered
82 : distribution at 0.6
83 : \plumedfile
84 : td_welltemp: TD_WELLTEMPERED BIASFACTOR=5
85 : td_uniform: TD_UNIFORM MINIMA=0.23 MAXIMA=0.6 SIGMA_MAXIMA=0.05
86 : td_combination: TD_PRODUCT_COMBINATION DISTRIBUTIONS=td_uniform,td_welltemp
87 : \endplumedfile
88 :
89 :
90 : In the following example the overall interval on which the
91 : target distribution is defined is from -4 to 4.
92 : We employ a product of a Gaussian distribution with two centers
93 : and distribution that is uniform on the interval -3 to 3 and
94 : then smoothly decays to zero outside that interval.
95 : The overall effect will then be to cut off the tails of the
96 : Gaussian distribution
97 : \plumedfile
98 : TD_GAUSSIAN ...
99 : CENTER1=-2.9 SIGMA1=1.0
100 : CENTER2=+2.9 SIGMA2=0.4
101 : LABEL=td_gauss
102 : ... TD_GAUSSIAN
103 :
104 : TD_UNIFORM ...
105 : MINIMA=-3.0 SIGMA_MINIMA=0.20
106 : MAXIMA=+3.0 SIGMA_MAXIMA=0.15
107 : LABEL=td_uni
108 : ... TD_UNIFORM
109 :
110 : td_pc: TD_PRODUCT_COMBINATION DISTRIBUTIONS=td_gauss,td_uni
111 : \endplumedfile
112 :
113 : */
114 : //+ENDPLUMEDOC
115 :
116 : class VesBias;
117 :
118 12 : class TD_ProductCombination: public TargetDistribution {
119 : private:
120 : std::vector<TargetDistribution*> distribution_pntrs_;
121 : std::vector<Grid*> grid_pntrs_;
122 : unsigned int ndist_;
123 : void setupAdditionalGrids(const std::vector<Value*>&, const std::vector<std::string>&, const std::vector<std::string>&, const std::vector<unsigned int>&);
124 : public:
125 : static void registerKeywords(Keywords&);
126 : explicit TD_ProductCombination(const ActionOptions& ao);
127 : void updateGrid();
128 : double getValue(const std::vector<double>&) const;
129 : //
130 : void linkVesBias(VesBias*);
131 : void linkAction(Action*);
132 : //
133 : void linkBiasGrid(Grid*);
134 : void linkBiasWithoutCutoffGrid(Grid*);
135 : void linkFesGrid(Grid*);
136 : //
137 : };
138 :
139 :
140 6456 : PLUMED_REGISTER_ACTION(TD_ProductCombination,"TD_PRODUCT_COMBINATION")
141 :
142 :
143 5 : void TD_ProductCombination::registerKeywords(Keywords& keys) {
144 5 : TargetDistribution::registerKeywords(keys);
145 20 : keys.add("compulsory","DISTRIBUTIONS","The labels of the target distribution actions to be used in the product combination.");
146 10 : keys.use("WELLTEMPERED_FACTOR");
147 10 : keys.use("SHIFT_TO_ZERO");
148 5 : }
149 :
150 :
151 4 : TD_ProductCombination::TD_ProductCombination(const ActionOptions& ao):
152 : PLUMED_VES_TARGETDISTRIBUTION_INIT(ao),
153 : distribution_pntrs_(0),
154 : grid_pntrs_(0),
155 4 : ndist_(0)
156 : {
157 4 : std::vector<std::string> targetdist_labels;
158 8 : parseVector("DISTRIBUTIONS",targetdist_labels);
159 :
160 4 : std::string error_msg = "";
161 12 : distribution_pntrs_ = VesTools::getPointersFromLabels<TargetDistribution*>(targetdist_labels,plumed.getActionSet(),error_msg);
162 4 : if(error_msg.size()>0) {plumed_merror("Error in keyword DISTRIBUTIONS of "+getName()+": "+error_msg);}
163 :
164 32 : for(unsigned int i=0; i<distribution_pntrs_.size(); i++) {
165 16 : if(distribution_pntrs_[i]->isDynamic()) {setDynamic();}
166 16 : if(distribution_pntrs_[i]->fesGridNeeded()) {setFesGridNeeded();}
167 16 : if(distribution_pntrs_[i]->biasGridNeeded()) {setBiasGridNeeded();}
168 : }
169 :
170 4 : ndist_ = distribution_pntrs_.size();
171 8 : grid_pntrs_.assign(ndist_,NULL);
172 4 : if(ndist_==0) {plumed_merror(getName()+ ": no distributions are given.");}
173 4 : if(ndist_==1) {plumed_merror(getName()+ ": giving only one distribution does not make sense.");}
174 : //
175 4 : checkRead();
176 4 : }
177 :
178 :
179 0 : double TD_ProductCombination::getValue(const std::vector<double>& argument) const {
180 0 : plumed_merror("getValue not implemented for TD_ProductCombination");
181 : return 0.0;
182 : }
183 :
184 :
185 4 : void TD_ProductCombination::setupAdditionalGrids(const std::vector<Value*>& arguments, const std::vector<std::string>& min, const std::vector<std::string>& max, const std::vector<unsigned int>& nbins) {
186 20 : for(unsigned int i=0; i<ndist_; i++) {
187 16 : distribution_pntrs_[i]->setupGrids(arguments,min,max,nbins);
188 16 : if(distribution_pntrs_[i]->getDimension()!=this->getDimension()) {
189 0 : plumed_merror(getName() + ": all target distribution must have the same dimension");
190 : }
191 8 : grid_pntrs_[i]=distribution_pntrs_[i]->getTargetDistGridPntr();
192 : }
193 4 : }
194 :
195 :
196 14 : void TD_ProductCombination::updateGrid() {
197 70 : for(unsigned int i=0; i<ndist_; i++) {
198 56 : distribution_pntrs_[i]->updateTargetDist();
199 : }
200 56 : std::vector<double> integration_weights = GridIntegrationWeights::getIntegrationWeights(getTargetDistGridPntr());
201 : double norm = 0.0;
202 23420 : for(Grid::index_t l=0; l<targetDistGrid().getSize(); l++) {
203 : double value = 1.0;
204 58515 : for(unsigned int i=0; i<ndist_; i++) {
205 46812 : value *= grid_pntrs_[i]->getValue(l);
206 : }
207 11703 : if(value<0.0 && !isTargetDistGridShiftedToZero()) {plumed_merror(getName()+": The target distribution function gives negative values. You should change the definition of the target distribution to avoid this. You can also use the SHIFT_TO_ZERO keyword to avoid this problem.");}
208 11703 : norm += integration_weights[l]*value;
209 11703 : targetDistGrid().setValue(l,value);
210 11703 : logTargetDistGrid().setValue(l,-std::log(value));
211 : }
212 :
213 14 : if(norm>0.0) {
214 14 : targetDistGrid().scaleAllValuesAndDerivatives(1.0/norm);
215 : }
216 0 : else if(!isTargetDistGridShiftedToZero()) {
217 0 : plumed_merror(getName()+": The target distribution function cannot be normalized proberly. You should change the definition of the target distribution to avoid this. You can also use the SHIFT_TO_ZERO keyword to avoid this problem.");
218 : }
219 14 : logTargetDistGrid().setMinToZero();
220 14 : }
221 :
222 :
223 1 : void TD_ProductCombination::linkVesBias(VesBias* vesbias_pntr_in) {
224 1 : TargetDistribution::linkVesBias(vesbias_pntr_in);
225 5 : for(unsigned int i=0; i<ndist_; i++) {
226 4 : distribution_pntrs_[i]->linkVesBias(vesbias_pntr_in);
227 : }
228 1 : }
229 :
230 :
231 0 : void TD_ProductCombination::linkAction(Action* action_pntr_in) {
232 0 : TargetDistribution::linkAction(action_pntr_in);
233 0 : for(unsigned int i=0; i<ndist_; i++) {
234 0 : distribution_pntrs_[i]->linkAction(action_pntr_in);
235 : }
236 0 : }
237 :
238 :
239 0 : void TD_ProductCombination::linkBiasGrid(Grid* bias_grid_pntr_in) {
240 0 : TargetDistribution::linkBiasGrid(bias_grid_pntr_in);
241 0 : for(unsigned int i=0; i<ndist_; i++) {
242 0 : distribution_pntrs_[i]->linkBiasGrid(bias_grid_pntr_in);
243 : }
244 0 : }
245 :
246 :
247 0 : void TD_ProductCombination::linkBiasWithoutCutoffGrid(Grid* bias_withoutcutoff_grid_pntr_in) {
248 0 : TargetDistribution::linkBiasWithoutCutoffGrid(bias_withoutcutoff_grid_pntr_in);
249 0 : for(unsigned int i=0; i<ndist_; i++) {
250 0 : distribution_pntrs_[i]->linkBiasWithoutCutoffGrid(bias_withoutcutoff_grid_pntr_in);
251 : }
252 0 : }
253 :
254 :
255 1 : void TD_ProductCombination::linkFesGrid(Grid* fes_grid_pntr_in) {
256 1 : TargetDistribution::linkFesGrid(fes_grid_pntr_in);
257 5 : for(unsigned int i=0; i<ndist_; i++) {
258 4 : distribution_pntrs_[i]->linkFesGrid(fes_grid_pntr_in);
259 : }
260 1 : }
261 :
262 :
263 : }
264 4839 : }
|