rb_construction.h

Go to the documentation of this file.

// rbOOmit: An implementation of the Certified Reduced Basis method.
// Copyright (C) 2009, 2010 David J. Knezevic

// This file is part of rbOOmit.

// rbOOmit is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.

// rbOOmit is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.

// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

#ifndef LIBMESH_RB_CONSTRUCTION_H
#define LIBMESH_RB_CONSTRUCTION_H

// rbOOmit includes
#include "libmesh/rb_construction_base.h"

// libMesh includes
#include "libmesh/linear_implicit_system.h"
#include "libmesh/dense_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dg_fem_context.h"
#include "libmesh/dirichlet_boundaries.h"

// C++ includes

namespace libMesh
{

class RBThetaExpansion;
class RBAssemblyExpansion;
class RBEvaluation;
class ElemAssembly;

class RBConstruction : public RBConstructionBase<LinearImplicitSystem>
{
public:

  RBConstruction (EquationSystems & es,
                  const std::string & name,
                  const unsigned int number);

  virtual ~RBConstruction ();

  typedef RBConstruction sys_type;

  virtual void solve_for_matrix_and_rhs (LinearSolver<Number> & input_solver,
                                         SparseMatrix<Number> & input_matrix,
                                         NumericVector<Number> & input_rhs);

  void set_rb_evaluation(RBEvaluation & rb_eval_in);

  RBEvaluation & get_rb_evaluation();

  bool is_rb_eval_initialized() const;

  RBThetaExpansion & get_rb_theta_expansion();

  void set_rb_assembly_expansion(RBAssemblyExpansion & rb_assembly_expansion_in);

  RBAssemblyExpansion & get_rb_assembly_expansion();

  sys_type & system () { return *this; }

  typedef RBConstructionBase<LinearImplicitSystem> Parent;

  virtual void clear () override;

  virtual std::string system_type () const override;

  virtual Real truth_solve(int plot_solution);

  virtual Real train_reduced_basis(const bool resize_rb_eval_data=true);

  void enrich_basis_from_rhs_terms(const bool resize_rb_eval_data=true);

  virtual Real compute_max_error_bound();

  const RBParameters & get_greedy_parameter(unsigned int i);

  void set_rel_training_tolerance(Real new_training_tolerance)
  {this->rel_training_tolerance = new_training_tolerance; }
  Real get_rel_training_tolerance() { return rel_training_tolerance; }

  void set_abs_training_tolerance(Real new_training_tolerance)
  {this->abs_training_tolerance = new_training_tolerance; }
  Real get_abs_training_tolerance() { return abs_training_tolerance; }

  void set_normalize_rb_bound_in_greedy(bool normalize_rb_bound_in_greedy_in)
  {this->normalize_rb_bound_in_greedy = normalize_rb_bound_in_greedy_in; }
  bool get_normalize_rb_bound_in_greedy() { return normalize_rb_bound_in_greedy; }

  unsigned int get_Nmax() const    { return Nmax; }
  virtual void set_Nmax(unsigned int Nmax);

  virtual void load_basis_function(unsigned int i);

  virtual void load_rb_solution();

  SparseMatrix<Number> * get_inner_product_matrix();

  SparseMatrix<Number> * get_non_dirichlet_inner_product_matrix();

  SparseMatrix<Number> * get_non_dirichlet_inner_product_matrix_if_avail();

  SparseMatrix<Number> * get_Aq(unsigned int q);

  SparseMatrix<Number> * get_non_dirichlet_Aq(unsigned int q);

  SparseMatrix<Number> * get_non_dirichlet_Aq_if_avail(unsigned int q);

  virtual void initialize_rb_construction(bool skip_matrix_assembly=false,
                                          bool skip_vector_assembly=false);

  NumericVector<Number> * get_Fq(unsigned int q);

  NumericVector<Number> * get_non_dirichlet_Fq(unsigned int q);

  NumericVector<Number> * get_non_dirichlet_Fq_if_avail(unsigned int q);

  NumericVector<Number> * get_output_vector(unsigned int n, unsigned int q_l);

  NumericVector<Number> * get_non_dirichlet_output_vector(unsigned int n, unsigned int q_l);

  virtual void get_all_matrices(std::map<std::string, SparseMatrix<Number> *> & all_matrices);

  virtual void get_all_vectors(std::map<std::string, NumericVector<Number> *> & all_vectors);

  virtual void get_output_vectors(std::map<std::string, NumericVector<Number> *> & all_vectors);

  virtual void assemble_affine_expansion(bool skip_matrix_assembly,
                                         bool skip_vector_assembly);

  void assemble_inner_product_matrix(SparseMatrix<Number> * input_matrix, bool apply_dof_constraints=true);

  void assemble_Aq_matrix(unsigned int q, SparseMatrix<Number> * input_matrix, bool apply_dof_constraints=true);

