libMesh
Public Member Functions | Public Attributes | Protected Attributes | List of all members
libMesh::RBParametrizedFunction Class Referenceabstract

A simple functor class that provides a RBParameter-dependent function. More...

#include <rb_parametrized_function.h>

Inheritance diagram for libMesh::RBParametrizedFunction:
[legend]

Public Member Functions

 RBParametrizedFunction ()
 Constructor. More...
 
 RBParametrizedFunction (RBParametrizedFunction &&)=default
 Special functions. More...
 
 RBParametrizedFunction (const RBParametrizedFunction &)=default
 
RBParametrizedFunctionoperator= (const RBParametrizedFunction &)=default
 
RBParametrizedFunctionoperator= (RBParametrizedFunction &&)=default
 
virtual ~RBParametrizedFunction ()
 
virtual unsigned int get_n_components () const =0
 Specify the number of components in this parametrized function. More...
 
virtual Number evaluate_comp (const RBParameters &mu, unsigned int comp, const Point &xyz, dof_id_type elem_id, unsigned int qp, subdomain_id_type subdomain_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
 Evaluate the parametrized function at the specified point for parameter mu. More...
 
virtual Number side_evaluate_comp (const RBParameters &mu, unsigned int comp, const Point &xyz, dof_id_type elem_id, unsigned int side_index, unsigned int qp, subdomain_id_type subdomain_id, boundary_id_type boundary_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
 Same as evaluate_comp() but for element sides. More...
 
virtual Number node_evaluate_comp (const RBParameters &mu, unsigned int comp, const Point &xyz, dof_id_type node_id, boundary_id_type boundary_id)
 Same as evaluate_comp() but for element nodes. More...
 
virtual std::vector< Numberevaluate (const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int qp, subdomain_id_type subdomain_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
 Evaluate the parametrized function at the specified point for parameter mu. More...
 
virtual std::vector< Numberside_evaluate (const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int side_index, unsigned int qp, subdomain_id_type subdomain_id, boundary_id_type boundary_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
 Same as evaluate() but for element sides. More...
 
virtual std::vector< Numbernode_evaluate (const RBParameters &mu, const Point &xyz, dof_id_type node_id, boundary_id_type boundary_id)
 Same as evaluate() but for element nodes. More...
 
virtual void vectorized_evaluate (const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
 Vectorized version of evaluate. More...
 
virtual void side_vectorized_evaluate (const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
 Same as vectorized_evaluate() but on element sides. More...
 
virtual void node_vectorized_evaluate (const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
 Same as vectorized_evaluate() but on element nodes. More...
 
virtual void preevaluate_parametrized_function_on_mesh (const RBParameters &mu, const std::unordered_map< dof_id_type, std::vector< Point >> &all_xyz, const std::unordered_map< dof_id_type, subdomain_id_type > &sbd_ids, const std::unordered_map< dof_id_type, std::vector< std::vector< Point >> > &all_xyz_perturb, const System &sys)
 Store the result of vectorized_evaluate. More...
 
virtual void preevaluate_parametrized_function_on_mesh_sides (const RBParameters &mu, const std::map< std::pair< dof_id_type, unsigned int >, std::vector< Point >> &side_all_xyz, const std::map< std::pair< dof_id_type, unsigned int >, subdomain_id_type > &sbd_ids, const std::map< std::pair< dof_id_type, unsigned int >, boundary_id_type > &side_boundary_ids, const std::map< std::pair< dof_id_type, unsigned int >, unsigned int > &side_types, const std::map< std::pair< dof_id_type, unsigned int >, std::vector< std::vector< Point >> > &side_all_xyz_perturb, const System &sys)
 Same as preevaluate_parametrized_function_on_mesh() except for mesh sides. More...
 
virtual void preevaluate_parametrized_function_on_mesh_nodes (const RBParameters &mu, const std::unordered_map< dof_id_type, Point > &all_xyz, const std::unordered_map< dof_id_type, boundary_id_type > &node_boundary_ids, const System &sys)
 Same as preevaluate_parametrized_function_on_mesh() except for mesh nodes. More...
 
virtual Number lookup_preevaluated_value_on_mesh (unsigned int comp, dof_id_type elem_id, unsigned int qp) const
 Look up the preevaluate values of the parametrized function for component comp, element elem_id, and quadrature point qp. More...
 
virtual Number lookup_preevaluated_side_value_on_mesh (unsigned int comp, dof_id_type elem_id, unsigned int side_index, unsigned int qp) const
 Look up the preevaluated values of the parametrized function for component comp, element elem_id, side_index, and quadrature point qp. More...
 
virtual Number lookup_preevaluated_node_value_on_mesh (unsigned int comp, dof_id_type node_id) const
 Look up the preevaluate values of the parametrized function for component comp, node node_id. More...
 
virtual void initialize_lookup_table ()
 If this parametrized function is defined based on a lookup table then we can call this function to initialize the table. More...
 
Number get_parameter_independent_data (const std::string &property_name, subdomain_id_type sbd_id) const
 Get the value stored in _parameter_independent_data associated with region_name and property_name. More...
 
const std::set< boundary_id_type > & get_parametrized_function_boundary_ids () const
 For RBParametrizedFunctions defined on element sides or nodes, we get/set the boundary IDs that this parametrized function is defined on. More...
 
void set_parametrized_function_boundary_ids (const std::set< boundary_id_type > &boundary_ids, bool is_nodal_boundary)
 
bool on_mesh_sides () const
 
bool on_mesh_nodes () const
 
virtual void get_spatial_indices (std::vector< std::vector< unsigned int >> &spatial_indices, const VectorizedEvalInput &v)
 In some cases a parametrized function is defined based on array data that we index into based on the spatial data from the mesh (e.g. More...
 
virtual void initialize_spatial_indices (const std::vector< std::vector< unsigned int >> &spatial_indices, const VectorizedEvalInput &v)
 The Online stage counterpart of get_spatial_indices(). More...
 
virtual void preevaluate_parametrized_function_cleanup ()
 Virtual function that performs cleanup after each "preevaluate parametrized function" evaluation. More...
 

Public Attributes

std::vector< std::vector< std::vector< Number > > > preevaluated_values
 Storage for pre-evaluated values. More...
 
std::unordered_map< dof_id_type, std::vector< unsigned int > > mesh_to_preevaluated_values_map
 Indexing into preevaluated_values for the case where the preevaluated values were obtained from evaluations at elements/quadrature points on a mesh. More...
 
std::map< std::pair< dof_id_type, unsigned int >, std::vector< unsigned int > > mesh_to_preevaluated_side_values_map
 Similar to the above except this map stores the data on element sides. More...
 
std::unordered_map< dof_id_type, unsigned intmesh_to_preevaluated_node_values_map
 Indexing into preevaluated_values for the case where the preevaluated values were obtained from evaluations at elements/quadrature points on a mesh. More...
 
bool requires_xyz_perturbations
 Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluate function values. More...
 
bool requires_all_elem_qp_data
 Boolean to indicate whether this parametrized function requires data from all qps on the current element at each qp location. More...
 
bool is_lookup_table
 Boolean to indicate if this parametrized function is defined based on a lookup table or not. More...
 
std::string lookup_table_param_name
 If this is a lookup table, then lookup_table_param_name specifies the parameter that is used to index into the lookup table. More...
 
Real fd_delta
 The finite difference step size in the case that this function in the case that this function uses finite differencing. More...
 

Protected Attributes

std::map< std::string, std::map< subdomain_id_type, Number > > _parameter_independent_data
 In some cases we need to store parameter-independent data which is related to this function but since it is parameter-indepedent should not be returned as part of evaluate(). More...
 
std::set< boundary_id_type_parametrized_function_boundary_ids
 In the case of an RBParametrizedFunction defined on element sides, this defines the set of boundary IDs that the function is defined on. More...
 
bool _is_nodal_boundary
 In the case that _parametrized_function_boundary_ids is not empty, then this parametrized function is defined on a mesh boundary. More...
 

Detailed Description

A simple functor class that provides a RBParameter-dependent function.

Author
David Knezevic
Date
2012 Provides a reduced basis parameterized function.

Definition at line 94 of file rb_parametrized_function.h.

Constructor & Destructor Documentation

◆ RBParametrizedFunction() [1/3]

libMesh::RBParametrizedFunction::RBParametrizedFunction ( )

Constructor.

Definition at line 47 of file rb_parametrized_function.C.

48 :
51 is_lookup_table(false),
52 fd_delta(1.e-6),
53 _is_nodal_boundary(false)
54 {}
bool is_lookup_table
Boolean to indicate if this parametrized function is defined based on a lookup table or not...
bool requires_all_elem_qp_data
Boolean to indicate whether this parametrized function requires data from all qps on the current elem...
bool requires_xyz_perturbations
Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluat...
bool _is_nodal_boundary
In the case that _parametrized_function_boundary_ids is not empty, then this parametrized function is...
Real fd_delta
The finite difference step size in the case that this function in the case that this function uses fi...

◆ RBParametrizedFunction() [2/3]

libMesh::RBParametrizedFunction::RBParametrizedFunction ( RBParametrizedFunction &&  )
default

Special functions.

  • This class can be default copy/move assigned/constructed.
  • The destructor is defaulted out-of-line.

◆ RBParametrizedFunction() [3/3]

libMesh::RBParametrizedFunction::RBParametrizedFunction ( const RBParametrizedFunction )
default

◆ ~RBParametrizedFunction()

libMesh::RBParametrizedFunction::~RBParametrizedFunction ( )
virtualdefault

Member Function Documentation

◆ evaluate()

std::vector< Number > libMesh::RBParametrizedFunction::evaluate ( const RBParameters mu,
const Point xyz,
dof_id_type  elem_id,
unsigned int  qp,
subdomain_id_type  subdomain_id,
const std::vector< Point > &  xyz_perturb,
const std::vector< Real > &  phi_i_qp 
)
virtual

Evaluate the parametrized function at the specified point for parameter mu.

If requires_xyz_perturbations==false, then xyz_perturb will not be used.

In this case we evaluate for all components.

Definition at line 109 of file rb_parametrized_function.C.

Referenced by evaluate_comp(), and vectorized_evaluate().

116 {
117  // This method should be overridden in subclasses, so we just give a not implemented error message here
118  libmesh_not_implemented();
119 
120  return std::vector<Number>();
121 }

◆ evaluate_comp()

Number libMesh::RBParametrizedFunction::evaluate_comp ( const RBParameters mu,
unsigned int  comp,
const Point xyz,
dof_id_type  elem_id,
unsigned int  qp,
subdomain_id_type  subdomain_id,
const std::vector< Point > &  xyz_perturb,
const std::vector< Real > &  phi_i_qp 
)
virtual

Evaluate the parametrized function at the specified point for parameter mu.

If requires_xyz_perturbations==false, then xyz_perturb will not be used.

In this case we return the value for component comp only, but the base class implementation simply calls the vector-returning evaluate() function below and returns the comp'th component, so derived classes should provide a more efficient routine or just call the vector-returning function instead.

Definition at line 59 of file rb_parametrized_function.C.

References evaluate().

67 {
68  std::vector<Number> values = evaluate(mu, xyz, elem_id, qp, subdomain_id, xyz_perturb, phi_i_qp);
69 
70  libmesh_error_msg_if(comp >= values.size(), "Error: Invalid value of comp");
71 
72  return values[comp];
73 }
virtual std::vector< Number > evaluate(const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int qp, subdomain_id_type subdomain_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
Evaluate the parametrized function at the specified point for parameter mu.

◆ get_n_components()

virtual unsigned int libMesh::RBParametrizedFunction::get_n_components ( ) const
pure virtual

Specify the number of components in this parametrized function.

A scalar-valued function has one component, a vector-valued function has more than one component.

Implemented in ShiftedGaussian, and Gxyz.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), node_vectorized_evaluate(), side_vectorized_evaluate(), and vectorized_evaluate().

◆ get_parameter_independent_data()

Number libMesh::RBParametrizedFunction::get_parameter_independent_data ( const std::string &  property_name,
subdomain_id_type  sbd_id 
) const

Get the value stored in _parameter_independent_data associated with region_name and property_name.

Definition at line 748 of file rb_parametrized_function.C.

References _parameter_independent_data.

750 {
751  return libmesh_map_find(libmesh_map_find(_parameter_independent_data, property_name), sbd_id);
752 }
std::map< std::string, std::map< subdomain_id_type, Number > > _parameter_independent_data
In some cases we need to store parameter-independent data which is related to this function but since...

◆ get_parametrized_function_boundary_ids()

const std::set< boundary_id_type > & libMesh::RBParametrizedFunction::get_parametrized_function_boundary_ids ( ) const

For RBParametrizedFunctions defined on element sides or nodes, we get/set the boundary IDs that this parametrized function is defined on.

Definition at line 771 of file rb_parametrized_function.C.

References _parametrized_function_boundary_ids.

Referenced by libMesh::RBEIMConstruction::initialize_qp_data(), libMesh::RBEIMEvaluation::node_distribute_bfs(), on_mesh_nodes(), on_mesh_sides(), and libMesh::RBEIMEvaluation::side_distribute_bfs().

772 {
774 }
std::set< boundary_id_type > _parametrized_function_boundary_ids
In the case of an RBParametrizedFunction defined on element sides, this defines the set of boundary I...

◆ get_spatial_indices()

void libMesh::RBParametrizedFunction::get_spatial_indices ( std::vector< std::vector< unsigned int >> &  spatial_indices,
const VectorizedEvalInput v 
)
virtual

In some cases a parametrized function is defined based on array data that we index into based on the spatial data from the mesh (e.g.

element, node, or side indices). We refer to the indices that we use to index into this array data as "spatial indices". This method sets spatial_indices based on the provided mesh-based indices.

Note that spatial_indices is defined as a doubly-nested vector so that we can handle the case where the spatial function evaluation requires us to have indices from all nodes of an element, since in that case we need a vector of indices (one per node) for each point. Other cases, such as when we define the parametrized function based on the element index only, only require a singly-nested vector which we handle as a special case of the doubly-nested vector.

This method is typically used in the Offline stage in order to generate and store the relevant spatial indices.

This method is a no-op by default, but it can be overridden in subclasses to provide the relevant behavior.

Definition at line 754 of file rb_parametrized_function.C.

Referenced by libMesh::RBEIMEvaluation::initialize_interpolation_points_spatial_indices().

756 {
757  // No-op by default
758 }

◆ initialize_lookup_table()

void libMesh::RBParametrizedFunction::initialize_lookup_table ( )
virtual

If this parametrized function is defined based on a lookup table then we can call this function to initialize the table.

This is a no-op by default, but it can be overridden in subclasses as needed.

Definition at line 743 of file rb_parametrized_function.C.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

744 {
745  // No-op by default, override in subclasses as needed
746 }

◆ initialize_spatial_indices()

void libMesh::RBParametrizedFunction::initialize_spatial_indices ( const std::vector< std::vector< unsigned int >> &  spatial_indices,
const VectorizedEvalInput v 
)
virtual

The Online stage counterpart of get_spatial_indices().

This method is used to initialize the spatial index data in this object so that we can evaluate it during an Online solve.

Definition at line 760 of file rb_parametrized_function.C.

Referenced by libMesh::RBEIMEvaluation::initialize_param_fn_spatial_indices().

762 {
763  // No-op by default
764 }

◆ lookup_preevaluated_node_value_on_mesh()

Number libMesh::RBParametrizedFunction::lookup_preevaluated_node_value_on_mesh ( unsigned int  comp,
dof_id_type  node_id 
) const
virtual

Look up the preevaluate values of the parametrized function for component comp, node node_id.

Definition at line 731 of file rb_parametrized_function.C.

References mesh_to_preevaluated_node_values_map, and preevaluated_values.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

733 {
734  unsigned int index =
735  libmesh_map_find(mesh_to_preevaluated_node_values_map, node_id);
736 
737  libmesh_error_msg_if(preevaluated_values.size() != 1, "Error: we expect only one parameter index");
738  libmesh_error_msg_if(index >= preevaluated_values[0].size(), "Error: invalid index");
739 
740  return preevaluated_values[0][index][comp];
741 }
std::unordered_map< dof_id_type, unsigned int > mesh_to_preevaluated_node_values_map
Indexing into preevaluated_values for the case where the preevaluated values were obtained from evalu...
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.

◆ lookup_preevaluated_side_value_on_mesh()

Number libMesh::RBParametrizedFunction::lookup_preevaluated_side_value_on_mesh ( unsigned int  comp,
dof_id_type  elem_id,
unsigned int  side_index,
unsigned int  qp 
) const
virtual

Look up the preevaluated values of the parametrized function for component comp, element elem_id, side_index, and quadrature point qp.

Definition at line 714 of file rb_parametrized_function.C.

References mesh_to_preevaluated_side_values_map, and preevaluated_values.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

718 {
719  const std::vector<unsigned int> & indices_at_qps =
720  libmesh_map_find(mesh_to_preevaluated_side_values_map, std::make_pair(elem_id,side_index));
721 
722  libmesh_error_msg_if(qp >= indices_at_qps.size(), "Error: invalid qp");
723 
724  unsigned int index = indices_at_qps[qp];
725  libmesh_error_msg_if(preevaluated_values.size() != 1, "Error: we expect only one parameter index");
726  libmesh_error_msg_if(index >= preevaluated_values[0].size(), "Error: invalid index");
727 
728  return preevaluated_values[0][index][comp];
729 }
std::map< std::pair< dof_id_type, unsigned int >, std::vector< unsigned int > > mesh_to_preevaluated_side_values_map
Similar to the above except this map stores the data on element sides.
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.

◆ lookup_preevaluated_value_on_mesh()

Number libMesh::RBParametrizedFunction::lookup_preevaluated_value_on_mesh ( unsigned int  comp,
dof_id_type  elem_id,
unsigned int  qp 
) const
virtual

Look up the preevaluate values of the parametrized function for component comp, element elem_id, and quadrature point qp.

Definition at line 698 of file rb_parametrized_function.C.

References mesh_to_preevaluated_values_map, and preevaluated_values.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

701 {
702  const std::vector<unsigned int> & indices_at_qps =
703  libmesh_map_find(mesh_to_preevaluated_values_map, elem_id);
704 
705  libmesh_error_msg_if(qp >= indices_at_qps.size(), "Error: invalid qp");
706 
707  unsigned int index = indices_at_qps[qp];
708  libmesh_error_msg_if(preevaluated_values.size() != 1, "Error: we expect only one parameter index");
709  libmesh_error_msg_if(index >= preevaluated_values[0].size(), "Error: invalid index");
710 
711  return preevaluated_values[0][index][comp];
712 }
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.
std::unordered_map< dof_id_type, std::vector< unsigned int > > mesh_to_preevaluated_values_map
Indexing into preevaluated_values for the case where the preevaluated values were obtained from evalu...

◆ node_evaluate()

std::vector< Number > libMesh::RBParametrizedFunction::node_evaluate ( const RBParameters mu,
const Point xyz,
dof_id_type  node_id,
boundary_id_type  boundary_id 
)
virtual

Same as evaluate() but for element nodes.

Definition at line 141 of file rb_parametrized_function.C.

Referenced by node_evaluate_comp(), and node_vectorized_evaluate().

145 {
146  // This method should be overridden in subclasses, so we just give a not implemented error message here
147  libmesh_not_implemented();
148 
149  return std::vector<Number>();
150 }

◆ node_evaluate_comp()

Number libMesh::RBParametrizedFunction::node_evaluate_comp ( const RBParameters mu,
unsigned int  comp,
const Point xyz,
dof_id_type  node_id,
boundary_id_type  boundary_id 
)
virtual

Same as evaluate_comp() but for element nodes.

Definition at line 95 of file rb_parametrized_function.C.

References node_evaluate().

100 {
101  std::vector<Number> values = node_evaluate(mu, xyz, node_id, boundary_id);
102 
103  libmesh_error_msg_if(comp >= values.size(), "Error: Invalid value of comp");
104 
105  return values[comp];
106 }
virtual std::vector< Number > node_evaluate(const RBParameters &mu, const Point &xyz, dof_id_type node_id, boundary_id_type boundary_id)
Same as evaluate() but for element nodes.

◆ node_vectorized_evaluate()

void libMesh::RBParametrizedFunction::node_vectorized_evaluate ( const std::vector< RBParameters > &  mus,
const VectorizedEvalInput v,
std::vector< std::vector< std::vector< Number >>> &  output 
)
virtual

Same as vectorized_evaluate() but on element nodes.

The base class implementation of this function loops over the input "mus" vector and calls node_evaluate() for each entry. The node_evaluate() function may be overridden in derived classes.

Definition at line 322 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::boundary_ids, get_n_components(), libMesh::index_range(), libMesh::make_range(), node_evaluate(), and libMesh::VectorizedEvalInput::node_ids.

Referenced by preevaluate_parametrized_function_on_mesh_nodes(), and libMesh::RBEIMEvaluation::rb_eim_solves().

325 {
326  LOG_SCOPE("node_vectorized_evaluate()", "RBParametrizedFunction");
327 
328  output.clear();
329  unsigned int n_points = v.all_xyz.size();
330 
331  // The number of components returned by this RBParametrizedFunction
332  auto n_components = this->get_n_components();
333 
334  // We first loop over all mus and all n_points, calling node_evaluate()
335  // for each and storing the results. It is easier to first
336  // pre-compute all the values before filling output, since, in the
337  // case of multi-sample RBParameters, the ordering of the loops is a
338  // bit complicated otherwise.
339  std::vector<std::vector<std::vector<Number>>> all_evals(mus.size());
340  for (auto mu_index : index_range(mus))
341  {
342  // Allocate enough space to store all points for the current mu
343  all_evals[mu_index].resize(n_points);
344  for (auto point_index : index_range(all_evals[mu_index]))
345  {
346  // Evaluate all samples for the current mu at the current interpolation point
347  all_evals[mu_index][point_index] =
348  this->node_evaluate(mus[mu_index],
349  v.all_xyz[point_index],
350  v.node_ids[point_index],
351  v.boundary_ids[point_index]);
352 
353  // The vector returned by node_evaluate() should contain:
354  // n_components * mus[mu_index].n_samples()
355  // entries. That is, for multi-sample RBParameters objects,
356  // the vector will be packed with entries as follows:
357  // [sample0_component0, sample0_component1, ..., sample0_componentN,
358  // sample1_component0, sample1_component1, ..., sample1_componentN,
359  // ...
360  // sampleM_component0, sampleM_component1, ..., sampleM_componentN]
361  auto n_samples = mus[mu_index].n_samples();
362  auto received_data = all_evals[mu_index][point_index].size();
363  libmesh_error_msg_if(received_data != n_components * n_samples,
364  "Recieved " << received_data <<
365  " evaluated values but expected to receive " << n_components * n_samples);
366  }
367  }
368 
369  // TODO: move this code for computing the total number of samples
370  // represented by a std::vector of RBParameters objects to a helper
371  // function.
372  unsigned int output_size = 0;
373  for (const auto & mu : mus)
374  output_size += mu.n_samples();
375 
376  output.resize(output_size);
377 
378  // We use traditional for-loops here (rather than range-based) so that we can declare and
379  // increment multiple loop counters all within the local scope of the for-loop.
380  for (auto [mu_index, output_index] = std::make_tuple(0u, 0u); mu_index < mus.size(); ++mu_index)
381  {
382  auto n_samples = mus[mu_index].n_samples();
383  for (auto mu_sample_idx = 0u; mu_sample_idx < n_samples; ++mu_sample_idx, ++output_index)
384  {
385  output[output_index].resize(n_points);
386  for (auto point_index : make_range(n_points))
387  {
388  output[output_index][point_index].resize(n_components);
389 
390  for (auto comp : make_range(n_components))
391  output[output_index][point_index][comp] = all_evals[mu_index][point_index][n_components*mu_sample_idx + comp];
392  }
393  }
394  }
395 }
virtual unsigned int get_n_components() const =0
Specify the number of components in this parametrized function.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
virtual std::vector< Number > node_evaluate(const RBParameters &mu, const Point &xyz, dof_id_type node_id, boundary_id_type boundary_id)
Same as evaluate() but for element nodes.

◆ on_mesh_nodes()

bool libMesh::RBParametrizedFunction::on_mesh_nodes ( ) const

◆ on_mesh_sides()

bool libMesh::RBParametrizedFunction::on_mesh_sides ( ) const

◆ operator=() [1/2]

RBParametrizedFunction& libMesh::RBParametrizedFunction::operator= ( const RBParametrizedFunction )
default

◆ operator=() [2/2]

RBParametrizedFunction& libMesh::RBParametrizedFunction::operator= ( RBParametrizedFunction &&  )
default

◆ preevaluate_parametrized_function_cleanup()

void libMesh::RBParametrizedFunction::preevaluate_parametrized_function_cleanup ( )
virtual

Virtual function that performs cleanup after each "preevaluate parametrized function" evaluation.

This function is a no-op by default, but it can be overridden in subclasses in order to do necessary cleanup, such as clearing cached data.

Definition at line 766 of file rb_parametrized_function.C.

Referenced by preevaluate_parametrized_function_on_mesh(), preevaluate_parametrized_function_on_mesh_nodes(), and preevaluate_parametrized_function_on_mesh_sides().

767 {
768  // No-op by default
769 }

◆ preevaluate_parametrized_function_on_mesh()

void libMesh::RBParametrizedFunction::preevaluate_parametrized_function_on_mesh ( const RBParameters mu,
const std::unordered_map< dof_id_type, std::vector< Point >> &  all_xyz,
const std::unordered_map< dof_id_type, subdomain_id_type > &  sbd_ids,
const std::unordered_map< dof_id_type, std::vector< std::vector< Point >> > &  all_xyz_perturb,
const System sys 
)
virtual

Store the result of vectorized_evaluate.

This is helpful during EIM training, since we can pre-evaluate and store the parameterized function for each training sample. If requires_xyz_perturbations==false, then all_xyz_perturb will not be used.

Definition at line 397 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::all_xyz_perturb, dim, libMesh::Elem::dim(), libMesh::FEMContext::elem_dimensions(), libMesh::VectorizedEvalInput::elem_ids, libMesh::MeshBase::elem_ref(), libMesh::VectorizedEvalInput::elem_types, libMesh::FEMContext::get_element_fe(), libMesh::System::get_mesh(), libMesh::if(), libMesh::index_range(), libMesh::VectorizedEvalInput::JxW_all_qp, mesh_to_preevaluated_values_map, libMesh::VectorizedEvalInput::phi_i_all_qp, libMesh::VectorizedEvalInput::phi_i_qp, libMesh::FEMContext::pre_fe_reinit(), preevaluate_parametrized_function_cleanup(), preevaluated_values, libMesh::VectorizedEvalInput::qps, requires_all_elem_qp_data, requires_xyz_perturbations, libMesh::VectorizedEvalInput::sbd_ids, libMesh::Elem::type(), and vectorized_evaluate().

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

402 {
404 
405  unsigned int n_points = 0;
406  for (const auto & xyz_pair : all_xyz)
407  {
408  const std::vector<Point> & xyz_vec = xyz_pair.second;
409  n_points += xyz_vec.size();
410  }
411 
412  VectorizedEvalInput v;
413  v.all_xyz.resize(n_points);
414  v.elem_ids.resize(n_points);
415  v.qps.resize(n_points);
416  v.sbd_ids.resize(n_points);
417  v.all_xyz_perturb.resize(n_points);
418  v.phi_i_qp.resize(n_points);
419  v.elem_types.resize(n_points);
420 
422  {
423  v.JxW_all_qp.resize(n_points);
424  v.phi_i_all_qp.resize(n_points);
425  }
426 
427  // Empty vector to be used when xyz perturbations are not required
428  std::vector<Point> empty_perturbs;
429 
430  // In order to compute phi_i_qp, we initialize a FEMContext
431  FEMContext con(sys);
432  for (auto dim : con.elem_dimensions())
433  {
434  auto fe = con.get_element_fe(/*var=*/0, dim);
435  fe->get_JxW();
436  fe->get_phi();
437  }
438 
439  unsigned int counter = 0;
440  for (const auto & [elem_id, xyz_vec] : all_xyz)
441  {
442  subdomain_id_type subdomain_id = libmesh_map_find(sbd_ids, elem_id);
443 
444  // The amount of data to be stored for each component
445  auto n_qp = xyz_vec.size();
446  mesh_to_preevaluated_values_map[elem_id].resize(n_qp);
447 
448  // Also initialize phi in order to compute phi_i_qp
449  const Elem & elem_ref = sys.get_mesh().elem_ref(elem_id);
450  con.pre_fe_reinit(sys, &elem_ref);
451 
452  auto elem_fe = con.get_element_fe(/*var=*/0, elem_ref.dim());
453  const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
454  const std::vector<Real> & JxW = elem_fe->get_JxW();
455 
456  elem_fe->reinit(&elem_ref);
457 
458  for (auto qp : index_range(xyz_vec))
459  {
460  mesh_to_preevaluated_values_map[elem_id][qp] = counter;
461 
462  v.all_xyz[counter] = xyz_vec[qp];
463  v.elem_ids[counter] = elem_id;
464  v.qps[counter] = qp;
465  v.sbd_ids[counter] = subdomain_id;
466  v.elem_types[counter] = elem_ref.type();
467 
468  v.phi_i_qp[counter].resize(phi.size());
469  for(auto i : index_range(phi))
470  v.phi_i_qp[counter][i] = phi[i][qp];
471 
473  {
474  const auto & qps_and_perturbs =
475  libmesh_map_find(all_xyz_perturb, elem_id);
476  libmesh_error_msg_if(qp >= qps_and_perturbs.size(), "Error: Invalid qp");
477 
478  v.all_xyz_perturb[counter] = qps_and_perturbs[qp];
479  }
480  else
481  {
482  v.all_xyz_perturb[counter] = empty_perturbs;
483  }
484 
486  {
487  // In this case we store data for all qps on this element
488  // at each point.
489  v.JxW_all_qp[counter].resize(JxW.size());
490  for(auto i : index_range(JxW))
491  v.JxW_all_qp[counter][i] = JxW[i];
492 
493  v.phi_i_all_qp[counter].resize(phi.size());
494  for(auto i : index_range(phi))
495  {
496  v.phi_i_all_qp[counter][i].resize(phi[i].size());
497  for(auto j : index_range(phi[i]))
498  v.phi_i_all_qp[counter][i][j] = phi[i][j];
499  }
500  }
501 
502  counter++;
503  }
504  }
505 
506  std::vector<RBParameters> mus {mu};
508 
510 }
unsigned int dim
TestClass subdomain_id_type
Based on the 4-byte comment warning above, this probably doesn&#39;t work with exodusII at all...
Definition: id_types.h:43
bool requires_all_elem_qp_data
Boolean to indicate whether this parametrized function requires data from all qps on the current elem...
virtual void preevaluate_parametrized_function_cleanup()
Virtual function that performs cleanup after each "preevaluate parametrized function" evaluation...
bool requires_xyz_perturbations
Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluat...
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.
std::unordered_map< dof_id_type, std::vector< unsigned int > > mesh_to_preevaluated_values_map
Indexing into preevaluated_values for the case where the preevaluated values were obtained from evalu...
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
virtual void vectorized_evaluate(const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
Vectorized version of evaluate.

◆ preevaluate_parametrized_function_on_mesh_nodes()

void libMesh::RBParametrizedFunction::preevaluate_parametrized_function_on_mesh_nodes ( const RBParameters mu,
const std::unordered_map< dof_id_type, Point > &  all_xyz,
const std::unordered_map< dof_id_type, boundary_id_type > &  node_boundary_ids,
const System sys 
)
virtual

Same as preevaluate_parametrized_function_on_mesh() except for mesh nodes.

Definition at line 661 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::boundary_ids, mesh_to_preevaluated_node_values_map, libMesh::VectorizedEvalInput::node_ids, node_vectorized_evaluate(), preevaluate_parametrized_function_cleanup(), and preevaluated_values.

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

665 {
667 
668  unsigned int n_points = all_xyz.size();
669 
670  VectorizedEvalInput v;
671  v.all_xyz.resize(n_points);
672  v.node_ids.resize(n_points);
673  v.boundary_ids.resize(n_points);
674 
675  // Empty vector to be used when xyz perturbations are not required
676  std::vector<Point> empty_perturbs;
677 
678  unsigned int counter = 0;
679  for (const auto & [node_id, p] : all_xyz)
680  {
681  boundary_id_type boundary_id = libmesh_map_find(node_boundary_ids, node_id);
682 
683  mesh_to_preevaluated_node_values_map[node_id] = counter;
684 
685  v.all_xyz[counter] = p;
686  v.node_ids[counter] = node_id;
687  v.boundary_ids[counter] = boundary_id;
688 
689  counter++;
690  }
691 
692  std::vector<RBParameters> mus {mu};
694 
696 }
std::unordered_map< dof_id_type, unsigned int > mesh_to_preevaluated_node_values_map
Indexing into preevaluated_values for the case where the preevaluated values were obtained from evalu...
int8_t boundary_id_type
Definition: id_types.h:51
virtual void preevaluate_parametrized_function_cleanup()
Virtual function that performs cleanup after each "preevaluate parametrized function" evaluation...
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.
virtual void node_vectorized_evaluate(const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
Same as vectorized_evaluate() but on element nodes.

◆ preevaluate_parametrized_function_on_mesh_sides()

void libMesh::RBParametrizedFunction::preevaluate_parametrized_function_on_mesh_sides ( const RBParameters mu,
const std::map< std::pair< dof_id_type, unsigned int >, std::vector< Point >> &  side_all_xyz,
const std::map< std::pair< dof_id_type, unsigned int >, subdomain_id_type > &  sbd_ids,
const std::map< std::pair< dof_id_type, unsigned int >, boundary_id_type > &  side_boundary_ids,
const std::map< std::pair< dof_id_type, unsigned int >, unsigned int > &  side_types,
const std::map< std::pair< dof_id_type, unsigned int >, std::vector< std::vector< Point >> > &  side_all_xyz_perturb,
const System sys 
)
virtual

Same as preevaluate_parametrized_function_on_mesh() except for mesh sides.

Definition at line 512 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::all_xyz_perturb, libMesh::VectorizedEvalInput::boundary_ids, libMesh::Elem::build_side_ptr(), dim, libMesh::Elem::dim(), libMesh::FEMContext::elem_dimensions(), libMesh::VectorizedEvalInput::elem_ids, libMesh::MeshBase::elem_ref(), libMesh::FEMContext::get_element_fe(), libMesh::System::get_mesh(), libMesh::FEMContext::get_side_fe(), libMesh::if(), libMesh::index_range(), mesh_to_preevaluated_side_values_map, libMesh::VectorizedEvalInput::phi_i_qp, libMesh::FEMContext::pre_fe_reinit(), preevaluate_parametrized_function_cleanup(), preevaluated_values, libMesh::VectorizedEvalInput::qps, requires_xyz_perturbations, libMesh::VectorizedEvalInput::sbd_ids, libMesh::VectorizedEvalInput::side_indices, and side_vectorized_evaluate().

Referenced by libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set().

519 {
521 
522  unsigned int n_points = 0;
523  for (const auto & xyz_pair : side_all_xyz)
524  {
525  const std::vector<Point> & xyz_vec = xyz_pair.second;
526  n_points += xyz_vec.size();
527  }
528 
529  VectorizedEvalInput v;
530  v.all_xyz.resize(n_points);
531  v.elem_ids.resize(n_points);
532  v.side_indices.resize(n_points);
533  v.qps.resize(n_points);
534  v.sbd_ids.resize(n_points);
535  v.boundary_ids.resize(n_points);
536  v.all_xyz_perturb.resize(n_points);
537  v.phi_i_qp.resize(n_points);
538 
539  // Empty vector to be used when xyz perturbations are not required
540  std::vector<Point> empty_perturbs;
541 
542  // In order to compute phi_i_qp, we initialize a FEMContext
543  FEMContext con(sys);
544  for (auto dim : con.elem_dimensions())
545  {
546  auto fe = con.get_element_fe(/*var=*/0, dim);
547  fe->get_phi();
548 
549  auto side_fe = con.get_side_fe(/*var=*/0, dim);
550  side_fe->get_phi();
551  }
552 
553  unsigned int counter = 0;
554  for (const auto & xyz_pair : side_all_xyz)
555  {
556  auto elem_side_pair = xyz_pair.first;
557  dof_id_type elem_id = elem_side_pair.first;
558  unsigned int side_index = elem_side_pair.second;
559 
560  const std::vector<Point> & xyz_vec = xyz_pair.second;
561 
562  subdomain_id_type subdomain_id = libmesh_map_find(sbd_ids, elem_side_pair);
563  boundary_id_type boundary_id = libmesh_map_find(side_boundary_ids, elem_side_pair);
564  unsigned int side_type = libmesh_map_find(side_types, elem_side_pair);
565 
566  // The amount of data to be stored for each component
567  auto n_qp = xyz_vec.size();
568  mesh_to_preevaluated_side_values_map[elem_side_pair].resize(n_qp);
569 
570  const Elem & elem_ref = sys.get_mesh().elem_ref(elem_id);
571  con.pre_fe_reinit(sys, &elem_ref);
572 
573  // side_type == 0 --> standard side
574  // side_type == 1 --> shellface
575  if (side_type == 0)
576  {
577  std::unique_ptr<const Elem> elem_side;
578  elem_ref.build_side_ptr(elem_side, side_index);
579 
580  auto side_fe = con.get_side_fe(/*var=*/0, elem_ref.dim());
581  side_fe->reinit(&elem_ref, side_index);
582 
583  const std::vector<std::vector<Real>> & phi = side_fe->get_phi();
584  for (auto qp : index_range(xyz_vec))
585  {
586  mesh_to_preevaluated_side_values_map[elem_side_pair][qp] = counter;
587 
588  v.all_xyz[counter] = xyz_vec[qp];
589  v.elem_ids[counter] = elem_side_pair.first;
590  v.side_indices[counter] = elem_side_pair.second;
591  v.qps[counter] = qp;
592  v.sbd_ids[counter] = subdomain_id;
593  v.boundary_ids[counter] = boundary_id;
594 
595  v.phi_i_qp[counter].resize(phi.size());
596  for(auto i : index_range(phi))
597  v.phi_i_qp[counter][i] = phi[i][qp];
598 
600  {
601  const auto & qps_and_perturbs =
602  libmesh_map_find(side_all_xyz_perturb, elem_side_pair);
603  libmesh_error_msg_if(qp >= qps_and_perturbs.size(), "Error: Invalid qp");
604 
605  v.all_xyz_perturb[counter] = qps_and_perturbs[qp];
606  }
607  else
608  {
609  v.all_xyz_perturb[counter] = empty_perturbs;
610  }
611  counter++;
612  }
613  }
614  else if (side_type == 1)
615  {
616  auto elem_fe = con.get_element_fe(/*var=*/0, elem_ref.dim());
617  const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
618 
619  elem_fe->reinit(&elem_ref);
620 
621  for (auto qp : index_range(xyz_vec))
622  {
623  mesh_to_preevaluated_side_values_map[elem_side_pair][qp] = counter;
624 
625  v.all_xyz[counter] = xyz_vec[qp];
626  v.elem_ids[counter] = elem_side_pair.first;
627  v.side_indices[counter] = elem_side_pair.second;
628  v.qps[counter] = qp;
629  v.sbd_ids[counter] = subdomain_id;
630  v.boundary_ids[counter] = boundary_id;
631 
632  v.phi_i_qp[counter].resize(phi.size());
633  for(auto i : index_range(phi))
634  v.phi_i_qp[counter][i] = phi[i][qp];
635 
637  {
638  const auto & qps_and_perturbs =
639  libmesh_map_find(side_all_xyz_perturb, elem_side_pair);
640  libmesh_error_msg_if(qp >= qps_and_perturbs.size(), "Error: Invalid qp");
641 
642  v.all_xyz_perturb[counter] = qps_and_perturbs[qp];
643  }
644  else
645  {
646  v.all_xyz_perturb[counter] = empty_perturbs;
647  }
648  counter++;
649  }
650  }
651  else
652  libmesh_error_msg ("Unrecognized side_type: " << side_type);
653  }
654 
655  std::vector<RBParameters> mus {mu};
657 
659 }
unsigned int dim
TestClass subdomain_id_type
Based on the 4-byte comment warning above, this probably doesn&#39;t work with exodusII at all...
Definition: id_types.h:43
virtual void side_vectorized_evaluate(const std::vector< RBParameters > &mus, const VectorizedEvalInput &v, std::vector< std::vector< std::vector< Number >>> &output)
Same as vectorized_evaluate() but on element sides.
int8_t boundary_id_type
Definition: id_types.h:51
virtual void preevaluate_parametrized_function_cleanup()
Virtual function that performs cleanup after each "preevaluate parametrized function" evaluation...
std::map< std::pair< dof_id_type, unsigned int >, std::vector< unsigned int > > mesh_to_preevaluated_side_values_map
Similar to the above except this map stores the data on element sides.
bool requires_xyz_perturbations
Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluat...
std::vector< std::vector< std::vector< Number > > > preevaluated_values
Storage for pre-evaluated values.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
uint8_t dof_id_type
Definition: id_types.h:67

◆ set_parametrized_function_boundary_ids()

void libMesh::RBParametrizedFunction::set_parametrized_function_boundary_ids ( const std::set< boundary_id_type > &  boundary_ids,
bool  is_nodal_boundary 
)

Definition at line 776 of file rb_parametrized_function.C.

References _is_nodal_boundary, and _parametrized_function_boundary_ids.

777 {
779  _is_nodal_boundary = is_nodal_boundary;
780 }
std::set< boundary_id_type > _parametrized_function_boundary_ids
In the case of an RBParametrizedFunction defined on element sides, this defines the set of boundary I...
bool _is_nodal_boundary
In the case that _parametrized_function_boundary_ids is not empty, then this parametrized function is...

◆ side_evaluate()

std::vector< Number > libMesh::RBParametrizedFunction::side_evaluate ( const RBParameters mu,
const Point xyz,
dof_id_type  elem_id,
unsigned int  side_index,
unsigned int  qp,
subdomain_id_type  subdomain_id,
boundary_id_type  boundary_id,
const std::vector< Point > &  xyz_perturb,
const std::vector< Real > &  phi_i_qp 
)
virtual

Same as evaluate() but for element sides.

Definition at line 124 of file rb_parametrized_function.C.

Referenced by side_evaluate_comp(), and side_vectorized_evaluate().

133 {
134  // This method should be overridden in subclasses, so we just give a not implemented error message here
135  libmesh_not_implemented();
136 
137  return std::vector<Number>();
138 }

◆ side_evaluate_comp()

Number libMesh::RBParametrizedFunction::side_evaluate_comp ( const RBParameters mu,
unsigned int  comp,
const Point xyz,
dof_id_type  elem_id,
unsigned int  side_index,
unsigned int  qp,
subdomain_id_type  subdomain_id,
boundary_id_type  boundary_id,
const std::vector< Point > &  xyz_perturb,
const std::vector< Real > &  phi_i_qp 
)
virtual

Same as evaluate_comp() but for element sides.

Definition at line 76 of file rb_parametrized_function.C.

References side_evaluate().

86 {
87  std::vector<Number> values = side_evaluate(mu, xyz, elem_id, side_index, qp, subdomain_id, boundary_id, xyz_perturb, phi_i_qp);
88 
89  libmesh_error_msg_if(comp >= values.size(), "Error: Invalid value of comp");
90 
91  return values[comp];
92 }
virtual std::vector< Number > side_evaluate(const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int side_index, unsigned int qp, subdomain_id_type subdomain_id, boundary_id_type boundary_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
Same as evaluate() but for element sides.

◆ side_vectorized_evaluate()

void libMesh::RBParametrizedFunction::side_vectorized_evaluate ( const std::vector< RBParameters > &  mus,
const VectorizedEvalInput v,
std::vector< std::vector< std::vector< Number >>> &  output 
)
virtual

Same as vectorized_evaluate() but on element sides.

The base class implementation of this function loops over the input "mus" vector and calls side_evaluate() for each entry. The side_evaluate() function may be overridden in derived classes.

Definition at line 236 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::all_xyz_perturb, libMesh::VectorizedEvalInput::boundary_ids, libMesh::VectorizedEvalInput::elem_ids, get_n_components(), libMesh::index_range(), libMesh::make_range(), libMesh::VectorizedEvalInput::phi_i_qp, libMesh::VectorizedEvalInput::qps, requires_xyz_perturbations, libMesh::VectorizedEvalInput::sbd_ids, side_evaluate(), and libMesh::VectorizedEvalInput::side_indices.

Referenced by preevaluate_parametrized_function_on_mesh_sides(), and libMesh::RBEIMEvaluation::rb_eim_solves().

239 {
240  LOG_SCOPE("side_vectorized_evaluate()", "RBParametrizedFunction");
241 
242  output.clear();
243  unsigned int n_points = v.all_xyz.size();
244 
245  libmesh_error_msg_if(v.sbd_ids.size() != n_points, "Error: invalid vector sizes");
246  libmesh_error_msg_if(requires_xyz_perturbations && (v.all_xyz_perturb.size() != n_points), "Error: invalid vector sizes");
247 
248  // Dummy vector to be used when xyz perturbations are not required
249  std::vector<Point> empty_perturbs;
250 
251  // The number of components returned by this RBParametrizedFunction
252  auto n_components = this->get_n_components();
253 
254  // We first loop over all mus and all n_points, calling side_evaluate()
255  // for each and storing the results. It is easier to first
256  // pre-compute all the values before filling output, since, in the
257  // case of multi-sample RBParameters, the ordering of the loops is a
258  // bit complicated otherwise.
259  std::vector<std::vector<std::vector<Number>>> all_evals(mus.size());
260  for (auto mu_index : index_range(mus))
261  {
262  // Allocate enough space to store all points for the current mu
263  all_evals[mu_index].resize(n_points);
264  for (auto point_index : index_range(all_evals[mu_index]))
265  {
266  // Evaluate all samples for the current mu at the current interpolation point
267  all_evals[mu_index][point_index] =
268  this->side_evaluate(mus[mu_index],
269  v.all_xyz[point_index],
270  v.elem_ids[point_index],
271  v.side_indices[point_index],
272  v.qps[point_index],
273  v.sbd_ids[point_index],
274  v.boundary_ids[point_index],
275  requires_xyz_perturbations ? v.all_xyz_perturb[point_index] : empty_perturbs,
276  v.phi_i_qp[point_index]);
277 
278  // The vector returned by side_evaluate() should contain:
279  // n_components * mus[mu_index].n_samples()
280  // entries. That is, for multi-sample RBParameters objects,
281  // the vector will be packed with entries as follows:
282  // [sample0_component0, sample0_component1, ..., sample0_componentN,
283  // sample1_component0, sample1_component1, ..., sample1_componentN,
284  // ...
285  // sampleM_component0, sampleM_component1, ..., sampleM_componentN]
286  auto n_samples = mus[mu_index].n_samples();
287  auto received_data = all_evals[mu_index][point_index].size();
288  libmesh_error_msg_if(received_data != n_components * n_samples,
289  "Recieved " << received_data <<
290  " evaluated values but expected to receive " << n_components * n_samples);
291  }
292  }
293 
294  // TODO: move this code for computing the total number of samples
295  // represented by a std::vector of RBParameters objects to a helper
296  // function.
297  unsigned int output_size = 0;
298  for (const auto & mu : mus)
299  output_size += mu.n_samples();
300 
301  output.resize(output_size);
302 
303  // We use traditional for-loops here (rather than range-based) so that we can declare and
304  // increment multiple loop counters all within the local scope of the for-loop.
305  for (auto [mu_index, output_index] = std::make_tuple(0u, 0u); mu_index < mus.size(); ++mu_index)
306  {
307  auto n_samples = mus[mu_index].n_samples();
308  for (auto mu_sample_idx = 0u; mu_sample_idx < n_samples; ++mu_sample_idx, ++output_index)
309  {
310  output[output_index].resize(n_points);
311  for (auto point_index : make_range(n_points))
312  {
313  output[output_index][point_index].resize(n_components);
314 
315  for (auto comp : make_range(n_components))
316  output[output_index][point_index][comp] = all_evals[mu_index][point_index][n_components*mu_sample_idx + comp];
317  }
318  }
319  }
320 }
virtual unsigned int get_n_components() const =0
Specify the number of components in this parametrized function.
virtual std::vector< Number > side_evaluate(const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int side_index, unsigned int qp, subdomain_id_type subdomain_id, boundary_id_type boundary_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
Same as evaluate() but for element sides.
bool requires_xyz_perturbations
Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluat...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ vectorized_evaluate()

void libMesh::RBParametrizedFunction::vectorized_evaluate ( const std::vector< RBParameters > &  mus,
const VectorizedEvalInput v,
std::vector< std::vector< std::vector< Number >>> &  output 
)
virtual

Vectorized version of evaluate.

If requires_xyz_perturbations==false, then all_xyz_perturb will not be used.

The base class implementation of this function loops over the input "mus" vector and calls evaluate() for each entry. The evaluate() function may be overridden in derived classes.

Definition at line 152 of file rb_parametrized_function.C.

References libMesh::VectorizedEvalInput::all_xyz, libMesh::VectorizedEvalInput::all_xyz_perturb, libMesh::VectorizedEvalInput::elem_ids, evaluate(), get_n_components(), libMesh::index_range(), libMesh::make_range(), libMesh::VectorizedEvalInput::phi_i_qp, libMesh::VectorizedEvalInput::qps, requires_xyz_perturbations, and libMesh::VectorizedEvalInput::sbd_ids.

Referenced by preevaluate_parametrized_function_on_mesh(), and libMesh::RBEIMEvaluation::rb_eim_solves().

155 {
156  LOG_SCOPE("vectorized_evaluate()", "RBParametrizedFunction");
157 
158  output.clear();
159  unsigned int n_points = v.all_xyz.size();
160 
161  libmesh_error_msg_if(v.sbd_ids.size() != n_points, "Error: invalid vector sizes");
162  libmesh_error_msg_if(requires_xyz_perturbations && (v.all_xyz_perturb.size() != n_points), "Error: invalid vector sizes");
163 
164  // Dummy vector to be used when xyz perturbations are not required
165  std::vector<Point> empty_perturbs;
166 
167  // The number of components returned by this RBParametrizedFunction
168  auto n_components = this->get_n_components();
169 
170  // We first loop over all mus and all n_points, calling evaluate()
171  // for each and storing the results. It is easier to first
172  // pre-compute all the values before filling output, since, in the
173  // case of multi-sample RBParameters, the ordering of the loops is a
174  // bit complicated otherwise.
175  std::vector<std::vector<std::vector<Number>>> all_evals(mus.size());
176  for (auto mu_index : index_range(mus))
177  {
178  // Allocate enough space to store all points for the current mu
179  all_evals[mu_index].resize(n_points);
180  for (auto point_index : index_range(all_evals[mu_index]))
181  {
182  // Evaluate all samples for the current mu at the current interpolation point
183  all_evals[mu_index][point_index] =
184  this->evaluate(mus[mu_index],
185  v.all_xyz[point_index],
186  v.elem_ids[point_index],
187  v.qps[point_index],
188  v.sbd_ids[point_index],
189  requires_xyz_perturbations ? v.all_xyz_perturb[point_index] : empty_perturbs,
190  v.phi_i_qp[point_index]);
191 
192  // The vector returned by evaluate() should contain:
193  // n_components * mus[mu_index].n_samples()
194  // entries. That is, for multi-sample RBParameters objects,
195  // the vector will be packed with entries as follows:
196  // [sample0_component0, sample0_component1, ..., sample0_componentN,
197  // sample1_component0, sample1_component1, ..., sample1_componentN,
198  // ...
199  // sampleM_component0, sampleM_component1, ..., sampleM_componentN]
200  auto n_samples = mus[mu_index].n_samples();
201  auto received_data = all_evals[mu_index][point_index].size();
202  libmesh_error_msg_if(received_data != n_components * n_samples,
203  "Recieved " << received_data <<
204  " evaluated values but expected to receive " << n_components * n_samples);
205  }
206  }
207 
208  // TODO: move this code for computing the total number of samples
209  // represented by a std::vector of RBParameters objects to a helper
210  // function.
211  unsigned int output_size = 0;
212  for (const auto & mu : mus)
213  output_size += mu.n_samples();
214 
215  output.resize(output_size);
216 
217  // We use traditional for-loops here (rather than range-based) so that we can declare and
218  // increment multiple loop counters all within the local scope of the for-loop.
219  for (auto [mu_index, output_index] = std::make_tuple(0u, 0u); mu_index < mus.size(); ++mu_index)
220  {
221  auto n_samples = mus[mu_index].n_samples();
222  for (auto mu_sample_idx = 0u; mu_sample_idx < n_samples; ++mu_sample_idx, ++output_index)
223  {
224  output[output_index].resize(n_points);
225  for (auto point_index : make_range(n_points))
226  {
227  output[output_index][point_index].resize(n_components);
228 
229  for (auto comp : make_range(n_components))
230  output[output_index][point_index][comp] = all_evals[mu_index][point_index][n_components*mu_sample_idx + comp];
231  }
232  }
233  }
234 }
virtual unsigned int get_n_components() const =0
Specify the number of components in this parametrized function.
bool requires_xyz_perturbations
Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluat...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
virtual std::vector< Number > evaluate(const RBParameters &mu, const Point &xyz, dof_id_type elem_id, unsigned int qp, subdomain_id_type subdomain_id, const std::vector< Point > &xyz_perturb, const std::vector< Real > &phi_i_qp)
Evaluate the parametrized function at the specified point for parameter mu.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

Member Data Documentation

◆ _is_nodal_boundary

bool libMesh::RBParametrizedFunction::_is_nodal_boundary
protected

In the case that _parametrized_function_boundary_ids is not empty, then this parametrized function is defined on a mesh boundary.

This boolean indicates if the mesh boundary under consideration is a set of sides, or a set of nodes.

Definition at line 446 of file rb_parametrized_function.h.

Referenced by on_mesh_nodes(), on_mesh_sides(), and set_parametrized_function_boundary_ids().

◆ _parameter_independent_data

std::map<std::string, std::map<subdomain_id_type, Number> > libMesh::RBParametrizedFunction::_parameter_independent_data
protected

In some cases we need to store parameter-independent data which is related to this function but since it is parameter-indepedent should not be returned as part of evaluate().

We index this data by "property name" –> subdomain_id –> value.

Definition at line 433 of file rb_parametrized_function.h.

Referenced by get_parameter_independent_data().

◆ _parametrized_function_boundary_ids

std::set<boundary_id_type> libMesh::RBParametrizedFunction::_parametrized_function_boundary_ids
protected

In the case of an RBParametrizedFunction defined on element sides, this defines the set of boundary IDs that the function is defined on.

Definition at line 439 of file rb_parametrized_function.h.

Referenced by get_parametrized_function_boundary_ids(), and set_parametrized_function_boundary_ids().

◆ fd_delta

Real libMesh::RBParametrizedFunction::fd_delta

The finite difference step size in the case that this function in the case that this function uses finite differencing.

Definition at line 422 of file rb_parametrized_function.h.

Referenced by libMesh::RBEIMConstruction::initialize_qp_data().

◆ is_lookup_table

bool libMesh::RBParametrizedFunction::is_lookup_table

Boolean to indicate if this parametrized function is defined based on a lookup table or not.

If it is defined based on a lookup table, then the evaluation functions will access a discrete parameter to determine the index to lookup.

Definition at line 410 of file rb_parametrized_function.h.

Referenced by libMesh::RBDataSerialization::add_rb_eim_evaluation_data_to_builder(), libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), libMesh::RBDataDeserialization::load_rb_eim_evaluation_data(), and libMesh::RBEIMConstruction::train_eim_approximation_with_greedy().

◆ lookup_table_param_name

std::string libMesh::RBParametrizedFunction::lookup_table_param_name

If this is a lookup table, then lookup_table_param_name specifies the parameter that is used to index into the lookup table.

Definition at line 416 of file rb_parametrized_function.h.

Referenced by libMesh::RBEIMConstruction::set_rb_construction_parameters().

◆ mesh_to_preevaluated_node_values_map

std::unordered_map<dof_id_type, unsigned int> libMesh::RBParametrizedFunction::mesh_to_preevaluated_node_values_map

Indexing into preevaluated_values for the case where the preevaluated values were obtained from evaluations at elements/quadrature points on a mesh.

The indexing here is: node_id –> point_index

Definition at line 387 of file rb_parametrized_function.h.

Referenced by lookup_preevaluated_node_value_on_mesh(), and preevaluate_parametrized_function_on_mesh_nodes().

◆ mesh_to_preevaluated_side_values_map

std::map<std::pair<dof_id_type,unsigned int>, std::vector<unsigned int> > libMesh::RBParametrizedFunction::mesh_to_preevaluated_side_values_map

Similar to the above except this map stores the data on element sides.

The indexing here is: (elem_id,side index) –> qp –> point_index

Definition at line 379 of file rb_parametrized_function.h.

Referenced by lookup_preevaluated_side_value_on_mesh(), and preevaluate_parametrized_function_on_mesh_sides().

◆ mesh_to_preevaluated_values_map

std::unordered_map<dof_id_type, std::vector<unsigned int> > libMesh::RBParametrizedFunction::mesh_to_preevaluated_values_map

Indexing into preevaluated_values for the case where the preevaluated values were obtained from evaluations at elements/quadrature points on a mesh.

The indexing here is: elem_id –> qp –> point_index Then preevaluated_values[0][point_index] provides the vector of component values at that point.

Definition at line 372 of file rb_parametrized_function.h.

Referenced by lookup_preevaluated_value_on_mesh(), and preevaluate_parametrized_function_on_mesh().

◆ preevaluated_values

std::vector<std::vector<std::vector<Number> > > libMesh::RBParametrizedFunction::preevaluated_values

◆ requires_all_elem_qp_data

bool libMesh::RBParametrizedFunction::requires_all_elem_qp_data

Boolean to indicate whether this parametrized function requires data from all qps on the current element at each qp location.

This can be necessary in certain cases, e.g. when the parametrized function depends on "element average" quantities.

Definition at line 403 of file rb_parametrized_function.h.

Referenced by libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), and preevaluate_parametrized_function_on_mesh().

◆ requires_xyz_perturbations

bool libMesh::RBParametrizedFunction::requires_xyz_perturbations

Boolean to indicate whether this parametrized function requires xyz perturbations in order to evaluate function values.

An example of where perturbations are required is when the parametrized function is based on finite difference approximations to derivatives.

Definition at line 395 of file rb_parametrized_function.h.

Referenced by preevaluate_parametrized_function_on_mesh(), preevaluate_parametrized_function_on_mesh_sides(), side_vectorized_evaluate(), and vectorized_evaluate().


The documentation for this class was generated from the following files: