Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "GeneralReporter.h"
11 : #include "OptimizationData.h"
12 :
13 : registerMooseObject("OptimizationApp", OptimizationData);
14 :
15 : template <typename T>
16 : InputParameters
17 2336 : OptimizationDataTempl<T>::validParams()
18 : {
19 2336 : InputParameters params = T::validParams();
20 :
21 2336 : params.addClassDescription(
22 : "Reporter to hold measurement and simulation data for optimization problems");
23 4672 : params.addParam<std::vector<Real>>(
24 : "measurement_values",
25 : "Measurement values collected from locations given by measurement_points");
26 4672 : params.addParam<std::vector<Point>>("measurement_points",
27 : "Point locations corresponding to each measurement value");
28 4672 : params.addParam<std::vector<Real>>("measurement_times",
29 : "Times corresponding to each measurement value");
30 :
31 4672 : params.addParam<FileName>("measurement_file",
32 : "CSV file with measurement value and coordinates (value, x, y, z).");
33 4672 : params.addParam<std::string>(
34 : "file_xcoord", "x", "x coordinate column name from measurement_file csv being read in.");
35 4672 : params.addParam<std::string>(
36 : "file_ycoord", "y", "y coordinate column name from csv file being read in.");
37 4672 : params.addParam<std::string>(
38 : "file_zcoord", "z", "z coordinate column name from csv file being read in.");
39 4672 : params.addParam<std::string>(
40 : "file_time", "time", "time column name from csv file being read in.");
41 4672 : params.addParam<std::vector<std::string>>(
42 : "file_variable_weights", {}, "variable weight column names from csv file being read in.");
43 4672 : params.addParam<std::string>(
44 : "file_value", "value", "measurement value column name from csv file being read in.");
45 :
46 4672 : params.addParam<std::vector<std::string>>(
47 : "variable_weight_names",
48 : "Vector of weight reporter names that will create a reporter to transfer weights into. The "
49 : "ordering of these weight reporter names corresponds to the ordering used in variable.");
50 4672 : params.addParam<std::vector<VariableName>>(
51 : "variable", "Vector of variable names to sample at measurement points.");
52 4672 : params.addParam<ReporterValueName>("objective_name",
53 : "Name of reporter value defining the objective.");
54 4672 : params.addParamNamesToGroup("measurement_points measurement_values measurement_times",
55 : "Input Measurement Data");
56 4672 : params.addParamNamesToGroup("measurement_file file_xcoord file_ycoord file_zcoord file_time "
57 : "file_value file_variable_weights",
58 : "File Measurement Data");
59 2336 : return params;
60 0 : }
61 :
62 : template <typename T>
63 1168 : OptimizationDataTempl<T>::OptimizationDataTempl(const InputParameters & parameters)
64 : : T(parameters),
65 2336 : _measurement_xcoord(this->template declareValueByName<std::vector<Real>>(
66 : "measurement_xcoord", REPORTER_MODE_REPLICATED)),
67 2336 : _measurement_ycoord(this->template declareValueByName<std::vector<Real>>(
68 : "measurement_ycoord", REPORTER_MODE_REPLICATED)),
69 2336 : _measurement_zcoord(this->template declareValueByName<std::vector<Real>>(
70 : "measurement_zcoord", REPORTER_MODE_REPLICATED)),
71 2336 : _measurement_time(this->template declareValueByName<std::vector<Real>>(
72 : "measurement_time", REPORTER_MODE_REPLICATED)),
73 2336 : _measurement_values(this->template declareValueByName<std::vector<Real>>(
74 : "measurement_values", REPORTER_MODE_REPLICATED)),
75 2336 : _simulation_values(this->template declareValueByName<std::vector<Real>>(
76 : "simulation_values", REPORTER_MODE_REPLICATED)),
77 1168 : _misfit_values(this->template declareValueByName<std::vector<Real>>("misfit_values",
78 : REPORTER_MODE_REPLICATED)),
79 3504 : _objective_val(this->isParamSetByUser("objective_name")
80 2434 : ? this->template declareValueByName<Real>(
81 : this->template getParam<ReporterValueName>("objective_name"),
82 : REPORTER_MODE_REPLICATED)
83 1914 : : this->template declareUnusedValue<Real>())
84 : {
85 : // read in data
86 2706 : if (this->isParamValid("measurement_file") && this->isParamValid("measurement_points"))
87 0 : mooseError("Input file can only define a single input for measurement data. Use only "
88 : "measurement_file or measurement_points, but never both");
89 2336 : else if (this->isParamValid("measurement_file"))
90 185 : readMeasurementsFromFile();
91 1966 : else if (this->isParamValid("measurement_points"))
92 465 : readMeasurementsFromInput();
93 :
94 1168 : _misfit_values.resize(_measurement_values.size());
95 :
96 2336 : if (this->isParamValid("variable"))
97 : {
98 1707 : std::vector<VariableName> var_names(
99 : this->template getParam<std::vector<VariableName>>("variable"));
100 1176 : for (const auto & name : var_names)
101 607 : _var_vec.push_back(&this->_fe_problem.getVariable(
102 607 : this->_tid, name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD));
103 569 : }
104 2336 : if (this->isParamValid("variable_weight_names"))
105 : {
106 165 : std::vector<std::string> weight_names(
107 : this->template getParam<std::vector<std::string>>("variable_weight_names"));
108 148 : for (const auto & name : weight_names)
109 : {
110 : if (_weight_names_weights_map.count(name) == 1)
111 : {
112 66 : _variable_weights.push_back(_weight_names_weights_map[name]);
113 : }
114 : else
115 : {
116 : // default is to create a new weight reporter and fill it with 1's
117 : // these will be overwritten by a reporter transfer.
118 27 : _variable_weights.push_back(
119 81 : &this->template declareValueByName<std::vector<Real>>(name, REPORTER_MODE_REPLICATED));
120 27 : _variable_weights.back()->assign(_measurement_xcoord.size(), 1);
121 : }
122 : }
123 55 : }
124 3520 : if (this->isParamValid("variable") && this->isParamValid("variable_weight_names") &&
125 : _variable_weights.size() != _var_vec.size())
126 : {
127 0 : this->paramError(
128 : "variable_weight_names",
129 : "The same number of names must be in both 'variable_weight_names' and 'variable'.");
130 : }
131 1168 : }
132 :
133 : template <typename T>
134 : void
135 45897 : OptimizationDataTempl<T>::execute()
136 : {
137 45897 : computeMisfit();
138 45897 : _objective_val = computeMisfitValue();
139 45897 : }
140 :
141 : template <typename T>
142 : void
143 45897 : OptimizationDataTempl<T>::computeMisfit()
144 : {
145 45897 : if (_var_vec.empty())
146 : return;
147 :
148 : // FIXME: This is basically copied from PointValue.
149 : // Implementation can be improved using the functionality in PointSamplerBase,
150 : // but this will require changes in MOOSE to work for reporters.
151 :
152 31389 : const std::size_t nvals = _measurement_values.size();
153 31389 : _simulation_values.resize(nvals, 0.0);
154 31389 : _misfit_values.resize(nvals);
155 :
156 31389 : errorCheckDataSize();
157 94223 : for (const auto var_index : make_range(_var_vec.size()))
158 : {
159 31417 : const auto & sys = _var_vec[var_index]->sys().system();
160 31417 : const auto vnum = _var_vec[var_index]->number();
161 : // A weight reporter is not automatically created and for those cases, we
162 : // set the weight to 1.
163 276 : std::vector<Real> weights(_variable_weights.empty()
164 31417 : ? std::vector<Real>(_measurement_xcoord.size(), 1)
165 : : (*_variable_weights[var_index]));
166 1284249 : for (const auto & i : make_range(nvals))
167 : {
168 1252832 : if (MooseUtils::absoluteFuzzyEqual(this->_t, _measurement_time[i]))
169 : {
170 : // If we are on the first var, make sure reset the simulation values so they aren't
171 : // accumulated on repeated timesteps
172 151238 : if (var_index == 0)
173 150846 : _simulation_values[i] = 0.0;
174 :
175 151238 : const Point point(_measurement_xcoord[i], _measurement_ycoord[i], _measurement_zcoord[i]);
176 151238 : const Real val = sys.point_value(vnum, point, false);
177 :
178 151238 : _simulation_values[i] += weights[i] * val;
179 151238 : _misfit_values[i] = _simulation_values[i] - _measurement_values[i];
180 : }
181 : }
182 : }
183 : }
184 :
185 : template <typename T>
186 : void
187 185 : OptimizationDataTempl<T>::readMeasurementsFromFile()
188 : {
189 370 : std::string xName = this->template getParam<std::string>("file_xcoord");
190 370 : std::string yName = this->template getParam<std::string>("file_ycoord");
191 370 : std::string zName = this->template getParam<std::string>("file_zcoord");
192 370 : std::string tName = this->template getParam<std::string>("file_time");
193 370 : std::string valueName = this->template getParam<std::string>("file_value");
194 370 : std::vector<std::string> weightNames =
195 : this->template getParam<std::vector<std::string>>("file_variable_weights");
196 :
197 : bool found_x = false;
198 : bool found_y = false;
199 : bool found_z = false;
200 : bool found_t = false;
201 : bool found_value = false;
202 :
203 370 : MooseUtils::DelimitedFileReader reader(this->template getParam<FileName>("measurement_file"));
204 185 : reader.read();
205 :
206 185 : auto const & names = reader.getNames();
207 185 : auto const & data = reader.getData();
208 :
209 : const std::size_t rows = data[0].size();
210 1522 : for (std::size_t i = 0; i < names.size(); ++i)
211 : {
212 : // make sure all data columns have the same length
213 1337 : if (data[i].size() != rows)
214 0 : this->paramError("file", "Mismatching column lengths in file");
215 :
216 1337 : if (names[i] == xName)
217 : {
218 185 : _measurement_xcoord = data[i];
219 : found_x = true;
220 : }
221 1152 : else if (names[i] == yName)
222 : {
223 185 : _measurement_ycoord = data[i];
224 : found_y = true;
225 : }
226 967 : else if (names[i] == zName)
227 : {
228 185 : _measurement_zcoord = data[i];
229 : found_z = true;
230 : }
231 782 : else if (names[i] == tName)
232 : {
233 73 : _measurement_time = data[i];
234 : found_t = true;
235 : }
236 709 : else if (names[i] == valueName)
237 : {
238 185 : _measurement_values = data[i];
239 : found_value = true;
240 : }
241 524 : else if (std::find(weightNames.begin(), weightNames.end(), names[i]) != weightNames.end())
242 : {
243 75 : _weight_names_weights_map.emplace(names[i],
244 300 : &(this->template declareValueByName<std::vector<Real>>(
245 : names[i], REPORTER_MODE_REPLICATED)));
246 75 : _weight_names_weights_map[names[i]]->assign(data[i].begin(), data[i].end());
247 : }
248 : }
249 :
250 : // check if all required columns were found
251 185 : if (!found_x)
252 0 : this->paramError(
253 : "measurement_file", "Column with name '", xName, "' missing from measurement file");
254 185 : if (!found_y)
255 0 : this->paramError(
256 : "measurement_file", "Column with name '", yName, "' missing from measurement file");
257 185 : if (!found_z)
258 0 : this->paramError(
259 : "measurement_file", "Column with name '", zName, "' missing from measurement file");
260 185 : if (!found_t)
261 112 : _measurement_time.assign(rows, 0);
262 185 : if (!found_value)
263 0 : this->paramError(
264 : "measurement_file", "Column with name '", valueName, "' missing from measurement file");
265 185 : if (_weight_names_weights_map.size() != weightNames.size())
266 : {
267 0 : std::string out("\n Measurement file column names: ");
268 0 : for (const auto & name : names)
269 0 : out += " " + name;
270 : out += "\n file_variable_weights names: ";
271 0 : for (const auto & name : weightNames)
272 0 : out += " " + name;
273 0 : this->paramError(
274 : "measurement_file",
275 : "Not all of the file_variable_weights names were found in the measurement_file.",
276 : out);
277 : }
278 370 : }
279 :
280 : template <typename T>
281 : void
282 465 : OptimizationDataTempl<T>::readMeasurementsFromInput()
283 : {
284 930 : if (!this->template getParam<std::vector<std::string>>("file_variable_weights").empty())
285 0 : this->paramError(
286 : "measurement_values",
287 : "file_variable_weights cannot be used with measurement data read from the input "
288 : "file, use measure_file input instead.");
289 :
290 3539 : for (const auto & p : this->template getParam<std::vector<Point>>("measurement_points"))
291 : {
292 2609 : _measurement_xcoord.push_back(p(0));
293 2609 : _measurement_ycoord.push_back(p(1));
294 2609 : _measurement_zcoord.push_back(p(2));
295 : }
296 :
297 930 : if (this->isParamValid("measurement_times"))
298 81 : _measurement_time = this->template getParam<std::vector<Real>>("measurement_times");
299 : else
300 438 : _measurement_time.assign(_measurement_xcoord.size(), 0.0);
301 :
302 930 : if (this->isParamValid("measurement_values"))
303 930 : _measurement_values = this->template getParam<std::vector<Real>>("measurement_values");
304 : else
305 0 : this->paramError("measurement_values", "Input file must contain measurement points and values");
306 465 : }
307 :
308 : template <typename T>
309 : void
310 31389 : OptimizationDataTempl<T>::errorCheckDataSize()
311 : {
312 31389 : const std::size_t nvals = _measurement_values.size();
313 31389 : std::string msg = "";
314 31389 : if (_measurement_xcoord.size() != nvals)
315 0 : msg += "x-coordinate data (" + std::to_string(_measurement_xcoord.size()) + "), ";
316 31389 : if (_measurement_ycoord.size() != nvals)
317 0 : msg += "y-coordinate data (" + std::to_string(_measurement_ycoord.size()) + "), ";
318 31389 : if (_measurement_zcoord.size() != nvals)
319 0 : msg += "z-coordinate data (" + std::to_string(_measurement_zcoord.size()) + "), ";
320 31389 : if (_measurement_time.size() != nvals)
321 0 : msg += "time data (" + std::to_string(_measurement_time.size()) + "), ";
322 31389 : if (!msg.empty())
323 0 : mooseError("Number of entries in ",
324 0 : std::string(msg.begin(), msg.end() - 2),
325 : " does not match number of entries in value data (",
326 : std::to_string(nvals),
327 : ").");
328 31389 : }
329 :
330 : template <typename T>
331 : Real
332 45897 : OptimizationDataTempl<T>::computeMisfitValue()
333 : {
334 : Real val = 0.0;
335 1566102 : for (auto & misfit : _misfit_values)
336 1520205 : val += misfit * misfit;
337 :
338 45897 : return val * 0.5;
339 : }
340 :
341 : template class OptimizationDataTempl<GeneralReporter>;
342 : template class OptimizationDataTempl<OptimizationReporterBase>;
|