  void assemble_Fq_vector(unsigned int q, NumericVector<Number> * input_vector, bool apply_dof_constraints=true);

  void add_scaled_Aq(Number scalar, unsigned int q_a,
                     SparseMatrix<Number> * input_matrix,
                     bool symmetrize);

  virtual void write_riesz_representors_to_files(const std::string & riesz_representors_dir,
                                                 const bool write_binary_residual_representors);

  virtual void read_riesz_representors_from_files(const std::string & riesz_representors_dir,
                                                  const bool write_binary_residual_representors);


  virtual void recompute_all_residual_terms(const bool compute_inner_products=true);

  virtual void process_parameters_file(const std::string & parameters_filename);

  void set_rb_construction_parameters(unsigned int n_training_samples_in,
                                      bool deterministic_training_in,
                                      unsigned 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,
                                      RBParameters mu_min_in,
                                      RBParameters mu_max_in,
                                      std::map<std::string, std::vector<Real>> discrete_parameter_values_in,
                                      std::map<std::string,bool> log_scaling);

  virtual void print_info();

  void print_basis_function_orthogonality();

  unsigned int get_delta_N() const { return delta_N; }

  void set_inner_product_assembly(ElemAssembly & inner_product_assembly_in);

  ElemAssembly & get_inner_product_assembly();

  void zero_constrained_dofs_on_vector(NumericVector<Number> & vector);

  static std::unique_ptr<DirichletBoundary> build_zero_dirichlet_boundary_object();

  void set_convergence_assertion_flag(bool flag);

  //----------- PUBLIC DATA MEMBERS -----------//

  std::vector<Real> training_error_bounds;

  std::unique_ptr<LinearSolver<Number>> inner_product_solver;

  LinearSolver<Number> * extra_linear_solver;

  std::unique_ptr<SparseMatrix<Number>> inner_product_matrix;

  std::vector<Number > truth_outputs;

  std::vector<std::vector<Number >> output_dual_innerprods;

  std::vector<std::unique_ptr<NumericVector<Number>>> Fq_representor;

  std::vector<Number> Fq_representor_innerprods;

  bool skip_residual_in_train_reduced_basis;

  bool exit_on_repeated_greedy_parameters;

  bool impose_internal_fluxes;

  bool skip_degenerate_sides;

  bool compute_RB_inner_product;

  bool store_non_dirichlet_operators;

  bool use_empty_rb_solve_in_greedy;

  bool Fq_representor_innerprods_computed;

protected:

  virtual void allocate_data_structures();

  virtual void truth_assembly();

  virtual std::unique_ptr<DGFEMContext> build_context();

  virtual SparseMatrix<Number> & get_matrix_for_output_dual_solves();

  virtual bool greedy_termination_test(Real abs_greedy_error,
                                       Real initial_greedy_error,
                                       int count);

  void update_greedy_param_list();

  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);

  virtual void set_context_solution_vec(NumericVector<Number> & vec);

  virtual void assemble_misc_matrices();

  virtual void assemble_all_affine_operators();

  virtual void assemble_all_affine_vectors();

  virtual void assemble_all_output_vectors();

  virtual void compute_output_dual_innerprods();

  virtual void compute_Fq_representor_innerprods(bool compute_inner_products=true);

  virtual void enrich_RB_space();

  virtual void update_system();

  virtual Real get_RB_error_bound();

  virtual void update_RB_system_matrices();

  virtual void update_residual_terms(bool compute_inner_products=true);

  virtual void init_context(FEMContext &) {}

  bool get_convergence_assertion_flag() const;

  void check_convergence(LinearSolver<Number> & input_solver);

  //----------- PROTECTED DATA MEMBERS -----------//

  unsigned int Nmax;

  unsigned int delta_N;

  bool output_dual_innerprods_computed;

  bool assert_convergence;

private:

  //----------- PRIVATE DATA MEMBERS -----------//

  RBEvaluation * rb_eval;

  RBAssemblyExpansion * rb_assembly_expansion;

  ElemAssembly * inner_product_assembly;

  std::vector<std::unique_ptr<SparseMatrix<Number>>> Aq_vector;

  std::vector<std::unique_ptr<NumericVector<Number>>> Fq_vector;

  std::vector<std::vector<std::unique_ptr<NumericVector<Number>>>> outputs_vector;

  std::vector<std::unique_ptr<SparseMatrix<Number>>> non_dirichlet_Aq_vector;
  std::vector<std::unique_ptr<NumericVector<Number>>> non_dirichlet_Fq_vector;
  std::vector<std::vector<std::unique_ptr<NumericVector<Number>>>> non_dirichlet_outputs_vector;
  std::unique_ptr<SparseMatrix<Number>> non_dirichlet_inner_product_matrix;

  Real rel_training_tolerance;
  Real abs_training_tolerance;

  bool normalize_rb_bound_in_greedy;

};

} // namespace libMesh

#endif // LIBMESH_RB_CONSTRUCTION_H