libMesh
transient_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_TRANSIENT_RB_CONSTRUCTION_H
21 #define LIBMESH_TRANSIENT_RB_CONSTRUCTION_H
22 
23 // rbOOmit includes
24 #include "libmesh/rb_construction.h"
25 #include "libmesh/transient_rb_evaluation.h"
26 #include "libmesh/rb_temporal_discretization.h"
27 
28 // libMesh includes
29 #include "libmesh/transient_system.h"
30 
31 // C++ includes
32 
33 namespace libMesh
34 {
35 
49 {
50 public:
51 
57  const std::string & name,
58  const unsigned int number);
59 
63  virtual ~TransientRBConstruction ();
64 
69 
74 
79  virtual void clear () libmesh_override;
80 
89  virtual void initialize_rb_construction(bool skip_matrix_assembly=false,
90  bool skip_vector_assembly=false) libmesh_override;
91 
95  virtual Real truth_solve(int write_interval) libmesh_override;
96 
104  virtual Real train_reduced_basis(const bool resize_rb_eval_data=true) libmesh_override;
105 
110  virtual void process_parameters_file (const std::string & parameters_filename) libmesh_override;
111 
115  virtual void print_info() libmesh_override;
116 
121  virtual bool greedy_termination_test(Real abs_greedy_error,
122  Real initial_greedy_error,
123  int count) libmesh_override;
124 
129  virtual void assemble_all_affine_operators() libmesh_override;
130 
134  virtual void assemble_misc_matrices() libmesh_override;
135 
139  void assemble_L2_matrix(SparseMatrix<Number> * input_matrix,
140  bool apply_dirichlet_bc=true);
141 
146  void assemble_mass_matrix(SparseMatrix<Number> * input_matrix);
147 
152  void add_scaled_mass_matrix(Number scalar,
153  SparseMatrix<Number> * input_matrix);
154 
159  void mass_matrix_scaled_matvec(Number scalar,
160  NumericVector<Number> & dest,
161  NumericVector<Number> & arg);
162 
166  void set_L2_assembly(ElemAssembly & L2_assembly_in);
167 
172 
176  void assemble_Mq_matrix(unsigned int q,
177  SparseMatrix<Number> * input_matrix,
178  bool apply_dirichlet_bc=true);
179 
183  SparseMatrix<Number> * get_M_q(unsigned int q);
184 
188  SparseMatrix<Number> * get_non_dirichlet_M_q(unsigned int q);
189 
193  virtual void get_all_matrices(std::map<std::string, SparseMatrix<Number> *> & all_matrices) libmesh_override;
194 
198  virtual void truth_assembly() libmesh_override;
199 
208  int get_max_truth_solves() const { return max_truth_solves; }
209  void set_max_truth_solves(int max_truth_solves_in) { this->max_truth_solves = max_truth_solves_in; }
210 
214  Real get_POD_tol() const { return POD_tol; }
215  void set_POD_tol(const Real POD_tol_in) { this->POD_tol = POD_tol_in; }
216 
221  void set_delta_N(const unsigned int new_delta_N) { this->delta_N = new_delta_N; }
222 
227  virtual void load_rb_solution() libmesh_override;
228 
236 
243 
248  virtual void write_riesz_representors_to_files(const std::string & riesz_representors_dir,
249  const bool write_binary_residual_representors) libmesh_override;
250 
255  virtual void read_riesz_representors_from_files(const std::string & riesz_representors_dir,
256  const bool write_binary_residual_representors) libmesh_override;
257 
258 
259  //----------- PUBLIC DATA MEMBERS -----------//
260 
265 
271 
275  std::vector<SparseMatrix<Number> *> M_q_vector;
276 
282  std::vector<SparseMatrix<Number> *> non_dirichlet_M_q_vector;
283 
288  std::vector<std::vector<Number>> truth_outputs_all_k;
289 
295 
302 
307  std::string init_filename;
308 
309 protected:
310 
315  virtual void allocate_data_structures() libmesh_override;
316 
321  virtual void assemble_affine_expansion(bool skip_matrix_assembly,
322  bool skip_vector_assembly) libmesh_override;
323 
329  virtual void initialize_truth();
330 
335  virtual SparseMatrix<Number> & get_matrix_for_output_dual_solves() libmesh_override;
336 
341  void add_IC_to_RB_space();
342 
348  virtual void enrich_RB_space() libmesh_override;
349 
353  virtual void update_system() libmesh_override;
354 
358  virtual void update_RB_system_matrices() libmesh_override;
359 
364  virtual void update_residual_terms(bool compute_inner_products) libmesh_override;
365 
372 
373  //----------- PROTECTED DATA MEMBERS -----------//
374 
381 
389 
394 
400 
401 private:
402 
403  //----------- PRIVATE DATA MEMBERS -----------//
404 
409 };
410 
411 } // namespace libMesh
412 
413 #endif // LIBMESH_TRANSIENT_RB_CONSTRUCTION_H
virtual void write_riesz_representors_to_files(const std::string &riesz_representors_dir, const bool write_binary_residual_representors) libmesh_override
Write out all the Riesz representor data to files.
virtual void update_system() libmesh_override
Update the system after enriching the RB space.
virtual void read_riesz_representors_from_files(const std::string &riesz_representors_dir, const bool write_binary_residual_representors) libmesh_override
Write out all the Riesz representor data to files.
This is the EquationSystems class.
int max_truth_solves
Maximum number of truth solves in the POD-Greedy.
void set_max_truth_solves(int max_truth_solves_in)
virtual void clear() libmesh_override
Clear all the data structures associated with the system.
virtual void load_rb_solution() libmesh_override
Load the RB solution from the current time-level into the libMesh solution vector.
std::vector< SparseMatrix< Number > * > non_dirichlet_M_q_vector
We sometimes also need a second set of M_q matrices that do not have the Dirichlet boundary condition...
virtual bool greedy_termination_test(Real abs_greedy_error, Real initial_greedy_error, int count) libmesh_override
Function that indicates when to terminate the Greedy basis training.
void assemble_mass_matrix(SparseMatrix< Number > *input_matrix)
Assemble the mass matrix at the current parameter and store it in input_matrix.
UniquePtr< SparseMatrix< Number > > L2_matrix
The L2 matrix.
SparseMatrix< Number > * get_M_q(unsigned int q)
Get a pointer to M_q.
virtual void enrich_RB_space() libmesh_override
Add a new basis functions to the RB space.
virtual void print_info() libmesh_override
Print out info that describes the current setup of this RBConstruction.
virtual SparseMatrix< Number > & get_matrix_for_output_dual_solves() libmesh_override
Override to return the L2 product matrix for output dual norm solves for transient state problems...
virtual void initialize_truth()
This function imposes a truth initial condition, defaults to zero initial condition if the flag nonze...
virtual void assemble_all_affine_operators() libmesh_override
Assemble and store all the affine operators.
TransientSystem< RBConstruction > Parent
The type of the parent.
virtual void initialize_rb_construction(bool skip_matrix_assembly=false, bool skip_vector_assembly=false) libmesh_override
Allocate all the data structures necessary for the construction stage of the RB method.
Numeric vector.
Definition: dof_map.h:66
This class provides a specific system class.
Real POD_tol
If positive, this tolerance determines the number of POD modes we add to the space on a call to enric...
std::vector< std::vector< Number > > truth_outputs_all_k
The truth outputs for all time-levels from the most recent truth_solve.
virtual void allocate_data_structures() libmesh_override
Helper function that actually allocates all the data structures required by this class.
void update_RB_initial_condition_all_N()
Compute the L2 projection of the initial condition onto the RB space for 1 <= N <= RB_size and store ...
The libMesh namespace provides an interface to certain functionality in the library.
std::string init_filename
The filename of the file containing the initial condition projected onto the truth mesh...
const std::string & name() const
Definition: system.h:1998
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual void assemble_misc_matrices() libmesh_override
Override to assemble the L2 matrix as well.
Generic sparse matrix.
Definition: dof_map.h:65
std::vector< SparseMatrix< Number > * > M_q_vector
Vector storing the Q_m matrices from the mass operator.
void assemble_L2_matrix(SparseMatrix< Number > *input_matrix, bool apply_dirichlet_bc=true)
Assemble the L2 matrix.
virtual void process_parameters_file(const std::string &parameters_filename) libmesh_override
Read in the parameters from file and set up the system accordingly.
void set_L2_assembly(ElemAssembly &L2_assembly_in)
Set the L2 object.
Number set_error_temporal_data()
Set column k (i.e.
virtual void truth_assembly() libmesh_override
Assemble the truth system in the transient linear case.
bool nonzero_initialization
Boolean flag to indicate whether we are using a non-zero initialization.
This class is part of the rbOOmit framework.
ElemAssembly * L2_assembly
Function pointer for assembling the L2 matrix.
unsigned int delta_N
The number of basis functions that we add at each greedy step.
TransientRBConstruction sys_type
The type of system.
TransientRBConstruction(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
void set_POD_tol(const Real POD_tol_in)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< NumericVector< Number > * > temporal_data
Dense matrix to store the data that we use for the temporal POD.
virtual void update_residual_terms(bool compute_inner_products) libmesh_override
Compute the terms that are combined `online&#39; to determine the dual norm of the residual.
unsigned int number() const
Definition: system.h:2006
ElemAssembly provides a per-element (interior and boundary) assembly functionality.
Definition: elem_assembly.h:38
UniquePtr< SparseMatrix< Number > > non_dirichlet_L2_matrix
The L2 matrix without Dirichlet conditions enforced.
void assemble_Mq_matrix(unsigned int q, SparseMatrix< Number > *input_matrix, bool apply_dirichlet_bc=true)
Assemble the q^th affine term of the mass matrix and store it in input_matrix.
void add_scaled_mass_matrix(Number scalar, SparseMatrix< Number > *input_matrix)
Add the scaled mass matrix (assembled for the current parameter) to input_matrix. ...
int get_max_truth_solves() const
Get/set max_truth_solves, the maximum number of RB truth solves we are willing to compute in the tran...
const NumericVector< Number > & get_error_temporal_data()
Get the column of temporal_data corresponding to the current time level.
Define a class that encapsulates the details of a "generalized Euler" temporal discretization to be u...
bool compute_truth_projection_error
Boolean flag that indicates whether we will compute the projection error for the truth solution into ...
Defines a dense vector for use in Finite Element-type computations.
void mass_matrix_scaled_matvec(Number scalar, NumericVector< Number > &dest, NumericVector< Number > &arg)
Perform a matrix-vector multiplication with the current mass matrix and store the result in dest...
virtual Real truth_solve(int write_interval) libmesh_override
Perform a truth solve at the current parameter.
virtual Real train_reduced_basis(const bool resize_rb_eval_data=true) libmesh_override
Train the reduced basis.
void set_delta_N(const unsigned int new_delta_N)
Set delta_N, the number of basis functions we add to the RB space from each POD.
void add_IC_to_RB_space()
Initialize RB space by adding the truth initial condition as the first RB basis function.
SparseMatrix< Number > * get_non_dirichlet_M_q(unsigned int q)
Get a pointer to non_dirichlet_M_q.
virtual void assemble_affine_expansion(bool skip_matrix_assembly, bool skip_vector_assembly) libmesh_override
Override assemble_affine_expansion to also initialize RB_ic_proj_rhs_all_N, if necessary.
virtual void update_RB_system_matrices() libmesh_override
Compute the reduced basis matrices for the current basis.
virtual void get_all_matrices(std::map< std::string, SparseMatrix< Number > * > &all_matrices) libmesh_override
Get a map that stores pointers to all of the matrices.
Real get_POD_tol() const
Get/set POD_tol.
DenseVector< Number > RB_ic_proj_rhs_all_N
The vector that stores the right-hand side for the initial condition projections. ...