libMesh
rb_construction.h
Go to the documentation of this file.
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_CONSTRUCTION_H
21 #define LIBMESH_RB_CONSTRUCTION_H
22 
23 // rbOOmit includes
24 #include "libmesh/rb_construction_base.h"
25 
26 // libMesh includes
27 #include "libmesh/linear_implicit_system.h"
28 #include "libmesh/dense_vector.h"
29 #include "libmesh/dense_matrix.h"
30 #include "libmesh/dg_fem_context.h"
31 #include "libmesh/dirichlet_boundaries.h"
32 
33 // C++ includes
34 
35 namespace libMesh
36 {
37 
38 class RBThetaExpansion;
39 class RBAssemblyExpansion;
40 class RBEvaluation;
41 class ElemAssembly;
42 
53 class RBConstruction : public RBConstructionBase<LinearImplicitSystem>
54 {
55 public:
56 
62  const std::string & name,
63  const unsigned int number);
64 
70  RBConstruction (RBConstruction &&) = default;
71  RBConstruction (const RBConstruction &) = delete;
72  RBConstruction & operator= (const RBConstruction &) = delete;
74  virtual ~RBConstruction ();
75 
80 
85  virtual void solve_for_matrix_and_rhs (LinearSolver<Number> & input_solver,
86  SparseMatrix<Number> & input_matrix,
87  NumericVector<Number> & input_rhs);
88 
92  void set_rb_evaluation(RBEvaluation & rb_eval_in);
93 
98  const RBEvaluation & get_rb_evaluation() const;
99 
103  bool is_rb_eval_initialized() const;
104 
111 
115  void set_rb_assembly_expansion(RBAssemblyExpansion & rb_assembly_expansion_in);
116 
121 
125  sys_type & system () { return *this; }
126 
131 
136  virtual void clear () override;
137 
141  virtual std::string system_type () const override;
142 
149  virtual Real truth_solve(int plot_solution);
150 
158  virtual Real train_reduced_basis(const bool resize_rb_eval_data=true);
159 
176  Real train_reduced_basis_with_greedy(const bool resize_rb_eval_data);
177 
183  void enrich_basis_from_rhs_terms(const bool resize_rb_eval_data=true);
184 
201 
207  virtual Real compute_max_error_bound();
208 
213  const RBParameters & get_greedy_parameter(unsigned int i);
214 
218  void set_rel_training_tolerance(Real new_training_tolerance)
219  {this->rel_training_tolerance = new_training_tolerance; }
221 
225  void set_abs_training_tolerance(Real new_training_tolerance)
226  {this->abs_training_tolerance = new_training_tolerance; }
228 
232  void set_normalize_rb_bound_in_greedy(bool normalize_rb_bound_in_greedy_in)
233  {this->normalize_rb_bound_in_greedy = normalize_rb_bound_in_greedy_in; }
235 
239  void set_RB_training_type(const std::string & RB_training_type_in);
240  const std::string & get_RB_training_type() const;
241 
246  unsigned int get_Nmax() const { return Nmax; }
247  virtual void set_Nmax(unsigned int Nmax);
248 
253  virtual void load_basis_function(unsigned int i);
254 
259  virtual void load_rb_solution();
260 
265  Real compute_residual_dual_norm_slow(const unsigned int N);
266 
275 
283 
290 
294  SparseMatrix<Number> * get_Aq(unsigned int q);
295 
299  SparseMatrix<Number> * get_non_dirichlet_Aq(unsigned int q);
300 
306 
316  virtual void initialize_rb_construction(bool skip_matrix_assembly=false,
317  bool skip_vector_assembly=false);
318 
322  NumericVector<Number> * get_Fq(unsigned int q);
323 
328 
334 
338  NumericVector<Number> * get_output_vector(unsigned int n, unsigned int q_l);
339 
343  NumericVector<Number> * get_non_dirichlet_output_vector(unsigned int n, unsigned int q_l);
344 
348  virtual void get_all_matrices(std::map<std::string, SparseMatrix<Number> *> & all_matrices);
349 
353  virtual void get_all_vectors(std::map<std::string, NumericVector<Number> *> & all_vectors);
354 
358  virtual void get_output_vectors(std::map<std::string, NumericVector<Number> *> & all_vectors);
359 
365  virtual void assemble_affine_expansion(bool skip_matrix_assembly,
366  bool skip_vector_assembly);
367 
371  void assemble_inner_product_matrix(SparseMatrix<Number> * input_matrix, bool apply_dof_constraints=true);
372 
376  void assemble_Aq_matrix(unsigned int q, SparseMatrix<Number> * input_matrix, bool apply_dof_constraints=true);
377 
381  void assemble_Fq_vector(unsigned int q, NumericVector<Number> * input_vector, bool apply_dof_constraints=true);
382 
387  void add_scaled_Aq(Number scalar, unsigned int q_a,
388  SparseMatrix<Number> * input_matrix,
389  bool symmetrize);
390 
396  virtual void write_riesz_representors_to_files(const std::string & riesz_representors_dir,
397  const bool write_binary_residual_representors);
398 
405  virtual void read_riesz_representors_from_files(const std::string & riesz_representors_dir,
406  const bool write_binary_residual_representors);
407 
408 
416  virtual void recompute_all_residual_terms(const bool compute_inner_products=true);
417 
422  virtual void process_parameters_file(const std::string & parameters_filename);
423 
428  void set_rb_construction_parameters(unsigned int n_training_samples_in,
429  bool deterministic_training_in,
430  int training_parameters_random_seed_in,
431  bool quiet_mode_in,
432  unsigned int Nmax_in,
433  Real rel_training_tolerance_in,
434  Real abs_training_tolerance_in,
435  bool normalize_rb_error_bound_in_greedy_in,
436  const std::string & RB_training_type_in,
437  const RBParameters & mu_min_in,
438  const RBParameters & mu_max_in,
439  const std::map<std::string, std::vector<Real>> & discrete_parameter_values_in,
440  const std::map<std::string,bool> & log_scaling,
441  std::map<std::string, std::vector<RBParameter>> * training_sample_list=nullptr);
442 
446  virtual void print_info() const;
447 
454 
462  unsigned int get_delta_N() const { return delta_N; }
463 
467  void set_inner_product_assembly(ElemAssembly & inner_product_assembly_in);
468 
473 
478  void set_energy_inner_product(const std::vector<Number> & energy_inner_product_coeffs_in);
479 
485 
489  virtual bool check_if_zero_truth_solve() const;
490 
491 #ifdef LIBMESH_ENABLE_DIRICHLET
492 
496  static std::unique_ptr<DirichletBoundary> build_zero_dirichlet_boundary_object();
497 #endif
498 
503  void set_convergence_assertion_flag(bool flag);
504 
508  bool get_preevaluate_thetas_flag() const;
509  void set_preevaluate_thetas_flag(bool flag);
510 
511 
512  //----------- PUBLIC DATA MEMBERS -----------//
513 
521  std::vector<Real> training_error_bounds;
522 
528  std::unique_ptr<LinearSolver<Number>> inner_product_solver;
529 
537 
541  std::unique_ptr<SparseMatrix<Number>> inner_product_matrix;
542 
547  std::vector<Number > truth_outputs;
548 
553  std::vector<std::vector<Number >> output_dual_innerprods;
554 
561  std::vector<std::unique_ptr<NumericVector<Number>>> Fq_representor;
562 
570  std::vector<Number> Fq_representor_innerprods;
571 
580 
588 
595 
603 
611 
618 
628 
635 
642 
643 protected:
644 
649  virtual void allocate_data_structures();
650 
655  virtual void truth_assembly();
656 
662  virtual std::unique_ptr<DGFEMContext> build_context();
663 
670 
675  virtual bool greedy_termination_test(Real abs_greedy_error,
676  Real initial_greedy_error,
677  int count);
678 
684 
693  ElemAssembly * elem_assembly,
694  SparseMatrix<Number> * input_matrix,
695  NumericVector<Number> * input_vector,
696  bool symmetrize=false,
697  bool apply_dof_constraints=true);
698 
708  (DGFEMContext & /*context*/) {}
709 
718  virtual void post_process_truth_solution() {}
719 
726 
733  virtual void assemble_misc_matrices();
734 
739  virtual void assemble_all_affine_operators();
740 
744  virtual void assemble_all_affine_vectors();
745 
749  virtual void assemble_all_output_vectors();
750 
754  virtual void compute_output_dual_innerprods();
755 
768  virtual void compute_Fq_representor_innerprods(bool compute_inner_products=true);
769 
774  virtual void enrich_RB_space();
775 
780  virtual void update_system();
781 
787  virtual Real get_RB_error_bound();
788 
792  virtual void update_RB_system_matrices();
793 
801  virtual void update_residual_terms(bool compute_inner_products=true);
802 
809  virtual void init_context(FEMContext &)
810  {
811  // Failing to rederive init_context() means your FE objects don't
812  // know what to compute.
813  libmesh_deprecated();
814  }
815 
820  bool get_convergence_assertion_flag() const;
821 
826  void check_convergence(LinearSolver<Number> & input_solver);
827 
831  unsigned int get_current_training_parameter_index() const;
832  void set_current_training_parameter_index(unsigned int index);
833 
837  const std::vector<Number> & get_evaluated_thetas(unsigned int training_parameter_index) const;
838 
839  /*
840  * Pre-evaluate the theta functions on the entire (local) training parameter set.
841  */
842  void preevaluate_thetas();
843 
849 
850  //----------- PROTECTED DATA MEMBERS -----------//
851 
856  unsigned int Nmax;
857 
862  unsigned int delta_N;
863 
870 
876 
877 private:
878 
879  //----------- PRIVATE DATA MEMBERS -----------//
880 
887 
893 
898 
904 
915  std::vector<Number> energy_inner_product_coeffs;
916 
920  std::vector<std::unique_ptr<SparseMatrix<Number>>> Aq_vector;
921 
926  std::vector<std::unique_ptr<NumericVector<Number>>> Fq_vector;
927 
932  std::vector<std::vector<std::unique_ptr<NumericVector<Number>>>> outputs_vector;
933 
939  std::vector<std::unique_ptr<SparseMatrix<Number>>> non_dirichlet_Aq_vector;
940  std::vector<std::unique_ptr<NumericVector<Number>>> non_dirichlet_Fq_vector;
941  std::vector<std::vector<std::unique_ptr<NumericVector<Number>>>> non_dirichlet_outputs_vector;
942  std::unique_ptr<SparseMatrix<Number>> non_dirichlet_inner_product_matrix;
943 
950 
956 
963  std::string RB_training_type;
964 
972  std::vector<std::unique_ptr<NumericVector<Number>>> _untransformed_basis_functions;
973 
978  std::unique_ptr<NumericVector<Number>> _untransformed_solution;
979 
984 
991 
996 
1002  std::vector<std::vector<Number>> _evaluated_thetas;
1003 };
1004 
1005 } // namespace libMesh
1006 
1007 #endif // LIBMESH_RB_CONSTRUCTION_H
virtual void read_riesz_representors_from_files(const std::string &riesz_representors_dir, const bool write_binary_residual_representors)
Read in all the Riesz representor data from files.
Real compute_residual_dual_norm_slow(const unsigned int N)
The slow (but simple, non-error prone) way to compute the residual dual norm.
bool normalize_rb_bound_in_greedy
This boolean indicates if we normalize the RB error in the greedy using RBEvaluation::get_error_bound...
virtual void initialize_rb_construction(bool skip_matrix_assembly=false, bool skip_vector_assembly=false)
Allocate all the data structures necessary for the construction stage of the RB method.
This is the EquationSystems class.
std::vector< Real > training_error_bounds
Vector storing the values of the error bound for each parameter in the training set — the parameter ...
bool get_normalize_rb_bound_in_greedy() const
virtual void truth_assembly()
Assemble the truth matrix and right-hand side for current_parameters.
void set_RB_training_type(const std::string &RB_training_type_in)
Get/set the string that determines the training type.
void zero_constrained_dofs_on_vector(NumericVector< Number > &vector) const
It is sometimes useful to be able to zero vector entries that correspond to constrained dofs...
void set_energy_inner_product(const std::vector< Number > &energy_inner_product_coeffs_in)
Specify the coefficients of the A_q operators to be used in the energy inner-product.
std::vector< std::vector< Number > > _evaluated_thetas
Storage of evaluated theta functions at a set of parameters.
std::vector< std::unique_ptr< NumericVector< Number > > > Fq_representor
Vector storing the residual representors associated with the right-hand side.
virtual void post_process_elem_matrix_and_vector(DGFEMContext &)
This function is called from add_scaled_matrix_and_vector() before each element matrix and vector are...
bool assert_convergence
A boolean flag to indicate whether to check for proper convergence after each solve.
virtual Real train_reduced_basis(const bool resize_rb_eval_data=true)
Train the reduced basis.
virtual void get_all_matrices(std::map< std::string, SparseMatrix< Number > *> &all_matrices)
Get a map that stores pointers to all of the matrices.
std::vector< std::unique_ptr< NumericVector< Number > > > _untransformed_basis_functions
In cases where we have dof transformations such as a change of coordinates at some nodes we need to s...
void set_preevaluate_thetas_flag(bool flag)
virtual Real truth_solve(int plot_solution)
Perform a "truth" solve, i.e.
RBConstruction(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
SparseMatrix< Number > * get_Aq(unsigned int q)
Get a pointer to Aq.
bool compute_RB_inner_product
Boolean flag to indicate whether we compute the RB_inner_product_matrix.
virtual void set_context_solution_vec(NumericVector< Number > &vec)
Set current_local_solution = vec so that we can access vec from FEMContext during assembly...
void enrich_basis_from_rhs_terms(const bool resize_rb_eval_data=true)
This function computes one basis function for each rhs term.
std::vector< std::vector< Number > > output_dual_innerprods
The vector storing the dual norm inner product terms for each output.
ElemAssembly & get_inner_product_assembly()
bool get_convergence_assertion_flag() const
Getter for the flag determining if convergence should be checked after each solve.
bool skip_degenerate_sides
In some cases meshes are intentionally created with degenerate sides as a way to represent, say, triangles using a hex-only mesh.
virtual bool greedy_termination_test(Real abs_greedy_error, Real initial_greedy_error, int count)
Function that indicates when to terminate the Greedy basis training.
RBConstructionBase< LinearImplicitSystem > Parent
The type of the parent.
void add_scaled_Aq(Number scalar, unsigned int q_a, SparseMatrix< Number > *input_matrix, bool symmetrize)
Add the scaled q^th affine matrix to input_matrix.
const RBParameters & get_greedy_parameter(unsigned int i)
Return the parameters chosen during the i^th step of the Greedy algorithm.
std::string RB_training_type
This string indicates the type of training that we will use.
SparseMatrix< Number > * get_non_dirichlet_inner_product_matrix()
Get the non-Dirichlet (or more generally no-constraints) version of the inner-product matrix...
NumericVector< Number > * get_non_dirichlet_Fq_if_avail(unsigned int q)
Get a pointer to non_dirichlet_Fq if it&#39;s available, otherwise get Fq.
virtual std::unique_ptr< DGFEMContext > build_context()
Builds a DGFEMContext object with enough information to do evaluations on each element.
std::unique_ptr< SparseMatrix< Number > > inner_product_matrix
The inner product matrix.
const std::string & get_RB_training_type() const
The libMesh namespace provides an interface to certain functionality in the library.
NumericVector< Number > * get_Fq(unsigned int q)
Get a pointer to Fq.
std::vector< Number > Fq_representor_innerprods
Vectors storing the residual representor inner products to be used in computing the residuals online...
unsigned int get_current_training_parameter_index() const
Get/set the current training parameter index.
Real get_rel_training_tolerance() const
void reset_preevaluate_thetas_completed()
Reset the _preevaluate_thetas_completed flag to false.
virtual void print_info() const
Print out info that describes the current setup of this RBConstruction.
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
static std::unique_ptr< DirichletBoundary > build_zero_dirichlet_boundary_object()
It&#39;s helpful to be able to generate a DirichletBoundary that stores a ZeroFunction in order to impose...
std::vector< Number > energy_inner_product_coeffs
We may optionally want to use the "energy inner-product" rather than the inner-product assembly speci...
virtual void init_context(FEMContext &)
Initialize the FEMContext prior to performing an element loop.
virtual Real compute_max_error_bound()
(i) Compute the a posteriori error bound for each set of parameters in the training set...
void print_basis_function_orthogonality() const
Print out a matrix that shows the orthogonality of the RB basis functions.
std::vector< std::unique_ptr< NumericVector< Number > > > non_dirichlet_Fq_vector
virtual void get_output_vectors(std::map< std::string, NumericVector< Number > *> &all_vectors)
Get a map that stores pointers to all of the vectors.
void set_rb_construction_parameters(unsigned int n_training_samples_in, bool deterministic_training_in, int training_parameters_random_seed_in, bool quiet_mode_in, unsigned int Nmax_in, Real rel_training_tolerance_in, Real abs_training_tolerance_in, bool normalize_rb_error_bound_in_greedy_in, const std::string &RB_training_type_in, const RBParameters &mu_min_in, const RBParameters &mu_max_in, const std::map< std::string, std::vector< Real >> &discrete_parameter_values_in, const std::map< std::string, bool > &log_scaling, std::map< std::string, std::vector< RBParameter >> *training_sample_list=nullptr)
Set the state of this RBConstruction object based on the arguments to this function.
void assemble_inner_product_matrix(SparseMatrix< Number > *input_matrix, bool apply_dof_constraints=true)
Assemble the inner product matrix and store it in input_matrix.
LinearSolver< Number > * extra_linear_solver
Also, we store a pointer to an extra linear solver.
virtual void recompute_all_residual_terms(const bool compute_inner_products=true)
This function computes all of the residual representors, can be useful when restarting a basis traini...
virtual void load_basis_function(unsigned int i)
Load the i^th RB function into the RBConstruction solution vector.
virtual void solve_for_matrix_and_rhs(LinearSolver< Number > &input_solver, SparseMatrix< Number > &input_matrix, NumericVector< Number > &input_rhs)
Assembles & solves the linear system A*x=b for the specified matrix input_matrix and right-hand side ...
bool impose_internal_fluxes
Boolean flag to indicate whether we impose "fluxes" (i.e.
bool Fq_representor_innerprods_computed
A boolean flag to indicate whether or not the Fq representor norms have already been computed — used...
bool use_empty_rb_solve_in_greedy
A boolean flag to indicate whether or not we initialize the Greedy algorithm by performing rb_solves ...
bool use_energy_inner_product
Boolean to indicate whether we&#39;re using the energy inner-product.
unsigned int number() const
Definition: system.h:2269
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the RBThetaExpansion object that that belongs to rb_eval.
bool get_preevaluate_thetas_flag() const
Get/set flag to pre-evaluate the theta functions.
void set_rb_assembly_expansion(RBAssemblyExpansion &rb_assembly_expansion_in)
Set the rb_assembly_expansion object.
std::unique_ptr< SparseMatrix< Number > > non_dirichlet_inner_product_matrix
void update_greedy_param_list()
Update the list of Greedily chosen parameters with current_parameters.
Real rel_training_tolerance
Relative and absolute tolerances for training reduced basis using the Greedy scheme.
unsigned int _current_training_parameter_index
The current training parameter index during reduced basis training.
unsigned int get_Nmax() const
Get/set Nmax, the maximum number of RB functions we are willing to compute.
RBAssemblyExpansion & get_rb_assembly_expansion()
SparseMatrix< Number > * get_inner_product_matrix()
Get a pointer to inner_product_matrix.
virtual void update_RB_system_matrices()
Compute the reduced basis matrices for the current basis.
virtual void update_residual_terms(bool compute_inner_products=true)
Compute the terms that are combined ‘online’ to determine the dual norm of the residual.
bool store_untransformed_basis
Boolean flag to indicate whether we store a second copy of the basis without constraints or dof trans...
void assemble_Fq_vector(unsigned int q, NumericVector< Number > *input_vector, bool apply_dof_constraints=true)
Assemble the q^th affine vector and store it in input_matrix.
RBConstruction sys_type
The type of system.
const std::vector< Number > & get_evaluated_thetas(unsigned int training_parameter_index) const
Return the evaluated theta functions at the given training parameter index.
virtual void update_system()
Update the system after enriching the RB space; this calls a series of functions to update the system...
Real get_abs_training_tolerance() const
This class extends FEMContext in order to provide extra data required to perform local element residu...
virtual SparseMatrix< Number > & get_matrix_for_output_dual_solves()
Return the matrix for the output residual dual norm solves.
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62
NumericVector< Number > * get_non_dirichlet_Fq(unsigned int q)
Get a pointer to non-Dirichlet Fq.
virtual Real get_RB_error_bound()
unsigned int get_delta_N() const
Get delta_N, the number of basis functions we add to the RB space per iteration of the greedy algorit...
SparseMatrix< Number > * get_non_dirichlet_Aq(unsigned int q)
Get a pointer to non_dirichlet_Aq.
void set_inner_product_assembly(ElemAssembly &inner_product_assembly_in)
Set the rb_assembly_expansion object.
std::vector< std::unique_ptr< SparseMatrix< Number > > > Aq_vector
Vector storing the Q_a matrices from the affine expansion.
void add_scaled_matrix_and_vector(Number scalar, ElemAssembly *elem_assembly, SparseMatrix< Number > *input_matrix, NumericVector< Number > *input_vector, bool symmetrize=false, bool apply_dof_constraints=true)
This function loops over the mesh and applies the specified interior and/or boundary assembly routine...
bool skip_residual_in_train_reduced_basis
Boolean flag to indicate if we skip residual calculations in train_reduced_basis. ...
void set_convergence_assertion_flag(bool flag)
Setter for the flag determining if convergence should be checked after each solve.
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:52
bool is_rb_eval_initialized() const
unsigned int delta_N
The number of basis functions that we add at each greedy step.
Real train_reduced_basis_with_greedy(const bool resize_rb_eval_data)
Train the reduced basis using the "Greedy algorithm.".
RBEvaluation * rb_eval
The current RBEvaluation object we are using to perform the Evaluation stage of the reduced basis met...
std::unique_ptr< NumericVector< Number > > _untransformed_solution
We also store a copy of the untransformed solution in order to create _untransformed_basis_functions...
This class is part of the rbOOmit framework.
Definition: rb_evaluation.h:50
virtual void post_process_truth_solution()
Similarly, provide an opportunity to post-process the truth solution after the solve is complete...
void set_current_training_parameter_index(unsigned int index)
void set_abs_training_tolerance(Real new_training_tolerance)
Get/set the absolute tolerance for the basis training.
bool store_non_dirichlet_operators
Boolean flag to indicate whether we store a second copy of each affine operator and vector which does...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void clear() override
Clear all the data structures associated with the system.
ElemAssembly * inner_product_assembly
Pointer to inner product assembly.
virtual std::string system_type() const override
void check_convergence(LinearSolver< Number > &input_solver)
Check if the linear solver reports convergence.
std::unique_ptr< LinearSolver< Number > > inner_product_solver
We store an extra linear solver object which we can optionally use for solving all systems in which t...
void set_normalize_rb_bound_in_greedy(bool normalize_rb_bound_in_greedy_in)
Get/set the boolean to indicate if we normalize the RB error in the greedy.
std::vector< std::unique_ptr< NumericVector< Number > > > Fq_vector
Vector storing the Q_f vectors in the affine decomposition of the right-hand side.
This class stores the set of ElemAssembly functor objects that define the "parameter-independent expa...
ElemAssembly provides a per-element (interior and boundary) assembly functionality.
Definition: elem_assembly.h:38
std::vector< Number > truth_outputs
Vector storing the truth output values from the most recent truth solve.
SparseMatrix< Number > * get_non_dirichlet_inner_product_matrix_if_avail()
Get the non-Dirichlet inner-product matrix if it&#39;s available, otherwise get the inner-product matrix ...
virtual void compute_Fq_representor_innerprods(bool compute_inner_products=true)
Compute the terms that are combined ‘online’ to determine the dual norm of the residual.
virtual void set_Nmax(unsigned int Nmax)
std::vector< std::unique_ptr< SparseMatrix< Number > > > non_dirichlet_Aq_vector
We may also need a second set of matrices/vectors that do not have the Dirichlet boundary conditions ...
void set_rel_training_tolerance(Real new_training_tolerance)
Get/set the relative tolerance for the basis training.
std::vector< std::vector< std::unique_ptr< NumericVector< Number > > > > non_dirichlet_outputs_vector
virtual void write_riesz_representors_to_files(const std::string &riesz_representors_dir, const bool write_binary_residual_representors)
Write out all the Riesz representor data to files.
virtual void assemble_all_output_vectors()
Assemble and store the output vectors.
bool _preevaluate_thetas_flag
Flag to indicate if we preevaluate the theta functions.
void assemble_Aq_matrix(unsigned int q, SparseMatrix< Number > *input_matrix, bool apply_dof_constraints=true)
Assemble the q^th affine matrix and store it in input_matrix.
This class is part of the rbOOmit framework.
virtual void process_parameters_file(const std::string &parameters_filename)
Read in from the file specified by parameters_filename and set the this system&#39;s member variables acc...
RBAssemblyExpansion * rb_assembly_expansion
This member holds the (parameter independent) assembly functors that define the "affine expansion" of...
This class is part of the rbOOmit framework.
NumericVector< Number > * get_output_vector(unsigned int n, unsigned int q_l)
Get a pointer to the n^th output.
RBConstruction & operator=(const RBConstruction &)=delete
virtual bool check_if_zero_truth_solve() const
RBEvaluation & get_rb_evaluation()
Get a reference to the RBEvaluation object.
virtual void assemble_misc_matrices()
Assemble and store all the inner-product matrix, the constraint matrix (for constrained problems) and...
virtual void get_all_vectors(std::map< std::string, NumericVector< Number > *> &all_vectors)
Get a map that stores pointers to all of the vectors.
const std::string & name() const
Definition: system.h:2261
virtual void load_rb_solution()
Load the RB solution from the most recent solve with rb_eval into this system&#39;s solution vector...
void set_rb_evaluation(RBEvaluation &rb_eval_in)
Set the RBEvaluation object.
virtual void compute_output_dual_innerprods()
Compute and store the dual norm of each output functional.
bool exit_on_repeated_greedy_parameters
Boolean flag to indicate whether we exit the greedy if we select the same parameters twice in a row...
bool output_dual_innerprods_computed
A boolean flag to indicate whether or not the output dual norms have already been computed — used to...
void train_reduced_basis_with_POD()
Train the reduced basis using Proper Orthogonal Decomposition (POD).
virtual void enrich_RB_space()
Add a new basis function to the RB space.
virtual void assemble_all_affine_vectors()
Assemble and store the affine RHS vectors.
NumericVector< Number > * get_non_dirichlet_output_vector(unsigned int n, unsigned int q_l)
Get a pointer to non-Dirichlet output vector.
virtual void allocate_data_structures()
Helper function that actually allocates all the data structures required by this class.
bool _preevaluate_thetas_completed
Flag to indicate if the preevaluate_thetas function has been called, since this allows us to avoid ca...
SparseMatrix< Number > * get_non_dirichlet_Aq_if_avail(unsigned int q)
Get a pointer to non_dirichlet_Aq if it&#39;s available, otherwise get Aq.
unsigned int Nmax
Maximum number of reduced basis functions we are willing to use.
virtual void assemble_affine_expansion(bool skip_matrix_assembly, bool skip_vector_assembly)
Assemble the matrices and vectors for this system.
std::vector< std::vector< std::unique_ptr< NumericVector< Number > > > > outputs_vector
The libMesh vectors that define the output functionals.
virtual void assemble_all_affine_operators()
Assemble and store all Q_a affine operators as well as the inner-product matrix.