Line data Source code
1 : // rbOOmit: An implementation of the Certified Reduced Basis method.
2 : // Copyright (C) 2009, 2010 David J. Knezevic
3 :
4 : // This file is part of rbOOmit.
5 :
6 : // rbOOmit is free software; you can redistribute it and/or
7 : // modify it under the terms of the GNU Lesser General Public
8 : // License as published by the Free Software Foundation; either
9 : // version 2.1 of the License, or (at your option) any later version.
10 :
11 : // rbOOmit is distributed in the hope that it will be useful,
12 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 : // Lesser General Public License for more details.
15 :
16 : // You should have received a copy of the GNU Lesser General Public
17 : // License along with this library; if not, write to the Free Software
18 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 :
20 : #ifndef LIBMESH_RB_PARAMETRIZED_FUNCTION_H
21 : #define LIBMESH_RB_PARAMETRIZED_FUNCTION_H
22 :
23 : // libMesh includes
24 : #include "libmesh/libmesh_common.h"
25 : #include "libmesh/enum_elem_type.h"
26 : #include "libmesh/enum_order.h"
27 :
28 : // C++ includes
29 : #include <unordered_map>
30 : #include <vector>
31 : #include <map>
32 : #include <set>
33 :
34 : namespace libMesh
35 : {
36 :
37 : class RBParameters;
38 : class Point;
39 : class System;
40 :
41 : namespace Parallel {
42 : class Communicator;
43 : }
44 :
45 : /**
46 : * Define a struct for the input to the "vectorized evaluate" functions below.
47 : * This encapsulates the arguments into a class to prevent having many function
48 : * arguments, and also makes it easier to make API changes in the future because
49 : * we can change these structs without changing the function arguments.
50 : */
51 : struct VectorizedEvalInput
52 : {
53 0 : VectorizedEvalInput () = default;
54 : VectorizedEvalInput (const VectorizedEvalInput &) = default;
55 : VectorizedEvalInput & operator= (const VectorizedEvalInput &) = default;
56 : VectorizedEvalInput (VectorizedEvalInput &&) = default;
57 : VectorizedEvalInput & operator= (VectorizedEvalInput &&) = default;
58 0 : virtual ~VectorizedEvalInput() = default;
59 :
60 : /**
61 : * Clear all the members.
62 : */
63 : void clear();
64 :
65 : /**
66 : * The members that define the inputs to the vectorized evaluate functions. Note
67 : * that some of these members may be unused, for example when we call the "interior"
68 : * vectorized evaluate function, we do not use node_ids.
69 : *
70 : * Some data below is the same for all points within an element, e.g. when we store
71 : * data at multiple qps per element the sbd_ids, elem_ids, JxW_all_qp, phi_i_all_qp
72 : * will store the same data repeated n_qp times per element. A possible optimization
73 : * for this would be to store this data based on element indices rather than qp
74 : * indices, and store an index per qp to index into the element-based data vector.
75 : * This optimization hasn't been implemented at this stage, but it could be added
76 : * at some point later.
77 : */
78 : std::vector<Point> all_xyz;
79 : std::vector<dof_id_type> elem_ids;
80 : std::vector<unsigned int> qps;
81 : std::vector<subdomain_id_type> sbd_ids;
82 : std::vector<std::vector<Point>> all_xyz_perturb;
83 : std::vector<std::vector<Real>> phi_i_qp;
84 : std::vector<unsigned int> side_indices;
85 : std::vector<boundary_id_type> boundary_ids;
86 : std::vector<dof_id_type> node_ids;
87 : std::vector<ElemType> elem_types;
88 : /**
89 : * The following containers are indexed by element id to avoid duplicated data.
90 : * The elements have a local indexing as elem_ids might not always be contiguous.
91 : */
92 : std::map<dof_id_type, unsigned int> elem_id_to_local_index;
93 : std::vector<std::vector<Real>> JxW_all_qp;
94 : std::vector<std::vector<std::vector<Real>>> phi_i_all_qp;
95 : std::vector<Point> dxyzdxi_elem_center;
96 : std::vector<Point> dxyzdeta_elem_center;
97 : std::vector<Order> qrule_orders;
98 :
99 : /**
100 : * Generic map that can be used to store any list of ids (elements, nodes, elemsets, subdomains...)
101 : * corresponding to a property (string) that can be used to apply a specific treatment to a group
102 : * of entities.
103 : */
104 : std::unordered_map<std::string, std::set<dof_id_type>> rb_property_map;
105 : };
106 :
107 : /**
108 : * A simple functor class that provides a RBParameter-dependent function.
109 : *
110 : * \author David Knezevic
111 : * \date 2012
112 : * \brief Provides a reduced basis parameterized function.
113 : */
114 0 : class RBParametrizedFunction
115 : {
116 : public:
117 :
118 : /**
119 : * Constructor.
120 : */
121 : RBParametrizedFunction();
122 :
123 : /**
124 : * Special functions.
125 : * - This class can be default copy/move assigned/constructed.
126 : * - The destructor is defaulted out-of-line.
127 : */
128 : RBParametrizedFunction (RBParametrizedFunction &&) = default;
129 : RBParametrizedFunction (const RBParametrizedFunction &) = default;
130 : RBParametrizedFunction & operator= (const RBParametrizedFunction &) = default;
131 : RBParametrizedFunction & operator= (RBParametrizedFunction &&) = default;
132 : virtual ~RBParametrizedFunction();
133 :
134 : /**
135 : * Specify the number of components in this parametrized function.
136 : * A scalar-valued function has one component, a vector-valued
137 : * function has more than one component.
138 : */
139 : virtual unsigned int get_n_components() const = 0;
140 :
141 : /**
142 : * Evaluate the parametrized function at the specified point for
143 : * parameter \p mu. If requires_xyz_perturbations==false, then
144 : * xyz_perturb will not be used.
145 : *
146 : * In this case we return the value for component \p comp only, but
147 : * the base class implementation simply calls the vector-returning
148 : * evaluate() function below and returns the comp'th component, so
149 : * derived classes should provide a more efficient routine or just call
150 : * the vector-returning function instead.
151 : */
152 : virtual Number evaluate_comp(const RBParameters & mu,
153 : unsigned int comp,
154 : const Point & xyz,
155 : dof_id_type elem_id,
156 : unsigned int qp,
157 : subdomain_id_type subdomain_id,
158 : const std::vector<Point> & xyz_perturb,
159 : const std::vector<Real> & phi_i_qp);
160 :
161 : /**
162 : * Same as evaluate_comp() but for element sides.
163 : */
164 : virtual Number side_evaluate_comp(const RBParameters & mu,
165 : unsigned int comp,
166 : const Point & xyz,
167 : dof_id_type elem_id,
168 : unsigned int side_index,
169 : unsigned int qp,
170 : subdomain_id_type subdomain_id,
171 : boundary_id_type boundary_id,
172 : const std::vector<Point> & xyz_perturb,
173 : const std::vector<Real> & phi_i_qp);
174 :
175 : /**
176 : * Same as evaluate_comp() but for element nodes.
177 : */
178 : virtual Number node_evaluate_comp(const RBParameters & mu,
179 : unsigned int comp,
180 : const Point & xyz,
181 : dof_id_type node_id,
182 : boundary_id_type boundary_id);
183 :
184 : /**
185 : * Evaluate the parametrized function at the specified point for
186 : * parameter \p mu. If requires_xyz_perturbations==false, then
187 : * xyz_perturb will not be used.
188 : *
189 : * In this case we evaluate for all components.
190 : */
191 : virtual std::vector<Number> evaluate(const RBParameters & mu,
192 : const Point & xyz,
193 : dof_id_type elem_id,
194 : unsigned int qp,
195 : subdomain_id_type subdomain_id,
196 : const std::vector<Point> & xyz_perturb,
197 : const std::vector<Real> & phi_i_qp);
198 :
199 : /**
200 : * Same as evaluate() but for element sides.
201 : */
202 : virtual std::vector<Number> side_evaluate(const RBParameters & mu,
203 : const Point & xyz,
204 : dof_id_type elem_id,
205 : unsigned int side_index,
206 : unsigned int qp,
207 : subdomain_id_type subdomain_id,
208 : boundary_id_type boundary_id,
209 : const std::vector<Point> & xyz_perturb,
210 : const std::vector<Real> & phi_i_qp);
211 :
212 : /**
213 : * Same as evaluate() but for element nodes.
214 : */
215 : virtual std::vector<Number> node_evaluate(const RBParameters & mu,
216 : const Point & xyz,
217 : dof_id_type node_id,
218 : boundary_id_type boundary_id);
219 :
220 : /**
221 : * Vectorized version of evaluate. If requires_xyz_perturbations==false, then all_xyz_perturb will not be used.
222 : *
223 : * The base class implementation of this function loops over the
224 : * input "mus" vector and calls evaluate() for each entry. The
225 : * evaluate() function may be overridden in derived classes.
226 : */
227 : virtual void vectorized_evaluate(const std::vector<RBParameters> & mus,
228 : const VectorizedEvalInput & v,
229 : std::vector<std::vector<std::vector<Number>>> & output);
230 :
231 : /**
232 : * Same as vectorized_evaluate() but on element sides.
233 : *
234 : * The base class implementation of this function loops over the
235 : * input "mus" vector and calls side_evaluate() for each entry. The
236 : * side_evaluate() function may be overridden in derived classes.
237 : */
238 : virtual void side_vectorized_evaluate(const std::vector<RBParameters> & mus,
239 : const VectorizedEvalInput & v,
240 : std::vector<std::vector<std::vector<Number>>> & output);
241 :
242 : /**
243 : * Same as vectorized_evaluate() but on element nodes.
244 : *
245 : * The base class implementation of this function loops over the
246 : * input "mus" vector and calls node_evaluate() for each entry. The
247 : * node_evaluate() function may be overridden in derived classes.
248 : */
249 : virtual void node_vectorized_evaluate(const std::vector<RBParameters> & mus,
250 : const VectorizedEvalInput & v,
251 : std::vector<std::vector<std::vector<Number>>> & output);
252 :
253 : /**
254 : * Store the result of vectorized_evaluate. This is helpful during EIM training,
255 : * since we can pre-evaluate and store the parameterized function for each training
256 : * sample. If requires_xyz_perturbations==false, then all_xyz_perturb will not be used.
257 : */
258 : virtual void preevaluate_parametrized_function_on_mesh(const RBParameters & mu,
259 : const std::unordered_map<dof_id_type, std::vector<Point>> & all_xyz,
260 : const std::unordered_map<dof_id_type, subdomain_id_type> & sbd_ids,
261 : const std::unordered_map<dof_id_type, std::vector<std::vector<Point>> > & all_xyz_perturb,
262 : const System & sys);
263 :
264 : /**
265 : * Same as preevaluate_parametrized_function_on_mesh() except for mesh sides.
266 : */
267 : virtual void preevaluate_parametrized_function_on_mesh_sides(const RBParameters & mu,
268 : const std::map<std::pair<dof_id_type,unsigned int>, std::vector<Point>> & side_all_xyz,
269 : const std::map<std::pair<dof_id_type,unsigned int>, subdomain_id_type> & sbd_ids,
270 : const std::map<std::pair<dof_id_type,unsigned int>, boundary_id_type> & side_boundary_ids,
271 : const std::map<std::pair<dof_id_type,unsigned int>, unsigned int> & side_types,
272 : const std::map<std::pair<dof_id_type,unsigned int>, std::vector<std::vector<Point>> > & side_all_xyz_perturb,
273 : const System & sys);
274 :
275 : /**
276 : * Same as preevaluate_parametrized_function_on_mesh() except for mesh nodes.
277 : */
278 : virtual void preevaluate_parametrized_function_on_mesh_nodes(const RBParameters & mu,
279 : const std::unordered_map<dof_id_type, Point> & all_xyz,
280 : const std::unordered_map<dof_id_type, boundary_id_type> & node_boundary_ids,
281 : const System & sys);
282 :
283 : /**
284 : * Look up the preevaluate values of the parametrized function for
285 : * component \p comp, element \p elem_id, and quadrature point \p qp.
286 : */
287 : virtual Number lookup_preevaluated_value_on_mesh(unsigned int comp,
288 : dof_id_type elem_id,
289 : unsigned int qp) const;
290 :
291 : /**
292 : * Look up the preevaluated values of the parametrized function for
293 : * component \p comp, element \p elem_id, \p side_index, and quadrature point \p qp.
294 : */
295 : virtual Number lookup_preevaluated_side_value_on_mesh(unsigned int comp,
296 : dof_id_type elem_id,
297 : unsigned int side_index,
298 : unsigned int qp) const;
299 :
300 : /**
301 : * Look up the preevaluate values of the parametrized function for
302 : * component \p comp, node \p node_id.
303 : */
304 : virtual Number lookup_preevaluated_node_value_on_mesh(unsigned int comp,
305 : dof_id_type node_id) const;
306 :
307 : /**
308 : * If this parametrized function is defined based on a lookup table then
309 : * we can call this function to initialize the table. This is a no-op by
310 : * default, but it can be overridden in subclasses as needed.
311 : */
312 : virtual void initialize_lookup_table();
313 :
314 : /**
315 : * Get the value stored in _parameter_independent_data associated with
316 : * \p region_name and \p property_name.
317 : */
318 : Number get_parameter_independent_data(const std::string & property_name,
319 : subdomain_id_type sbd_id) const;
320 :
321 : /**
322 : * For RBParametrizedFunctions defined on element sides or nodes, we get/set the boundary
323 : * IDs that this parametrized function is defined on.
324 : */
325 : const std::set<boundary_id_type> & get_parametrized_function_boundary_ids() const;
326 : void set_parametrized_function_boundary_ids(const std::set<boundary_id_type> & boundary_ids, bool is_nodal_boundary);
327 :
328 : /**
329 : * @return true if this parametrized function is defined on mesh sides.
330 : */
331 : bool on_mesh_sides() const;
332 :
333 : /**
334 : * @return true if this parametrized function is defined on mesh nodes.
335 : */
336 : bool on_mesh_nodes() const;
337 :
338 : /**
339 : * In some cases a parametrized function is defined based on array data that
340 : * we index into based on the spatial data from the mesh (e.g. element, node,
341 : * or side indices). We refer to the indices that we use to index into this
342 : * array data as "spatial indices". This method sets \p spatial_indices based
343 : * on the provided mesh-based indices.
344 : *
345 : * Note that \p spatial_indices is defined as a doubly-nested vector so that we
346 : * can handle the case where the spatial function evaluation requires us to have
347 : * indices from all nodes of an element, since in that case we need a vector of
348 : * indices (one per node) for each point. Other cases, such as when we define
349 : * the parametrized function based on the element index only, only require a
350 : * singly-nested vector which we handle as a special case of the doubly-nested
351 : * vector.
352 : *
353 : * This method is typically used in the Offline stage in order to generate
354 : * and store the relevant spatial indices.
355 : *
356 : * This method is a no-op by default, but it can be overridden in subclasses
357 : * to provide the relevant behavior.
358 : */
359 : virtual void get_spatial_indices(std::vector<std::vector<unsigned int>> & spatial_indices,
360 : const VectorizedEvalInput & v);
361 :
362 : /**
363 : * The Online stage counterpart of get_spatial_indices(). This method
364 : * is used to initialize the spatial index data in this object so that
365 : * we can evaluate it during an Online solve.
366 : */
367 : virtual void initialize_spatial_indices(const std::vector<std::vector<unsigned int>> & spatial_indices,
368 : const VectorizedEvalInput & v);
369 :
370 : /**
371 : * Virtual function that performs cleanup after each
372 : * "preevaluate parametrized function" evaluation. This function
373 : * is a no-op by default, but it can be overridden in subclasses
374 : * in order to do necessary cleanup, such as clearing cached data.
375 : */
376 : virtual void preevaluate_parametrized_function_cleanup();
377 :
378 : /**
379 : * Function that returns a reference to the rb_property_map stored in the RBParametrizedFunction.
380 : */
381 : const std::unordered_map<std::string, std::set<dof_id_type>> & get_rb_property_map() const;
382 :
383 : /**
384 : * Function that adds a property to the RBParametrizedFunction rb_property_map.
385 : * The function checks that there is no duplicated properties before adding the property to the map.
386 : */
387 : void add_rb_property_map_entry(std::string & property_name, std::set<dof_id_type> & entity_ids);
388 :
389 : /**
390 : * Virtual function that can be overridden in RBParametrizedFunction subclasses to store
391 : * properties in the rb_property_map.
392 : */
393 : virtual void add_interpolation_data_to_rb_property_map(
394 : const Parallel::Communicator & comm,
395 : std::unordered_map<std::string, std::set<dof_id_type>> & rb_property_map,
396 : dof_id_type elem_id);
397 :
398 : /**
399 : * Storage for pre-evaluated values. The indexing is given by:
400 : * parameter index --> point index --> component index --> value.
401 : */
402 : std::vector<std::vector<std::vector<Number>>> preevaluated_values;
403 :
404 : /**
405 : * Indexing into preevaluated_values for the case where the preevaluated values
406 : * were obtained from evaluations at elements/quadrature points on a mesh.
407 : * The indexing here is:
408 : * elem_id --> qp --> point_index
409 : * Then preevaluated_values[0][point_index] provides the vector of component values at
410 : * that point.
411 : */
412 : std::unordered_map<dof_id_type, std::vector<unsigned int>> mesh_to_preevaluated_values_map;
413 :
414 : /**
415 : * Similar to the above except this map stores the data on element sides.
416 : * The indexing here is:
417 : * (elem_id,side index) --> qp --> point_index
418 : */
419 : std::map<std::pair<dof_id_type,unsigned int>, std::vector<unsigned int>> mesh_to_preevaluated_side_values_map;
420 :
421 : /**
422 : * Indexing into preevaluated_values for the case where the preevaluated values
423 : * were obtained from evaluations at elements/quadrature points on a mesh.
424 : * The indexing here is:
425 : * node_id --> point_index
426 : */
427 : std::unordered_map<dof_id_type, unsigned int> mesh_to_preevaluated_node_values_map;
428 :
429 : /**
430 : * Boolean to indicate whether this parametrized function requires xyz perturbations
431 : * in order to evaluate function values. An example of where perturbations are
432 : * required is when the parametrized function is based on finite difference
433 : * approximations to derivatives.
434 : */
435 : bool requires_xyz_perturbations;
436 :
437 : /**
438 : * Boolean to indicate whether this parametrized function requires data from
439 : * all qps on the current element at each qp location. This can be necessary
440 : * in certain cases, e.g. when the parametrized function depends on "element
441 : * average" quantities.
442 : */
443 : bool requires_all_elem_qp_data;
444 :
445 : /**
446 : * Boolean to indicate whether this parametrized function requires data
447 : * from the center on the current element.
448 : */
449 : bool requires_all_elem_center_data;
450 :
451 : /**
452 : * Boolean to indicate if this parametrized function is defined based on a lookup
453 : * table or not. If it is defined based on a lookup table, then the evaluation
454 : * functions will access a discrete parameter to determine the index to lookup.
455 : */
456 : bool is_lookup_table;
457 :
458 : /**
459 : * If this is a lookup table, then lookup_table_param_name specifies the parameter
460 : * that is used to index into the lookup table.
461 : */
462 : std::string lookup_table_param_name;
463 :
464 : /**
465 : * The finite difference step size in the case that this function in the case
466 : * that this function uses finite differencing.
467 : */
468 : Real fd_delta;
469 :
470 : protected:
471 :
472 : /**
473 : * In some cases we need to store parameter-independent data which is related
474 : * to this function but since it is parameter-indepedent should not be returned
475 : * as part of evaluate().
476 : *
477 : * We index this data by "property name" --> subdomain_id --> value.
478 : */
479 : std::map<std::string, std::map<subdomain_id_type, Number>> _parameter_independent_data;
480 :
481 : /**
482 : * In the case of an RBParametrizedFunction defined on element sides, this defines
483 : * the set of boundary IDs that the function is defined on.
484 : */
485 : std::set<boundary_id_type> _parametrized_function_boundary_ids;
486 :
487 : /**
488 : * In the case that _parametrized_function_boundary_ids is not empty, then this
489 : * parametrized function is defined on a mesh boundary. This boolean indicates
490 : * if the mesh boundary under consideration is a set of sides, or a set of nodes.
491 : */
492 : bool _is_nodal_boundary;
493 :
494 : /**
495 : * Generic property map used to store data during mesh pre-evaluation. It should be initialized
496 : * from a vectorized_evaluate() call during the mesh pre-evaluation. Then this RBParametrizedFunction
497 : * rb_property_map is used to fill the RBEIMEvaluation VectorizedEvalInput rb_property_map with the same properties but
498 : * interpolation point data only instead of the entire mesh.
499 : */
500 : std::unordered_map<std::string, std::set<dof_id_type>> _rb_property_map;
501 : };
502 :
503 : }
504 :
505 : #endif // LIBMESH_RB_PARAMETRIZED_FUNCTION_H
|
