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 :
25 : #include "GridIntegrationWeights.h"
26 : #include "VesTools.h"
27 :
28 : #include "core/ActionRegister.h"
29 : #include "core/PlumedMain.h"
30 : #include "core/Value.h"
31 : #include "tools/File.h"
32 : #include "tools/Grid.h"
33 :
34 :
35 : namespace PLMD {
36 : namespace ves {
37 :
38 : //+PLUMEDOC VES_UTILS VES_OUTPUT_TARGET_DISTRIBUTION
39 : /*
40 : Output target distribution to file.
41 :
42 : This action can be used to output target distributions to a grid file,
43 : for example to see how they look like before using them in a VES bias.
44 : This action only support static target distributions.
45 :
46 : This action is normally used through the \ref driver.
47 :
48 :
49 : \par Examples
50 :
51 : In the following input we define a target distribution that is uniform for
52 : argument 1 and a Gaussian for argument 2 and then output it to a file
53 : called targetdist-1.data.
54 : \plumedfile
55 : t1_1: TD_UNIFORM MINIMA=-4.0 MAXIMA=+4.0
56 : t1_2: TD_GAUSSIAN CENTER1=-2.0 SIGMA1=0.5
57 : t1: TD_PRODUCT_DISTRIBUTION DISTRIBUTIONS=t1_1,t1_2
58 :
59 : VES_OUTPUT_TARGET_DISTRIBUTION ...
60 : GRID_MIN=-4.0,-4.0
61 : GRID_MAX=+4.0,+4.0
62 : GRID_BINS=100,100
63 : TARGET_DISTRIBUTION=t1
64 : TARGETDIST_FILE=targetdist-1.data
65 : LOG_TARGETDIST_FILE=targetdist-1.log.data
66 : FMT_GRIDS=%11.6f
67 : ... VES_OUTPUT_TARGET_DISTRIBUTION
68 : \endplumedfile
69 :
70 : This input should be run through the driver by using a command similar to the
71 : following one where the trajectory/configuration file configuration.gro is needed to
72 : trick the code to exit correctly.
73 : \verbatim
74 : plumed driver --plumed plumed.dat --igro configuration.gro
75 : \endverbatim
76 :
77 : */
78 : //+ENDPLUMEDOC
79 :
80 :
81 : class OutputTargetDistribution :
82 : public Action {
83 : public:
84 : explicit OutputTargetDistribution(const ActionOptions&);
85 0 : void calculate() override {}
86 0 : void apply() override {}
87 : static void registerKeywords(Keywords& keys);
88 : };
89 :
90 :
91 : PLUMED_REGISTER_ACTION(OutputTargetDistribution,"VES_OUTPUT_TARGET_DISTRIBUTION")
92 :
93 122 : void OutputTargetDistribution::registerKeywords(Keywords& keys) {
94 122 : Action::registerKeywords(keys);
95 122 : keys.add("compulsory","GRID_MIN","the lower bounds for the grid");
96 122 : keys.add("compulsory","GRID_MAX","the upper bounds for the grid");
97 122 : keys.add("compulsory","GRID_BINS","the number of bins used for the grid.");
98 122 : keys.add("optional","GRID_PERIODICITY","specify if the individual arguments should be made periodic (YES) or not (NO). By default all arguments are taken as not periodic.");
99 122 : keys.add("compulsory","TARGETDIST_FILE","filename of the file for writing the target distribution");
100 122 : keys.add("optional","LOG_TARGETDIST_FILE","filename of the file for writing the log of the target distribution");
101 122 : keys.add("compulsory","TARGET_DISTRIBUTION","the target distribution to be used.");
102 122 : keys.add("optional","FMT_GRIDS","the numerical format of the target distribution grids written to file. By default it is %14.9f");
103 122 : keys.addFlag("DO_1D_PROJECTIONS",false,"Also output the one-dimensional marginal distributions for multi-dimensional target distribution.");
104 122 : }
105 :
106 120 : OutputTargetDistribution::OutputTargetDistribution(const ActionOptions&ao):
107 120 : Action(ao) {
108 :
109 : std::string targetdist_fname;
110 240 : parse("TARGETDIST_FILE",targetdist_fname);
111 : std::string log_targetdist_fname;
112 120 : parse("LOG_TARGETDIST_FILE",log_targetdist_fname);
113 120 : if(targetdist_fname==log_targetdist_fname) {
114 0 : plumed_merror("error in " + getName() + ":TARGETDIST_FILE and LOG_TARGETDIST_FILE cannot be the same");
115 : }
116 :
117 : std::vector<unsigned int> grid_bins;
118 240 : parseVector("GRID_BINS",grid_bins);
119 120 : unsigned int nargs = grid_bins.size();
120 :
121 120 : std::vector<std::string> grid_min(nargs);
122 120 : parseVector("GRID_MIN",grid_min);
123 120 : std::vector<std::string> grid_max(nargs);
124 120 : parseVector("GRID_MAX",grid_max);
125 :
126 120 : std::vector<std::string> grid_periodicity(nargs);
127 240 : parseVector("GRID_PERIODICITY",grid_periodicity);
128 120 : if(grid_periodicity.size()==0) {
129 224 : grid_periodicity.assign(nargs,"NO");
130 : }
131 :
132 120 : std::string fmt_grids="%14.9f";
133 120 : parse("FMT_GRIDS",fmt_grids);
134 :
135 120 : bool do_1d_proj = false;
136 120 : parseFlag("DO_1D_PROJECTIONS",do_1d_proj);
137 120 : if(do_1d_proj && nargs==1) {
138 0 : plumed_merror("doesn't make sense to use the DO_1D_PROJECTIONS keyword for a one-dimensional distribution");
139 : }
140 :
141 120 : plumed_massert(grid_min.size()==nargs,"mismatch between number of values given for grid parameters");
142 120 : plumed_massert(grid_max.size()==nargs,"mismatch between number of values given for grid parameters");
143 120 : plumed_massert(grid_periodicity.size()==nargs,"mismatch between number of values given for grid parameters");
144 :
145 : std::string targetdist_label;
146 120 : parse("TARGET_DISTRIBUTION",targetdist_label);
147 120 : checkRead();
148 : //
149 120 : std::vector<std::unique_ptr<Value>> arguments(nargs);
150 287 : for(unsigned int i=0; i < nargs; i++) {
151 : std::string is;
152 167 : Tools::convert(i+1,is);
153 167 : if(nargs==1) {
154 : is="";
155 : }
156 167 : arguments[i]= Tools::make_unique<Value>(nullptr,"arg"+is,false);
157 167 : if(grid_periodicity[i]=="YES") {
158 11 : arguments[i]->setDomain(grid_min[i],grid_max[i]);
159 156 : } else if(grid_periodicity[i]=="NO") {
160 156 : arguments[i]->setNotPeriodic();
161 : } else {
162 0 : plumed_merror("wrong value given in GRID_PERIODICITY, either specify YES or NO");
163 : }
164 : }
165 :
166 120 : std::string error_msg = "";
167 120 : TargetDistribution* targetdist_pntr = VesTools::getPointerFromLabel<TargetDistribution*>(targetdist_label,plumed.getActionSet(),error_msg);
168 120 : if(error_msg.size()>0) {
169 0 : plumed_merror("Error in keyword TARGET_DISTRIBUTION of "+getName()+": "+error_msg);
170 : }
171 : //
172 120 : if(targetdist_pntr->isDynamic()) {
173 0 : plumed_merror(getName() + " only works for static target distributions");
174 : }
175 120 : targetdist_pntr->setupGrids(Tools::unique2raw(arguments),grid_min,grid_max,grid_bins);
176 120 : targetdist_pntr->updateTargetDist();
177 : Grid* targetdist_grid_pntr = targetdist_pntr->getTargetDistGridPntr();
178 : Grid* log_targetdist_grid_pntr = targetdist_pntr->getLogTargetDistGridPntr();
179 :
180 :
181 120 : double sum_grid = TargetDistribution::integrateGrid(targetdist_grid_pntr);
182 120 : log.printf(" target distribution integrated over the grid: %16.12f\n",sum_grid);
183 120 : log.printf(" (%30.16e)\n",sum_grid);
184 : //
185 120 : OFile ofile;
186 120 : ofile.link(*this);
187 120 : ofile.enforceBackup();
188 120 : ofile.open(targetdist_fname);
189 : targetdist_grid_pntr->setOutputFmt(fmt_grids);
190 120 : targetdist_grid_pntr->writeToFile(ofile);
191 120 : ofile.close();
192 120 : if(log_targetdist_fname.size()>0) {
193 23 : OFile ofile2;
194 23 : ofile2.link(*this);
195 23 : ofile2.enforceBackup();
196 23 : ofile2.open(log_targetdist_fname);
197 : log_targetdist_grid_pntr->setOutputFmt(fmt_grids);
198 23 : log_targetdist_grid_pntr->writeToFile(ofile2);
199 23 : ofile2.close();
200 23 : }
201 :
202 120 : if(do_1d_proj) {
203 12 : for(unsigned int i=0; i<nargs; i++) {
204 8 : std::vector<std::string> arg1d(1);
205 8 : arg1d[0] = arguments[i]->getName();
206 8 : Grid marginal_grid = targetdist_pntr->getMarginal(arg1d);
207 : //
208 : std::string suffix;
209 8 : Tools::convert(i+1,suffix);
210 8 : suffix = "proj-" + suffix;
211 8 : std::string marginal_fname = FileBase::appendSuffix(targetdist_fname,"."+suffix);
212 : //
213 8 : OFile ofile3;
214 8 : ofile3.link(*this);
215 8 : ofile3.enforceBackup();
216 8 : ofile3.open(marginal_fname);
217 : marginal_grid.setOutputFmt(fmt_grids);
218 8 : marginal_grid.writeToFile(ofile3);
219 16 : }
220 : }
221 :
222 480 : }
223 :
224 :
225 :
226 :
227 :
228 : }
229 : }
|