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 1549 : OptimizationDataTempl<T>::validParams()
18 : {
19 1549 : InputParameters params = T::validParams();
20 :
21 1549 : params.addClassDescription(
22 : "Reporter to hold measurement and simulation data for optimization problems");
23 3098 : params.addParam<std::vector<Real>>(
24 : "measurement_values",
25 : "Measurement values collected from locations given by measurement_points");
26 3098 : params.addParam<std::vector<Point>>("measurement_points",
27 : "Point locations corresponding to each measurement value");
28 3098 : params.addParam<std::vector<Real>>("measurement_times",
29 : "Times corresponding to each measurement value");
30 :
31 3098 : params.addParam<FileName>("measurement_file",
32 : "CSV file with measurement value and coordinates (value, x, y, z).");
33 3098 : params.addParam<std::string>(
34 : "file_xcoord", "x", "x coordinate column name from measurement_file csv being read in.");
35 3098 : params.addParam<std::string>(
36 : "file_ycoord", "y", "y coordinate column name from csv file being read in.");
37 3098 : params.addParam<std::string>(
38 : "file_zcoord", "z", "z coordinate column name from csv file being read in.");
39 3098 : params.addParam<std::string>(
40 : "file_time", "time", "time column name from csv file being read in.");
41 3098 : params.addParam<std::vector<std::string>>(
42 : "file_variable_weights", {}, "variable weight column names from csv file being read in.");
43 3098 : params.addParam<std::string>(
44 : "file_value", "value", "measurement value column name from csv file being read in.");
45 :
46 3098 : 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 3098 : params.addParam<std::vector<VariableName>>(
51 : "variable", "Vector of variable names to sample at measurement points.");
52 3098 : params.addParam<ReporterValueName>("objective_name",
53 : "Name of reporter value defining the objective.");
54 3098 : params.addParamNamesToGroup("measurement_points measurement_values measurement_times",
55 : "Input Measurement Data");
56 3098 : params.addParamNamesToGroup("measurement_file file_xcoord file_ycoord file_zcoord file_time "
57 : "file_value file_variable_weights",
58 : "File Measurement Data");
59 1549 : return params;
60 0 : }
61 :
62 : template <typename T>
63 774 : OptimizationDataTempl<T>::OptimizationDataTempl(const InputParameters & parameters)
64 : : T(parameters),
65 1548 : _measurement_xcoord(this->template declareValueByName<std::vector<Real>>(
66 : "measurement_xcoord", REPORTER_MODE_REPLICATED)),
67 1548 : _measurement_ycoord(this->template declareValueByName<std::vector<Real>>(
68 : "measurement_ycoord", REPORTER_MODE_REPLICATED)),
69 1548 : _measurement_zcoord(this->template declareValueByName<std::vector<Real>>(
70 : "measurement_zcoord", REPORTER_MODE_REPLICATED)),
71 1548 : _measurement_time(this->template declareValueByName<std::vector<Real>>(
72 : "measurement_time", REPORTER_MODE_REPLICATED)),
73 1548 : _measurement_values(this->template declareValueByName<std::vector<Real>>(
74 : "measurement_values", REPORTER_MODE_REPLICATED)),
75 1548 : _simulation_values(this->template declareValueByName<std::vector<Real>>(
76 : "simulation_values", REPORTER_MODE_REPLICATED)),
77 774 : _misfit_values(this->template declareValueByName<std::vector<Real>>("misfit_values",
78 : REPORTER_MODE_REPLICATED)),
79 2322 : _objective_val(this->isParamSetByUser("objective_name")
80 1638 : ? this->template declareValueByName<Real>(
81 : this->template getParam<ReporterValueName>("objective_name"),
82 : REPORTER_MODE_REPLICATED)
83 1260 : : this->template declareUnusedValue<Real>())
84 : {
85 : // read in data
86 1828 : 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 1548 : else if (this->isParamValid("measurement_file"))
90 140 : readMeasurementsFromFile();
91 1268 : else if (this->isParamValid("measurement_points"))
92 287 : readMeasurementsFromInput();
93 :
94 774 : _misfit_values.resize(_measurement_values.size());
95 :
96 1548 : if (this->isParamValid("variable"))
97 : {
98 1155 : std::vector<VariableName> var_names(
99 : this->template getParam<std::vector<VariableName>>("variable"));
100 788 : for (const auto & name : var_names)
101 403 : _var_vec.push_back(&this->_fe_problem.getVariable(
102 403 : this->_tid, name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD));
103 385 : }
104 1548 : if (this->isParamValid("variable_weight_names"))
105 : {
106 99 : std::vector<std::string> weight_names(
107 : this->template getParam<std::vector<std::string>>("variable_weight_names"));
108 84 : for (const auto & name : weight_names)
109 : {
110 : if (_weight_names_weights_map.count(name) == 1)
111 : {
112 33 : _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 18 : _variable_weights.push_back(
119 54 : &this->template declareValueByName<std::vector<Real>>(name, REPORTER_MODE_REPLICATED));
120 18 : _variable_weights.back()->assign(_measurement_xcoord.size(), 1);
121 : }
122 : }
123 33 : }
124 2345 : 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 774 : }
132 :
133 : template <typename T>
134 : void
135 33893 : OptimizationDataTempl<T>::execute()
136 : {
137 33893 : computeMisfit();
138 33893 : _objective_val = computeMisfitValue();
139 33893 : }
140 :
141 : template <typename T>
142 : void
143 33893 : OptimizationDataTempl<T>::computeMisfit()
144 : {
145 33893 : 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 22011 : const std::size_t nvals = _measurement_values.size();
153 22011 : _simulation_values.resize(nvals, 0.0);
154 22011 : _misfit_values.resize(nvals);
155 :
156 22011 : errorCheckDataSize();
157 66065 : for (const auto var_index : make_range(_var_vec.size()))
158 : {
159 22027 : const auto & sys = _var_vec[var_index]->sys().system();
160 22027 : 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 180 : std::vector<Real> weights(_variable_weights.empty()
164 22027 : ? std::vector<Real>(_measurement_xcoord.size(), 1)
165 : : (*_variable_weights[var_index]));
166 875025 : for (const auto & i : make_range(nvals))
167 : {
168 852998 : 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 124112 : if (var_index == 0)
173 123888 : _simulation_values[i] = 0.0;
174 :
175 124112 : const Point point(_measurement_xcoord[i], _measurement_ycoord[i], _measurement_zcoord[i]);
176 124112 : const Real val = sys.point_value(vnum, point, false);
177 :
178 124112 : _simulation_values[i] += weights[i] * val;
179 124112 : _misfit_values[i] = _simulation_values[i] - _measurement_values[i];
180 : }
181 : }
182 : }
183 : }
184 :
185 : template <typename T>
186 : void
187 140 : OptimizationDataTempl<T>::readMeasurementsFromFile()
188 : {
189 280 : std::string xName = this->template getParam<std::string>("file_xcoord");
190 280 : std::string yName = this->template getParam<std::string>("file_ycoord");
191 280 : std::string zName = this->template getParam<std::string>("file_zcoord");
192 280 : std::string tName = this->template getParam<std::string>("file_time");
193 280 : std::string valueName = this->template getParam<std::string>("file_value");
194 280 : 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 280 : MooseUtils::DelimitedFileReader reader(this->template getParam<FileName>("measurement_file"));
204 140 : reader.read();
205 :
206 140 : auto const & names = reader.getNames();
207 140 : auto const & data = reader.getData();
208 :
209 : const std::size_t rows = data[0].size();
210 1090 : for (std::size_t i = 0; i < names.size(); ++i)
211 : {
212 : // make sure all data columns have the same length
213 950 : if (data[i].size() != rows)
214 0 : this->paramError("file", "Mismatching column lengths in file");
215 :
216 950 : if (names[i] == xName)
217 : {
218 140 : _measurement_xcoord = data[i];
219 : found_x = true;
220 : }
221 810 : else if (names[i] == yName)
222 : {
223 140 : _measurement_ycoord = data[i];
224 : found_y = true;
225 : }
226 670 : else if (names[i] == zName)
227 : {
228 140 : _measurement_zcoord = data[i];
229 : found_z = true;
230 : }
231 530 : else if (names[i] == tName)
232 : {
233 45 : _measurement_time = data[i];
234 : found_t = true;
235 : }
236 485 : else if (names[i] == valueName)
237 : {
238 140 : _measurement_values = data[i];
239 : found_value = true;
240 : }
241 345 : else if (std::find(weightNames.begin(), weightNames.end(), names[i]) != weightNames.end())
242 : {
243 39 : _weight_names_weights_map.emplace(names[i],
244 156 : &(this->template declareValueByName<std::vector<Real>>(
245 : names[i], REPORTER_MODE_REPLICATED)));
246 39 : _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 140 : if (!found_x)
252 0 : this->paramError(
253 : "measurement_file", "Column with name '", xName, "' missing from measurement file");
254 140 : if (!found_y)
255 0 : this->paramError(
256 : "measurement_file", "Column with name '", yName, "' missing from measurement file");
257 140 : if (!found_z)
258 0 : this->paramError(
259 : "measurement_file", "Column with name '", zName, "' missing from measurement file");
260 140 : if (!found_t)
261 95 : _measurement_time.assign(rows, 0);
262 140 : if (!found_value)
263 0 : this->paramError(
264 : "measurement_file", "Column with name '", valueName, "' missing from measurement file");
265 140 : 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 280 : }
279 :
280 : template <typename T>
281 : void
282 287 : OptimizationDataTempl<T>::readMeasurementsFromInput()
283 : {
284 574 : 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 2288 : for (const auto & p : this->template getParam<std::vector<Point>>("measurement_points"))
291 : {
292 1714 : _measurement_xcoord.push_back(p(0));
293 1714 : _measurement_ycoord.push_back(p(1));
294 1714 : _measurement_zcoord.push_back(p(2));
295 : }
296 :
297 574 : if (this->isParamValid("measurement_times"))
298 54 : _measurement_time = this->template getParam<std::vector<Real>>("measurement_times");
299 : else
300 269 : _measurement_time.assign(_measurement_xcoord.size(), 0.0);
301 :
302 574 : if (this->isParamValid("measurement_values"))
303 574 : _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 287 : }
307 :
308 : template <typename T>
309 : void
310 22011 : OptimizationDataTempl<T>::errorCheckDataSize()
311 : {
312 22011 : const std::size_t nvals = _measurement_values.size();
313 22011 : std::string msg = "";
314 22011 : if (_measurement_xcoord.size() != nvals)
315 0 : msg += "x-coordinate data (" + std::to_string(_measurement_xcoord.size()) + "), ";
316 22011 : if (_measurement_ycoord.size() != nvals)
317 0 : msg += "y-coordinate data (" + std::to_string(_measurement_ycoord.size()) + "), ";
318 22011 : if (_measurement_zcoord.size() != nvals)
319 0 : msg += "z-coordinate data (" + std::to_string(_measurement_zcoord.size()) + "), ";
320 22011 : if (_measurement_time.size() != nvals)
321 0 : msg += "time data (" + std::to_string(_measurement_time.size()) + "), ";
322 22011 : 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 22011 : }
329 :
330 : template <typename T>
331 : Real
332 33893 : OptimizationDataTempl<T>::computeMisfitValue()
333 : {
334 : Real val = 0.0;
335 1104875 : for (auto & misfit : _misfit_values)
336 1070982 : val += misfit * misfit;
337 :
338 33893 : return val * 0.5;
339 : }
340 :
341 : template class OptimizationDataTempl<GeneralReporter>;
342 : template class OptimizationDataTempl<OptimizationReporterBase>;
|
