20 #ifndef LIBMESH_RB_EIM_EVALUATION_H 21 #define LIBMESH_RB_EIM_EVALUATION_H 24 #include "libmesh/point.h" 25 #include "libmesh/rb_theta_expansion.h" 26 #include "libmesh/rb_parametrized.h" 27 #include "libmesh/parallel_object.h" 28 #include "libmesh/dense_matrix.h" 29 #include "libmesh/dense_vector.h" 30 #include "libmesh/rb_parametrized_function.h" 42 class RBParametrizedFunction;
77 typedef std::map<dof_id_type, std::vector<std::vector<Number>>>
QpDataMap;
82 typedef std::map<std::pair<dof_id_type,unsigned int>, std::vector<std::vector<Number>>>
SideQpDataMap;
87 typedef std::map<dof_id_type, std::vector<Number>>
NodeDataMap;
92 virtual void clear()
override;
128 void rb_eim_solves(
const std::vector<RBParameters> & mus,
unsigned int N);
208 std::vector<Number> & values);
216 unsigned int side_index,
218 std::vector<Number> & values);
248 unsigned int side_index,
274 std::vector<Number> & values)
const;
281 unsigned int side_index,
283 std::vector<Number> & values)
const;
292 unsigned int var)
const;
301 unsigned int qp)
const;
308 unsigned int side_index,
310 unsigned int qp)
const;
320 unsigned int var)
const;
437 const std::vector<Point> & perturbs,
438 const std::vector<Real> & phi_i_qp,
440 const std::vector<Real> & JxW_all_qp,
441 const std::vector<std::vector<Real>> & phi_i_all_qp);
456 unsigned int side_index,
460 const std::vector<Point> & perturbs,
461 const std::vector<Real> & phi_i_qp);
499 bool write_binary_basis_functions =
true);
512 const std::string & directory_name =
"offline_data",
513 bool read_binary_basis_functions =
true);
532 const std::string & directory_name =
"offline_data");
574 Number error_indicator_rhs,
611 bool write_binary_basis_functions);
618 bool write_binary_basis_functions);
625 bool write_binary_basis_functions);
632 const std::string & directory_name,
633 bool read_binary_basis_functions);
640 const std::string & directory_name,
641 bool read_binary_basis_functions);
648 const std::string & directory_name,
649 bool read_binary_basis_functions);
819 #endif // LIBMESH_RB_EIM_EVALUATION_H void set_rb_eim_solutions(const std::vector< DenseVector< Number >> &rb_eim_solutions)
Set _rb_eim_solutions.
virtual void clear() override
Clear this object.
VectorizedEvalInput _vec_eval_input
We store the EIM interpolation point data in this object.
ElemType
Defines an enum for geometric element types.
void add_interpolation_points_boundary_id(boundary_id_type b_id)
void write_out_node_basis_functions(const std::string &directory_name, bool write_binary_basis_functions)
Method that writes out element node EIM basis functions.
std::map< std::pair< dof_id_type, unsigned int >, std::vector< std::vector< Number > > > SideQpDataMap
Type of the data structure used to map from (elem id, side index) -> [n_vars][n_qp] data...
static Number get_parametrized_function_node_local_value(const NodeDataMap &pf, dof_id_type node_id, unsigned int comp)
Same as get_parametrized_function_values_at_qps() except for node data.
const SideQpDataMap & get_side_basis_function(unsigned int i) const
Get a reference to the i^th side basis function.
void set_eim_error_indicator_active(bool is_active)
Activate/decative the error indicator in EIM solves.
const DenseMatrix< Number > & get_interpolation_matrix() const
Get the EIM interpolation matrix.
void add_basis_function(const QpDataMap &bf)
Add bf to our EIM basis.
void set_interpolation_matrix_entry(unsigned int i, unsigned int j, Number value)
Set entry of the EIM interpolation matrix.
std::vector< RBParameters > _rb_eim_solves_mus
The parameters and the number of basis functions that were used in the most recent call to rb_eim_sol...
std::map< dof_id_type, std::vector< Number > > NodeDataMap
Type of the data structure used to map from (node id) -> [n_vars] data.
subdomain_id_type get_interpolation_points_subdomain_id(unsigned int index) const
void distribute_bfs(const System &sys)
Helper function that distributes the entries of _local_eim_basis_functions to their respective proces...
void get_eim_basis_function_side_values_at_qps(unsigned int basis_function_index, dof_id_type elem_id, unsigned int side_index, unsigned int var, std::vector< Number > &values) const
Same as get_eim_basis_function_values_at_qps() except for side data.
const std::set< unsigned int > & scale_components_in_enrichment() const
Get _scale_components_in_enrichment.
void get_eim_basis_function_values_at_qps(unsigned int basis_function_index, dof_id_type elem_id, unsigned int var, std::vector< Number > &values) const
Fill up values with the basis function values for basis function basis_function_index and variable va...
void rb_eim_solves(const std::vector< RBParameters > &mus, unsigned int N)
Perform rb_eim_solves at each mu in mus and store the results in _rb_eim_solutions.
std::set< unsigned int > _eim_vars_to_project_and_write
This set specifies which EIM variables will be projected and written out in write_out_projected_basis...
Number get_eim_basis_function_value(unsigned int basis_function_index, dof_id_type elem_id, unsigned int comp, unsigned int qp) const
Same as above, except that we just return the value at the qp^th quadrature point.
A simple functor class that provides a RBParameter-dependent function.
unsigned int get_n_basis_functions() const
Return the current number of EIM basis functions.
void add_node_basis_function(const NodeDataMap &node_bf)
Add node_bf to our EIM basis.
RBEIMEvaluation(const Parallel::Communicator &comm)
Constructor.
unsigned int get_n_interpolation_points_spatial_indices() const
_interpolation_points_spatial_indices is optional data, so we need to be able to check how many _inte...
void add_interpolation_points_elem_id(dof_id_type elem_id)
std::vector< Number > get_rb_eim_solutions_entries(unsigned int index) const
Return entry index for each solution in _rb_eim_solutions.
RBEIMEvaluation & operator=(const RBEIMEvaluation &)=delete
void read_in_node_basis_functions(const System &sys, const std::string &directory_name, bool read_binary_basis_functions)
Method that reads in element node EIM basis functions.
void add_interpolation_data(Point p, unsigned int comp, dof_id_type elem_id, subdomain_id_type subdomain_id, unsigned int qp, const std::vector< Point > &perturbs, const std::vector< Real > &phi_i_qp, ElemType elem_type, const std::vector< Real > &JxW_all_qp, const std::vector< std::vector< Real >> &phi_i_all_qp)
Add interpolation data associated with a new basis function.
const std::vector< std::vector< Real > > & get_interpolation_points_phi_i_all_qp(unsigned int index) const
const std::vector< DenseVector< Number > > & get_rb_eim_solutions() const
Return the EIM solution coefficients from the most recent call to rb_eim_solves().
std::unique_ptr< RBParametrizedFunction > _parametrized_function
Store the parametrized function that will be approximated by this EIM system.
const Parallel::Communicator & comm() const
const std::vector< Real > & get_rb_eim_error_indicators() const
Return the EIM error indicator values from the most recent call to rb_eim_solves().
const std::vector< Point > & get_interpolation_points_xyz_perturbations(unsigned int index) const
unsigned int get_interpolation_points_side_index(unsigned int index) const
bool get_preserve_rb_eim_solutions() const
Get _preserve_rb_eim_solutions.
unsigned int get_n_interpolation_points() const
Return the number of interpolation points.
void decrement_vector(QpDataMap &v, const DenseVector< Number > &coeffs)
Subtract coeffs[i]*basis_function[i] from v.
const NodeDataMap & get_node_basis_function(unsigned int i) const
Get a reference to the i^th node basis function.
The libMesh namespace provides an interface to certain functionality in the library.
unsigned int get_interpolation_points_comp(unsigned int index) const
const std::vector< DenseVector< Number > > & get_eim_solutions_for_training_set() const
Return a const reference to the EIM solutions for the parameters in the training set.
void initialize_param_fn_spatial_indices()
The Online counterpart of initialize_interpolation_points_spatial_indices().
void set_parametrized_function(std::unique_ptr< RBParametrizedFunction > pf)
Set the parametrized function that we will approximate using the Empirical Interpolation Method...
void add_interpolation_points_xyz_perturbations(const std::vector< Point > &perturbs)
bool _is_eim_error_indicator_active
Indicate if the EIM error indicator is active in RB EIM solves.
void set_preserve_rb_eim_solutions(bool preserve_rb_eim_solutions)
Set _preserve_rb_eim_solutions.
std::vector< SideQpDataMap > _local_side_eim_basis_functions
The EIM basis functions on element sides.
void initialize_eim_theta_objects()
Build a vector of RBTheta objects that accesses the components of the RB_solution member variable of ...
void add_side_basis_function(const SideQpDataMap &side_bf)
Add side_bf to our EIM basis.
static Number get_parametrized_function_value(const Parallel::Communicator &comm, const QpDataMap &pf, dof_id_type elem_id, unsigned int comp, unsigned int qp)
Same as above, except that we just return the value at the qp^th quadrature point.
void side_decrement_vector(SideQpDataMap &v, const DenseVector< Number > &coeffs)
Same as decrement_vector() except for Side data.
void side_distribute_bfs(const System &sys)
Same as distribute_bfs() except for side data.
Number get_eim_basis_function_node_local_value(unsigned int basis_function_index, dof_id_type node_id, unsigned int var) const
Same as get_eim_basis_function_values_at_qps() except for node data.
std::vector< std::unique_ptr< RBTheta > > & get_eim_theta_objects()
const std::vector< Real > & get_interpolation_points_JxW_all_qp(unsigned int index) const
std::vector< std::unique_ptr< RBTheta > > _rb_eim_theta_objects
The vector of RBTheta objects that are created to point to this RBEIMEvaluation.
static void get_parametrized_function_values_at_qps(const QpDataMap &pf, dof_id_type elem_id, unsigned int comp, std::vector< Number > &values)
Fill up values by evaluating the parametrized function pf for all quadrature points on element elem_i...
void side_gather_bfs()
Same as gather_bfs() except for side data.
std::vector< NodeDataMap > _local_node_eim_basis_functions
The EIM basis functions on element nodes (e.g.
const QpDataMap & get_basis_function(unsigned int i) const
Get a reference to the i^th basis function.
void read_in_basis_functions(const System &sys, const std::string &directory_name="offline_data", bool read_binary_basis_functions=true)
Read in all the basis functions from file.
boundary_id_type get_interpolation_points_boundary_id(unsigned int index) const
void print_local_eim_basis_functions() const
Print the contents of _local_eim_basis_functions to libMesh::out.
void resize_data_structures(const unsigned int Nmax)
Resize the data structures for storing data associated with this object.
std::vector< std::vector< unsigned int > > _interpolation_points_spatial_indices
Here we store the spatial indices that were initialized by initialize_spatial_indices_at_interp_pts()...
void set_n_basis_functions(unsigned int n_bfs)
Set the number of basis functions.
Manages consistently variables, degrees of freedom, and coefficient vectors.
void node_distribute_bfs(const System &sys)
Same as distribute_bfs() except for node data.
void read_in_interior_basis_functions(const System &sys, const std::string &directory_name, bool read_binary_basis_functions)
Method that reads in element interior EIM basis functions.
void write_out_interior_basis_functions(const std::string &directory_name, bool write_binary_basis_functions)
Method that writes out element interior EIM basis functions.
std::map< dof_id_type, std::vector< std::vector< Number > > > QpDataMap
Type of the data structure used to map from (elem id) -> [n_vars][n_qp] data.
Number get_eim_basis_function_side_value(unsigned int basis_function_index, dof_id_type elem_id, unsigned int side_index, unsigned int comp, unsigned int qp) const
Same as get_eim_basis_function_value() except for side data.
const std::vector< unsigned int > & get_interpolation_points_spatial_indices(unsigned int index) const
const DenseVector< Number > & get_error_indicator_interpolation_row() const
Get/set _extra_points_interpolation_matrix.
ElemType get_interpolation_points_elem_type(unsigned int index) const
void add_side_interpolation_data(Point p, unsigned int comp, dof_id_type elem_id, unsigned int side_index, subdomain_id_type subdomain_id, boundary_id_type boundary_id, unsigned int qp, const std::vector< Point > &perturbs, const std::vector< Real > &phi_i_qp)
Add interpolation data associated with a new basis function.
void add_interpolation_points_side_index(unsigned int side_index)
void node_gather_bfs()
Same as gather_bfs() except for node data.
void add_node_interpolation_data(Point p, unsigned int comp, dof_id_type node_id, boundary_id_type boundary_id)
Add interpolation data associated with a new basis function.
void add_interpolation_points_phi_i_all_qp(const std::vector< std::vector< Real >> &phi_i_all_qp)
void write_out_projected_basis_functions(System &sys, const std::string &directory_name="offline_data")
Project all basis functions using project_qp_data_map_onto_system() and then write out the resulting ...
void node_decrement_vector(NodeDataMap &v, const DenseVector< Number > &coeffs)
Same as decrement_vector() except for node data.
void add_interpolation_points_phi_i_qp(const std::vector< Real > &phi_i_qp)
const std::set< unsigned int > & get_eim_vars_to_project_and_write() const
Get _eim_vars_to_project_and_write.
virtual ~RBEIMEvaluation()
bool _preserve_rb_eim_solutions
Boolean to indicate if we skip updating _rb_eim_solutions in rb_eim_solves().
std::vector< DenseVector< Number > > _eim_solutions_for_training_set
Storage for EIM solutions from the training set.
An object whose state is distributed along a set of processors.
static void get_parametrized_function_side_values_at_qps(const SideQpDataMap &pf, dof_id_type elem_id, unsigned int side_index, unsigned int comp, std::vector< Number > &values)
Same as get_parametrized_function_values_at_qps() except for side data.
virtual std::unique_ptr< RBTheta > build_eim_theta(unsigned int index)
Build a theta object corresponding to EIM index index.
Point get_interpolation_points_xyz(unsigned int index) const
Get the data associated with EIM interpolation points.
dof_id_type get_interpolation_points_node_id(unsigned int index) const
virtual bool use_eim_error_indicator() const
Virtual function to indicate if we use the EIM error indicator in this case.
void add_interpolation_points_subdomain_id(subdomain_id_type sbd_id)
std::vector< Real > _rb_eim_error_indicators
If we're using the EIM error indicator, then we store the error indicator values corresponding to _rb...
DenseVector< Number > rb_eim_solve(DenseVector< Number > &EIM_rhs)
Calculate the EIM approximation for the given right-hand side vector EIM_rhs.
void add_interpolation_points_node_id(dof_id_type node_id)
void add_interpolation_points_xyz(Point p)
Set the data associated with EIM interpolation points.
unsigned int _rb_eim_solves_N
std::pair< Real, Real > get_eim_error_indicator(Number error_indicator_rhs, const DenseVector< Number > &eim_solution, const DenseVector< Number > &eim_rhs)
Evaluates the EIM error indicator based on error_indicator_rhs, eim_solution, and _error_indicator_in...
void read_in_side_basis_functions(const System &sys, const std::string &directory_name, bool read_binary_basis_functions)
Method that reads in element side EIM basis functions.
Number get_eim_basis_function_node_value(unsigned int basis_function_index, dof_id_type node_id, unsigned int var) const
Same as get_eim_basis_function_value() except for node data.
const VectorizedEvalInput & get_vec_eval_input() const
Get the VectorizedEvalInput data.
void add_interpolation_points_JxW_all_qp(const std::vector< Real > &JxW_all_qp)
static Number get_parametrized_function_side_value(const Parallel::Communicator &comm, const SideQpDataMap &pf, dof_id_type elem_id, unsigned int side_index, unsigned int comp, unsigned int qp)
Same as get_parametrized_function_value() except for side data.
void project_qp_data_map_onto_system(System &sys, const QpDataMap &bf_data, unsigned int var)
Project variable var of bf_data into the solution vector of System.
This class is part of the rbOOmit framework.
std::vector< DenseVector< Number > > _rb_eim_solutions
The EIM solution coefficients from the most recent call to rb_eim_solves().
void set_error_indicator_interpolation_row(const DenseVector< Number > &error_indicator_row)
std::vector< unsigned int > _interpolation_points_comp
In the case of a "vector-valued" EIM, this vector determines which component of the parameterized fun...
std::set< unsigned int > _scale_components_in_enrichment
This set that specifies which EIM variables will be scaled during EIM enrichment so that their maximu...
dof_id_type get_interpolation_points_elem_id(unsigned int index) const
void add_interpolation_points_comp(unsigned int comp)
void add_interpolation_points_qp(unsigned int qp)
std::vector< QpDataMap > _local_eim_basis_functions
The EIM basis functions.
RBParametrizedFunction & get_parametrized_function()
Get a reference to the parametrized function.
void add_interpolation_points_spatial_indices(const std::vector< unsigned int > &spatial_indices)
const std::vector< Real > & get_interpolation_points_phi_i_qp(unsigned int index) const
This class enables evaluation of an Empirical Interpolation Method (EIM) approximation.
void write_out_basis_functions(const std::string &directory_name="offline_data", bool write_binary_basis_functions=true)
Write out all the basis functions to file.
void initialize_interpolation_points_spatial_indices()
Initialize _interpolation_points_spatial_indices.
A Point defines a location in LIBMESH_DIM dimensional Real space.
unsigned int get_interpolation_points_qp(unsigned int index) const
void write_out_side_basis_functions(const std::string &directory_name, bool write_binary_basis_functions)
Method that writes out element side EIM basis functions.
DenseVector< Number > _error_indicator_interpolation_row
Here we store an extra row of the interpolation matrix which is used to compute the EIM error indicat...
void add_interpolation_points_elem_type(ElemType elem_type)
void gather_bfs()
Helper function that gathers the contents of _local_eim_basis_functions to processor 0 in preparation...
DenseMatrix< Number > _interpolation_matrix
Dense matrix that stores the lower triangular interpolation matrix that can be used.
static Number get_parametrized_function_node_value(const Parallel::Communicator &comm, const NodeDataMap &pf, dof_id_type node_id, unsigned int comp)
Same as get_parametrized_function_value() except for node data.