libMesh
rb_scm_construction.C
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 // Configuration data
21 #include "libmesh/libmesh_config.h"
22 
23 // Currently, the RBSCMConstruction should only be available
24 // if SLEPc support is enabled.
25 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK)
26 
27 #include "libmesh/rb_scm_construction.h"
28 #include "libmesh/rb_construction.h"
29 #include "libmesh/rb_scm_evaluation.h"
30 
31 #include "libmesh/libmesh_logging.h"
32 #include "libmesh/numeric_vector.h"
33 #include "libmesh/sparse_matrix.h"
34 #include "libmesh/equation_systems.h"
35 #include "libmesh/getpot.h"
36 #include "libmesh/parallel.h"
37 #include "libmesh/dof_map.h"
38 // For creating a directory
39 #include <sys/types.h>
40 #include <sys/stat.h>
41 #include <errno.h>
42 
43 namespace libMesh
44 {
45 
47  const std::string & name_in,
48  const unsigned int number_in)
49  : Parent(es, name_in, number_in),
50  SCM_training_tolerance(0.5),
51  RB_system_name(""),
52  rb_scm_eval(libmesh_nullptr)
53 {
54 
55  // set assemble_before_solve flag to false
56  // so that we control matrix assembly.
57  assemble_before_solve = false;
58 
59  // We symmetrize all operators hence use symmetric solvers
61 
62 }
63 
65 {
66  this->clear();
67 }
68 
69 
71 {
72  Parent::clear();
73 }
74 
76 {
77  rb_scm_eval = &rb_scm_eval_in;
78 }
79 
81 {
82  if (!rb_scm_eval)
83  libmesh_error_msg("Error: RBSCMEvaluation object hasn't been initialized yet");
84 
85  return *rb_scm_eval;
86 }
87 
89 {
91 }
92 
93 void RBSCMConstruction::process_parameters_file(const std::string & parameters_filename)
94 {
95  // First read in data from parameters_filename
96  GetPot infile(parameters_filename);
97  const unsigned int n_training_samples = infile("n_training_samples",1);
98  const bool deterministic_training = infile("deterministic_training",false);
99 
100  // Read in training_parameters_random_seed value. This is used to
101  // seed the RNG when picking the training parameters. By default the
102  // value is -1, which means use std::time to seed the RNG.
103  unsigned int training_parameters_random_seed_in = static_cast<int>(-1);
104  training_parameters_random_seed_in = infile("training_parameters_random_seed",
105  training_parameters_random_seed_in);
106  set_training_random_seed(training_parameters_random_seed_in);
107 
108  // SCM Greedy termination tolerance
109  const Real SCM_training_tolerance_in = infile("SCM_training_tolerance", SCM_training_tolerance);
110  set_SCM_training_tolerance(SCM_training_tolerance_in);
111 
112  // Initialize the parameter ranges and the parameters themselves
113  unsigned int n_continuous_parameters = infile.vector_variable_size("parameter_names");
114  RBParameters mu_min_in;
115  RBParameters mu_max_in;
116  for (unsigned int i=0; i<n_continuous_parameters; i++)
117  {
118  // Read in the parameter names
119  std::string param_name = infile("parameter_names", "NONE", i);
120 
121  {
122  Real min_val = infile(param_name, 0., 0);
123  mu_min_in.set_value(param_name, min_val);
124  }
125 
126  {
127  Real max_val = infile(param_name, 0., 1);
128  mu_max_in.set_value(param_name, max_val);
129  }
130  }
131 
132  std::map<std::string, std::vector<Real>> discrete_parameter_values_in;
133 
134  unsigned int n_discrete_parameters = infile.vector_variable_size("discrete_parameter_names");
135  for (unsigned int i=0; i<n_discrete_parameters; i++)
136  {
137  std::string param_name = infile("discrete_parameter_names", "NONE", i);
138 
139  unsigned int n_vals_for_param = infile.vector_variable_size(param_name);
140  std::vector<Real> vals_for_param(n_vals_for_param);
141  for (std::size_t j=0; j<vals_for_param.size(); j++)
142  vals_for_param[j] = infile(param_name, 0., j);
143 
144  discrete_parameter_values_in[param_name] = vals_for_param;
145  }
146 
147  initialize_parameters(mu_min_in, mu_max_in, discrete_parameter_values_in);
148 
149  std::map<std::string,bool> log_scaling;
150  const RBParameters & mu = get_parameters();
152  RBParameters::const_iterator it_end = mu.end();
153  unsigned int i=0;
154  for ( ; it != it_end; ++it)
155  {
156  std::string param_name = it->first;
157  log_scaling[param_name] = static_cast<bool>(infile("log_scaling", 0, i));
158  i++;
159  }
160 
162  this->get_parameters_max(),
163  n_training_samples,
164  log_scaling,
165  deterministic_training); // use deterministic parameters
166 }
167 
169 {
170  // Print out info that describes the current setup
171  libMesh::out << std::endl << "RBSCMConstruction parameters:" << std::endl;
172  libMesh::out << "system name: " << this->name() << std::endl;
173  libMesh::out << "SCM Greedy tolerance: " << get_SCM_training_tolerance() << std::endl;
174  if (rb_scm_eval)
175  {
176  libMesh::out << "A_q operators attached: " << get_rb_theta_expansion().get_n_A_terms() << std::endl;
177  }
178  else
179  {
180  libMesh::out << "RBThetaExpansion member is not set yet" << std::endl;
181  }
182  libMesh::out << "Number of parameters: " << get_n_params() << std::endl;
185  for ( ; it != it_end; ++it)
186  {
187  std::string param_name = it->first;
188  libMesh::out << "Parameter " << param_name
189  << ": Min = " << get_parameter_min(param_name)
190  << ", Max = " << get_parameter_max(param_name) << std::endl;
191  }
193  libMesh::out << "n_training_samples: " << get_n_training_samples() << std::endl;
194  libMesh::out << std::endl;
195 }
196 
198 {
199  // Clear SCM data vectors
200  rb_scm_eval->B_min.clear();
201  rb_scm_eval->B_max.clear();
202  rb_scm_eval->C_J.clear();
204  for (std::size_t i=0; i<rb_scm_eval->SCM_UB_vectors.size(); i++)
205  rb_scm_eval->SCM_UB_vectors[i].clear();
206  rb_scm_eval->SCM_UB_vectors.clear();
207 
208  // Resize the bounding box vectors
209  rb_scm_eval->B_min.resize(get_rb_theta_expansion().get_n_A_terms());
210  rb_scm_eval->B_max.resize(get_rb_theta_expansion().get_n_A_terms());
211 }
212 
213 void RBSCMConstruction::add_scaled_symm_Aq(unsigned int q_a, Number scalar)
214 {
215  LOG_SCOPE("add_scaled_symm_Aq()", "RBSCMConstruction");
216  // Load the operators from the RBConstruction
217  EquationSystems & es = this->get_equation_systems();
219  rb_system.add_scaled_Aq(scalar, q_a, matrix_A.get(), true);
220 }
221 
223 {
224  // Load the operators from the RBConstruction
225  EquationSystems & es = this->get_equation_systems();
227 
228  matrix_B->zero();
229  matrix_B->close();
230  matrix_B->add(1.,*rb_system.get_inner_product_matrix());
231 }
232 
234 {
235  LOG_SCOPE("perform_SCM_greedy()", "RBSCMConstruction");
236 
237  // initialize rb_scm_eval's parameters
239 
240  // Get a list of constrained dofs from rb_system
241  std::set<dof_id_type> constrained_dofs_set;
242  EquationSystems & es = this->get_equation_systems();
244 
245  for (dof_id_type i=0; i<rb_system.n_dofs(); i++)
246  {
247  if (rb_system.get_dof_map().is_constrained_dof(i))
248  {
249  constrained_dofs_set.insert(i);
250  }
251  }
252 
253  // Use these constrained dofs to identify which dofs we want to "get rid of"
254  // (i.e. condense) in our eigenproblems.
255  this->initialize_condensed_dofs(constrained_dofs_set);
256 
257  // Copy the inner product matrix over from rb_system to be used as matrix_B
258  load_matrix_B();
259 
261 
263  // This loads the new parameter into current_parameters
264  enrich_C_J(0);
265 
266  unsigned int SCM_iter=0;
267  while (true)
268  {
269  // matrix_A is reinitialized for the current parameters
270  // on each call to evaluate_stability_constant
272 
273  std::pair<unsigned int,Real> SCM_error_pair = compute_SCM_bounds_on_training_set();
274 
275  libMesh::out << "SCM iteration " << SCM_iter
276  << ", max_SCM_error = " << SCM_error_pair.second << std::endl;
277 
278  if (SCM_error_pair.second < SCM_training_tolerance)
279  {
280  libMesh::out << std::endl << "SCM tolerance of " << SCM_training_tolerance << " reached."
281  << std::endl << std::endl;
282  break;
283  }
284 
285  // If we need another SCM iteration, then enrich C_J
286  enrich_C_J(SCM_error_pair.first);
287 
288  libMesh::out << std::endl << "-----------------------------------" << std::endl << std::endl;
289 
290  SCM_iter++;
291  }
292 }
293 
295 {
296  LOG_SCOPE("compute_SCM_bounding_box()", "RBSCMConstruction");
297 
298  // Resize the bounding box vectors
299  rb_scm_eval->B_min.resize(get_rb_theta_expansion().get_n_A_terms());
300  rb_scm_eval->B_max.resize(get_rb_theta_expansion().get_n_A_terms());
301 
302  for (unsigned int q=0; q<get_rb_theta_expansion().get_n_A_terms(); q++)
303  {
304  matrix_A->zero();
305  add_scaled_symm_Aq(q, 1.);
306 
307  // Compute B_min(q)
308  eigen_solver->set_position_of_spectrum(SMALLEST_REAL);
310 
311  solve();
312  // TODO: Assert convergence for eigensolver
313 
314  unsigned int nconv = get_n_converged();
315  if (nconv != 0)
316  {
317  std::pair<Real, Real> eval = get_eigenpair(0);
318 
319  // ensure that the eigenvalue is real
320  libmesh_assert_less (eval.second, TOLERANCE);
321 
322  rb_scm_eval->set_B_min(q, eval.first);
323  libMesh::out << std::endl << "B_min("<<q<<") = " << rb_scm_eval->get_B_min(q) << std::endl;
324  }
325  else
326  libmesh_error_msg("Eigen solver for computing B_min did not converge");
327 
328  // Compute B_max(q)
329  eigen_solver->set_position_of_spectrum(LARGEST_REAL);
331 
332  solve();
333  // TODO: Assert convergence for eigensolver
334 
335  nconv = get_n_converged();
336  if (nconv != 0)
337  {
338  std::pair<Real, Real> eval = get_eigenpair(0);
339 
340  // ensure that the eigenvalue is real
341  libmesh_assert_less (eval.second, TOLERANCE);
342 
343  rb_scm_eval->set_B_max(q,eval.first);
344  libMesh::out << "B_max("<<q<<") = " << rb_scm_eval->get_B_max(q) << std::endl;
345  }
346  else
347  libmesh_error_msg("Eigen solver for computing B_max did not converge");
348  }
349 }
350 
352 {
353  LOG_SCOPE("evaluate_stability_constant()", "RBSCMConstruction");
354 
355  // Get current index of C_J
356  const unsigned int j = rb_scm_eval->C_J.size()-1;
357 
358  eigen_solver->set_position_of_spectrum(SMALLEST_REAL);
359 
360  // We assume B is set in system assembly
361  // For coercive problems, B is set to the inner product matrix
362  // For non-coercive time-dependent problems, B is set to the mass matrix
363 
364  // Set matrix A corresponding to mu_star
365  matrix_A->zero();
366  for (unsigned int q=0; q<get_rb_theta_expansion().get_n_A_terms(); q++)
367  {
369  }
370 
372  solve();
373  // TODO: Assert convergence for eigensolver
374 
375  unsigned int nconv = get_n_converged();
376  if (nconv != 0)
377  {
378  std::pair<Real, Real> eval = get_eigenpair(0);
379 
380  // ensure that the eigenvalue is real
381  libmesh_assert_less (eval.second, TOLERANCE);
382 
383  // Store the coercivity constant corresponding to mu_star
385  libMesh::out << std::endl << "Stability constant for C_J("<<j<<") = "
386  << rb_scm_eval->get_C_J_stability_constraint(j) << std::endl << std::endl;
387 
388  // Compute and store the vector y = (y_1, \ldots, y_Q) for the
389  // eigenvector currently stored in eigen_system.solution.
390  // We use this later to compute the SCM upper bounds.
392 
393  for (unsigned int q=0; q<get_rb_theta_expansion().get_n_A_terms(); q++)
394  {
395  Real norm_Aq2 = libmesh_real( Aq_inner_product(q, *solution, *solution) );
396 
397  rb_scm_eval->set_SCM_UB_vector(j,q,norm_Aq2/norm_B2);
398  }
399  }
400  else
401  libmesh_error_msg("Error: Eigensolver did not converge in evaluate_stability_constant");
402 }
403 
405  const NumericVector<Number> & w) const
406 {
407  matrix_B->vector_mult(*inner_product_storage_vector, w);
408 
410 }
411 
413  const NumericVector<Number> & v,
414  const NumericVector<Number> & w)
415 {
416  if (q >= get_rb_theta_expansion().get_n_A_terms())
417  libmesh_error_msg("Error: We must have q < Q_a in Aq_inner_product.");
418 
419  matrix_A->zero();
420  add_scaled_symm_Aq(q, 1.);
421  matrix_A->vector_mult(*inner_product_storage_vector, w);
422 
424 }
425 
427 {
428  LOG_SCOPE("compute_SCM_bounds_on_training_set()", "RBSCMConstruction");
429 
430  // Now compute the maximum bound error over training_parameters
431  unsigned int new_C_J_index = 0;
432  Real max_SCM_error = 0.;
433 
435  for (unsigned int i=0; i<get_local_n_training_samples(); i++)
436  {
437  set_params_from_training_set(first_index+i);
439  Real LB = rb_scm_eval->get_SCM_LB();
440  Real UB = rb_scm_eval->get_SCM_UB();
441 
442  Real error_i = SCM_greedy_error_indicator(LB, UB);
443 
444  if (error_i > max_SCM_error)
445  {
446  max_SCM_error = error_i;
447  new_C_J_index = i;
448  }
449  }
450 
451  numeric_index_type global_index = first_index + new_C_J_index;
452  std::pair<numeric_index_type,Real> error_pair(global_index, max_SCM_error);
453  get_global_max_error_pair(this->comm(),error_pair);
454 
455  return error_pair;
456 }
457 
458 void RBSCMConstruction::enrich_C_J(unsigned int new_C_J_index)
459 {
460  LOG_SCOPE("enrich_C_J()", "RBSCMConstruction");
461 
463 
464  rb_scm_eval->C_J.push_back(get_parameters());
465 
466  libMesh::out << std::endl << "SCM: Added mu = (";
467 
470  for ( ; it != it_end; ++it)
471  {
472  if (it != get_parameters().begin())
473  libMesh::out << ",";
474  std::string param_name = it->first;
475  RBParameters C_J_params = rb_scm_eval->C_J[rb_scm_eval->C_J.size()-1];
476  libMesh::out << C_J_params.get_value(param_name);
477  }
478  libMesh::out << ")" << std::endl;
479 
480  // Finally, resize C_J_stability_vector and SCM_UB_vectors
481  rb_scm_eval->C_J_stability_vector.push_back(0.);
482 
483  std::vector<Real> zero_vector(get_rb_theta_expansion().get_n_A_terms());
484  rb_scm_eval->SCM_UB_vectors.push_back(zero_vector);
485 }
486 
487 
488 } // namespace libMesh
489 
490 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
T libmesh_real(T a)
std::vector< Real > B_min
B_min, B_max define the bounding box.
numeric_index_type get_local_n_training_samples() const
Get the total number of training samples local to this processor.
const RBParameters & get_parameters_max() const
Get an RBParameters object that specifies the maximum allowable value for each parameter.
virtual Real get_SCM_UB()
Evaluate single SCM upper bound.
This is the EquationSystems class.
Number Aq_inner_product(unsigned int q, const NumericVector< Number > &v, const NumericVector< Number > &w)
Compute the inner product between two vectors using matrix Aq.
Real get_B_min(unsigned int i) const
Get B_min and B_max.
virtual void clear() libmesh_override
Clear all the data structures associated with the system.
UniquePtr< NumericVector< Number > > inner_product_storage_vector
We keep an extra temporary vector that is useful for performing inner products (avoids unnecessary me...
void set_SCM_UB_vector(unsigned int j, unsigned int q, Real y_q)
Set entries of SCM_UB_vector, which stores the vector y, corresponding to the minimizing eigenvectors...
unsigned int get_n_A_terms() const
Get Q_a, the number of terms in the affine expansion for the bilinear form.
Real get_value(const std::string &param_name) const
Get the value of the specific parameter.
Definition: rb_parameters.C:45
void print_discrete_parameter_values() const
Print out all the discrete parameter values.
unsigned int get_n_params() const
Get the number of parameters.
UniquePtr< SparseMatrix< Number > > matrix_B
A second system matrix for generalized eigenvalue problems.
Definition: eigen_system.h:155
virtual void solve() libmesh_override
Override to solve the condensed eigenproblem with the dofs in local_non_condensed_dofs_vector strippe...
virtual void enrich_C_J(unsigned int new_C_J_index)
Enrich C_J by adding the element of SCM_training_samples that has the largest gap between alpha_LB an...
const class libmesh_nullptr_t libmesh_nullptr
unsigned int get_n_converged() const
Definition: eigen_system.h:124
virtual void evaluate_stability_constant()
Compute the stability constant for current_parameters by solving a generalized eigenvalue problem ove...
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.
static const Real TOLERANCE
static void get_global_max_error_pair(const Parallel::Communicator &communicator, std::pair< numeric_index_type, Real > &error_pair)
Static function to return the error pair (index,error) that is corresponds to the largest error on al...
virtual ~RBSCMConstruction()
Destructor.
The libMesh namespace provides an interface to certain functionality in the library.
Real get_B_max(unsigned int i) const
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
std::vector< Real > B_max
RBSCMEvaluation * rb_scm_eval
The current RBSCMEvaluation object we are using to perform the Evaluation stage of the SCM...
UniquePtr< SparseMatrix< Number > > matrix_A
The system matrix for standard eigenvalue problems.
Definition: eigen_system.h:150
const std::string & name() const
Definition: system.h:1998
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
std::vector< std::vector< Real > > SCM_UB_vectors
This matrix stores the infimizing vectors y_1( ),...,y_Q_a( ), for each in C_J, which are used in co...
virtual void set_params_from_training_set_and_broadcast(unsigned int index)
Load the specified training parameter and then broadcast to all processors.
std::vector< Real > C_J_stability_vector
Vector storing the (truth) stability values at the parameters in C_J.
Real SCM_training_tolerance
Tolerance which controls when to terminate the SCM Greedy.
numeric_index_type get_n_training_samples() const
Get the total number of training samples.
virtual Real get_SCM_LB()
Evaluate single SCM lower bound.
virtual void print_info()
Print out info that describes the current setup of this RBSCMConstruction.
void set_rb_scm_evaluation(RBSCMEvaluation &rb_scm_eval_in)
Set the RBSCMEvaluation object.
Real get_C_J_stability_constraint(unsigned int j) const
Get stability constraints (i.e.
bool is_constrained_dof(const dof_id_type dof) const
Definition: dof_map.h:1734
virtual void load_matrix_B()
Copy over the matrix to store in matrix_B, usually this is the mass or inner-product matrix...
const DofMap & get_dof_map() const
Definition: system.h:2030
Real get_parameter_max(const std::string &param_name) const
Get maximum allowable value of parameter param_name.
dof_id_type numeric_index_type
Definition: id_types.h:92
std::string RB_system_name
The name of the associated RB system.
void set_parameters(const RBParameters &params)
Set the current parameters to params.
std::vector< RBParameters > C_J
Vector storing the greedily selected parameters during SCM training.
Number B_inner_product(const NumericVector< Number > &v, const NumericVector< Number > &w) const
Compute the inner product between two vectors using the system&#39;s matrix_B.
void set_params_from_training_set(unsigned int index)
Set parameters to the RBParameters stored in index index of the training set.
void set_B_min(unsigned int i, Real B_min_val)
Set B_min and B_max.
SparseMatrix< Number > * get_inner_product_matrix()
Get a pointer to inner_product_matrix.
void set_C_J_stability_constraint(unsigned int j, Real stability_constraint_in)
Set stability constraints (i.e.
virtual Real SCM_greedy_error_indicator(Real LB, Real UB)
Helper function which provides an error indicator to be used in the SCM greedy.
virtual void attach_deflation_space()
Attach the deflation space defined by the specified vector, can be useful in solving constrained eige...
const_iterator end() const
Get a constant iterator to the end of this RBParameters object.
Definition: rb_parameters.C:85
virtual void set_eigensolver_properties(int)
This function is called before truth eigensolves in compute_SCM_bounding_box and evaluate_stability_c...
RBSCMConstruction(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Constructor.
virtual void initialize_training_parameters(const RBParameters &mu_min, const RBParameters &mu_max, unsigned int n_training_parameters, std::map< std::string, bool > log_param_scale, bool deterministic=true)
Initialize the parameter ranges and indicate whether deterministic or random training parameters shou...
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
void initialize_condensed_dofs(const std::set< dof_id_type > &global_condensed_dofs_set=std::set< dof_id_type >())
Loop over the dofs on each processor to initialize the list of non-condensed dofs.
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:42
const RBParameters & get_parameters_min() const
Get an RBParameters object that specifies the minimum allowable value for each parameter.
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the RBThetaExpansion object.
virtual std::pair< unsigned int, Real > compute_SCM_bounds_on_training_set()
Compute upper and lower bounds for each SCM training point.
virtual void process_parameters_file(const std::string &parameters_filename)
Read in the parameters from file specified by parameters_filename and set the this system&#39;s member va...
virtual void resize_SCM_vectors()
Clear and resize the SCM data vectors.
void set_SCM_training_tolerance(Real SCM_training_tolerance_in)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
numeric_index_type get_first_local_training_index() const
Get the first local index of the training parameters.
void set_training_random_seed(unsigned int seed)
Set the seed that is used to randomly generate training parameters.
virtual void perform_SCM_greedy()
Perform the SCM greedy algorithm to develop a lower bound over the training set.
void set_B_max(unsigned int i, Real B_max_val)
virtual void add_scaled_symm_Aq(unsigned int q_a, Number scalar)
Add the scaled symmetrized affine matrix from the associated RBSystem to matrix_A.
UniquePtr< EigenSolver< Number > > eigen_solver
The EigenSolver, defining which interface, i.e solver package to use.
Definition: eigen_system.h:161
void set_value(const std::string &param_name, Real value)
Set the value of the specified parameter.
Definition: rb_parameters.C:57
virtual void compute_SCM_bounding_box()
Compute the SCM bounding box.
const Parallel::Communicator & comm() const
OStreamProxy out
RBSCMEvaluation & get_rb_scm_evaluation()
Get a reference to the RBSCMEvaluation object.
const EquationSystems & get_equation_systems() const
Definition: system.h:712
void set_eigenproblem_type(EigenProblemType ept)
Sets the type of the current eigen problem.
Definition: eigen_system.C:78
virtual T dot(const NumericVector< T > &v) const =0
This class is part of the rbOOmit framework.
dof_id_type n_dofs() const
Definition: system.C:148
const T_sys & get_system(const std::string &name) const
const_iterator begin() const
Get a constant iterator to beginning of this RBParameters object.
Definition: rb_parameters.C:80
const RBParameters & get_parameters() const
Get the current parameters.
bool assemble_before_solve
Flag which tells the system to whether or not to call the user assembly function during each call to ...
Definition: system.h:1477
std::map< std::string, Real >::const_iterator const_iterator
Definition: rb_parameters.h:57
void initialize_parameters(const RBParameters &mu_min_in, const RBParameters &mu_max_in, const std::map< std::string, std::vector< Real >> &discrete_parameter_values)
Initialize the parameter ranges and set current_parameters.
virtual std::pair< Real, Real > get_eigenpair(dof_id_type i) libmesh_override
Override get_eigenpair() to retrieve the eigenpair for the condensed eigensolve.
virtual void clear()
Clear all the data structures associated with the system.
Real get_SCM_training_tolerance() const
Get/set SCM_training_tolerance: tolerance for SCM greedy.
This class is part of the rbOOmit framework.
uint8_t dof_id_type
Definition: id_types.h:64
Real get_parameter_min(const std::string &param_name) const
Get minimum allowable value of parameter param_name.