libMesh
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Attributes | List of all members
libMesh::RBConstructionBase< Base > Class Template Reference

This class is part of the rbOOmit framework. More...

#include <rb_construction_base.h>

Inheritance diagram for libMesh::RBConstructionBase< Base >:
[legend]

Public Types

typedef RBConstructionBase< Base > sys_type
 The type of system. More...
 

Public Member Functions

 RBConstructionBase (EquationSystems &es, const std::string &name, const unsigned int number)
 Constructor. More...
 
 RBConstructionBase (RBConstructionBase &&)=default
 Special functions. More...
 
RBConstructionBaseoperator= (RBConstructionBase &&)=delete
 
 RBConstructionBase (const RBConstructionBase &)=delete
 
RBConstructionBaseoperator= (const RBConstructionBase &)=delete
 
virtual ~RBConstructionBase ()
 
sys_typesystem ()
 
virtual void clear ()
 Clear all the data structures associated with the system. More...
 
void set_quiet_mode (bool quiet_mode_in)
 Set the quiet_mode flag. More...
 
bool is_quiet () const
 Is the system in quiet mode? More...
 
numeric_index_type get_n_training_samples () const
 Get the number of global training samples. More...
 
numeric_index_type get_local_n_training_samples () const
 Get the total number of training samples local to this processor. More...
 
numeric_index_type get_first_local_training_index () const
 Get the first local index of the training parameters. More...
 
numeric_index_type get_last_local_training_index () const
 Get the last local index of the training parameters. More...
 
virtual void initialize_training_parameters (const RBParameters &mu_min, const RBParameters &mu_max, const unsigned int n_global_training_samples, const std::map< std::string, bool > &log_param_scale, const bool deterministic=true)
 Initialize the parameter ranges and indicate whether deterministic or random training parameters should be used and whether or not we want the parameters to be scaled logarithmically. More...
 
virtual void load_training_set (const std::map< std::string, std::vector< RBParameter >> &new_training_set)
 Overwrite the training parameters with new_training_set. More...
 
void set_training_parameter_values (const std::string &param_name, const std::vector< RBParameter > &values)
 Overwrite the local training samples for param_name using values. More...
 
void broadcast_parameters (const unsigned int proc_id)
 Broadcasts parameters from processor proc_id to all processors. More...
 
void set_training_random_seed (int seed)
 Set the seed that is used to randomly generate training parameters. More...
 
void set_deterministic_training_parameter_name (const std::string &name)
 In some cases we only want to allow discrete parameter values, instead of parameters that may take any value in a specified interval. More...
 
const std::string & get_deterministic_training_parameter_name () const
 Get the name of the parameter that we will generate deterministic training parameters for. More...
 
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. More...
 
void initialize_parameters (const RBParametrized &rb_parametrized)
 Initialize the parameter ranges and set current_parameters. More...
 
unsigned int get_n_params () const
 Get the number of parameters. More...
 
unsigned int get_n_continuous_params () const
 Get the number of continuous parameters. More...
 
unsigned int get_n_discrete_params () const
 Get the number of discrete parameters. More...
 
std::set< std::string > get_parameter_names () const
 Get a set that stores the parameter names. More...
 
const RBParametersget_parameters () const
 Get the current parameters. More...
 
bool set_parameters (const RBParameters &params)
 Set the current parameters to params The parameters are checked for validity; an error is thrown if the number of parameters or samples is different than expected. More...
 
const RBParametersget_parameters_min () const
 Get an RBParameters object that specifies the minimum allowable value for each parameter. More...
 
const RBParametersget_parameters_max () const
 Get an RBParameters object that specifies the maximum allowable value for each parameter. More...
 
Real get_parameter_min (const std::string &param_name) const
 Get minimum allowable value of parameter param_name. More...
 
Real get_parameter_max (const std::string &param_name) const
 Get maximum allowable value of parameter param_name. More...
 
void print_parameters () const
 Print the current parameters. More...
 
void write_parameter_data_to_files (const std::string &continuous_param_file_name, const std::string &discrete_param_file_name, const bool write_binary_data)
 Write out the parameter ranges to files. More...
 
void read_parameter_data_from_files (const std::string &continuous_param_file_name, const std::string &discrete_param_file_name, const bool read_binary_data)
 Read in the parameter ranges from files. More...
 
bool is_discrete_parameter (const std::string &mu_name) const
 Is parameter mu_name discrete? More...
 
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values () const
 Get a const reference to the discrete parameter values. More...
 
void print_discrete_parameter_values () const
 Print out all the discrete parameter values. More...
 

Static Public Member Functions

static std::pair< std::size_t, std::size_t > generate_training_parameters_random (const Parallel::Communicator &communicator, const std::map< std::string, bool > &log_param_scale, std::map< std::string, std::vector< RBParameter >> &local_training_parameters_in, const unsigned int n_global_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, const int training_parameters_random_seed=-1, const bool serial_training_set=false)
 Static helper function for generating a randomized set of parameters. More...
 
static std::pair< std::size_t, std::size_t > generate_training_parameters_deterministic (const Parallel::Communicator &communicator, const std::map< std::string, bool > &log_param_scale, std::map< std::string, std::vector< RBParameter >> &local_training_parameters_in, const unsigned int n_global_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, const bool serial_training_set=false)
 Static helper function for generating a deterministic set of parameters. More...
 
static Real get_closest_value (Real value, const std::vector< Real > &list_of_values)
 
static std::string get_info ()
 Gets a string containing the reference information. More...
 
static void print_info (std::ostream &out_stream=libMesh::out)
 Prints the reference information, by default to libMesh::out. More...
 
static unsigned int n_objects ()
 Prints the number of outstanding (created, but not yet destroyed) objects. More...
 
static void enable_print_counter_info ()
 Methods to enable/disable the reference counter output from print_info() More...
 
static void disable_print_counter_info ()
 

Public Attributes

bool verbose_mode
 Public boolean to toggle verbose mode. More...
 

Protected Types

typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 Data structure to log the information. More...
 

Protected Member Functions

virtual void init_data ()
 Initializes the member data fields associated with the system, so that, e.g., assemble() may be used. More...
 
RBParameters get_params_from_training_set (unsigned int global_index)
 Return the RBParameters in index global_index of the global training set. More...
 
void set_params_from_training_set (unsigned int global_index)
 Set parameters to the RBParameters stored in index global_index of the global training set. More...
 
virtual void set_params_from_training_set_and_broadcast (unsigned int global_index)
 Load the specified training parameter and then broadcast to all processors. More...
 
void increment_constructor_count (const std::string &name) noexcept
 Increments the construction counter. More...
 
void increment_destructor_count (const std::string &name) noexcept
 Increments the destruction counter. More...
 

Static Protected Member Functions

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 all processors. More...
 

Protected Attributes

bool quiet_mode
 Flag to indicate whether we print out extra information during the Offline stage. More...
 
bool serial_training_set
 This boolean flag indicates whether or not the training set should be the same on all processors. More...
 
std::unique_ptr< NumericVector< Number > > inner_product_storage_vector
 We keep an extra temporary vector that is useful for performing inner products (avoids unnecessary memory allocation/deallocation). More...
 

Static Protected Attributes

static Counts _counts
 Actually holds the data. More...
 
static Threads::atomic< unsigned int_n_objects
 The number of objects. More...
 
static Threads::spin_mutex _mutex
 Mutual exclusion object to enable thread-safe reference counting. More...
 
static bool _enable_print_counter = true
 Flag to control whether reference count information is printed when print_info is called. More...
 

Private Attributes

bool _training_parameters_initialized
 Boolean flag to indicate whether or not the parameter ranges have been initialized. More...
 
std::map< std::string, std::vector< RBParameter > > _training_parameters
 The training samples for each parameter. More...
 
numeric_index_type _first_local_index
 The first sample-vector index from the global vector which is stored in the _training_parameters on this processor. More...
 
numeric_index_type _n_local_training_samples
 
numeric_index_type _n_global_training_samples
 
int _training_parameters_random_seed
 If < 0, use std::time() * processor_id() to seed the random number generator for the training parameters (default). More...
 

Detailed Description

template<class Base>
class libMesh::RBConstructionBase< Base >

This class is part of the rbOOmit framework.

This is the base class for the Construction stage of the certified reduced basis (RB) method. We template the Base class so that we can derive from the appropriate libMesh System type (e.g. LinearImplicitSystem for standard reduced basis, EigenSystem for SCM) at compile time.

Author
David J. Knezevic
Date
2009

Definition at line 55 of file rb_construction_base.h.

Member Typedef Documentation

◆ Counts

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information.

The log is identified by the class name.

Definition at line 119 of file reference_counter.h.

◆ sys_type

template<class Base>
typedef RBConstructionBase<Base> libMesh::RBConstructionBase< Base >::sys_type

The type of system.

Definition at line 82 of file rb_construction_base.h.

Constructor & Destructor Documentation

◆ RBConstructionBase() [1/3]

template<class Base >
libMesh::RBConstructionBase< Base >::RBConstructionBase ( EquationSystems es,
const std::string &  name,
const unsigned int  number 
)

Constructor.

Definition at line 95 of file rb_construction_base.C.

98  : Base(es, name_in, number_in),
99  quiet_mode(true),
100  serial_training_set(false),
105  _training_parameters_random_seed(-1) // by default, use std::time to seed RNG
106 {
107 }
bool quiet_mode
Flag to indicate whether we print out extra information during the Offline stage. ...
numeric_index_type _first_local_index
The first sample-vector index from the global vector which is stored in the _training_parameters on t...
numeric_index_type _n_local_training_samples
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...
int _training_parameters_random_seed
If < 0, use std::time() * processor_id() to seed the random number generator for the training paramet...
numeric_index_type _n_global_training_samples
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ RBConstructionBase() [2/3]

template<class Base>
libMesh::RBConstructionBase< Base >::RBConstructionBase ( RBConstructionBase< Base > &&  )
default

Special functions.

  • This class has the same restrictions as the union of its potential base classes (currently LinearImplicitSystem and EigenSystem).
  • Destructor is defaulted out-of-line.

◆ RBConstructionBase() [3/3]

template<class Base>
libMesh::RBConstructionBase< Base >::RBConstructionBase ( const RBConstructionBase< Base > &  )
delete

◆ ~RBConstructionBase()

template<class Base >
libMesh::RBConstructionBase< Base >::~RBConstructionBase ( )
virtualdefault

Member Function Documentation

◆ broadcast_parameters()

template<class Base >
void libMesh::RBConstructionBase< Base >::broadcast_parameters ( const unsigned int  proc_id)

Broadcasts parameters from processor proc_id to all processors.

This broadcasts the RBParameters object from .get_parameters(), and then sets it on all processors with .set_parameters().

Definition at line 712 of file rb_construction_base.C.

713 {
714  libmesh_assert_less (proc_id, this->n_processors());
715 
716  // create a copy of the current parameters
717  RBParameters current_parameters = get_parameters();
718  libmesh_error_msg_if(current_parameters.n_samples()!=1,
719  "Only single-sample RBParameter objects can be broadcast.");
720 
721  // Serialize the current_parameters to current_parameters_vector in order to broadcast.
722  // We handle multiple samples and vector values.
723  // However, the vector values are assumed to remain the same size across samples.
724  const std::size_t nparams = current_parameters.n_parameters();
725  const std::size_t nsamples = current_parameters.n_samples();
726 
727  // First we get the sizes of all the parameter value vectors.
728  std::vector<std::size_t> param_value_sizes;
729  param_value_sizes.reserve(nparams);
730  for (const auto & pr : current_parameters)
731  param_value_sizes.push_back(pr.second[0].size());
732 
733  // Broadcast the sizes vector and reserve memory.
734  this->comm().broadcast(param_value_sizes, proc_id);
735  std::size_t buffsize = std::accumulate(param_value_sizes.cbegin(), param_value_sizes.cend(), 0ul);
736  std::vector<Real> serialized_parameters;
737  serialized_parameters.reserve(buffsize);
738 
739  // Then we serialize the parameters/sample/value vectors into a single vector.
740  for (const auto & pr : current_parameters)
741  {
742  for (const auto sample_idx : make_range(nsamples))
743  serialized_parameters.insert(serialized_parameters.end(),
744  pr.second[sample_idx].cbegin(),
745  pr.second[sample_idx].cend());
746  }
747 
748  // Do the broadcasts.
749  this->comm().broadcast(serialized_parameters, proc_id);
750 
751  // Deserialize into the copy of the RBParameters object.
752  std::size_t param_idx = 0;
753  auto val_idx = serialized_parameters.cbegin();
754  for (const auto & pr : current_parameters)
755  {
756  const std::size_t param_value_size = param_value_sizes[param_idx];
757  for (const auto sample_idx: make_range(nsamples))
758  {
759  auto end_val_idx = std::next(val_idx,param_value_size);
760  RBParameter sample_val(val_idx, end_val_idx);
761  current_parameters.set_value(pr.first, sample_idx, sample_val);
762  val_idx = end_val_idx;
763  }
764  ++param_idx;
765  }
766 
767  // Overwrite the parameters globally.
768  set_parameters(current_parameters);
769 }
std::vector< Real > RBParameter
Typedef for an individual RB parameter.
Definition: rb_parameters.h:39
const RBParameters & get_parameters() const
Get the current parameters.
bool set_parameters(const RBParameters &params)
Set the current parameters to params The parameters are checked for validity; an error is thrown if t...
static const unsigned int next[3]
A lookup table for the increment modulo 3 operation, for iterating through the three nodes per elemen...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134

◆ clear()

template<class Base >
void libMesh::RBConstructionBase< Base >::clear ( )
virtual

Clear all the data structures associated with the system.

Reimplemented from libMesh::RBParametrized.

Reimplemented in libMesh::RBConstruction, libMesh::RBEIMConstruction, libMesh::TransientSystem< RBConstruction >, libMesh::RBSCMConstruction, and libMesh::TransientRBConstruction.

Definition at line 113 of file rb_construction_base.C.

Referenced by libMesh::RBEIMConstruction::clear().

114 {
115  // clear the parent data
116  Base::clear();
118  _training_parameters.clear();
119 }
virtual void clear()
Clear all the data structures associated with the system.
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

101 {
102  _enable_print_counter = false;
103  return;
104 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 94 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

95 {
96  _enable_print_counter = true;
97  return;
98 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ generate_training_parameters_deterministic()

template<class Base >
std::pair< std::size_t, std::size_t > libMesh::RBConstructionBase< Base >::generate_training_parameters_deterministic ( const Parallel::Communicator communicator,
const std::map< std::string, bool > &  log_param_scale,
std::map< std::string, std::vector< RBParameter >> &  local_training_parameters_in,
const unsigned int  n_global_training_samples_in,
const RBParameters min_parameters,
const RBParameters max_parameters,
const bool  serial_training_set = false 
)
static

Static helper function for generating a deterministic set of parameters.

Only works with 1 or 2 parameters (as defined by the lengths of min/max parameters vectors), otherwise throws an error. The parameter n_global_training_samples_in is the total number of parameters to generate, and they will be split across all the processors (unless serial_training_set=true) in the local_training_parameters_in map.

Returns
a pair of {first_local_index,last_local_index}

Definition at line 564 of file rb_construction_base.C.

571 {
572  libmesh_assert_equal_to ( min_parameters.n_parameters(), max_parameters.n_parameters() );
573  const unsigned int num_params = min_parameters.n_parameters();
574 
575  if (num_params == 0)
576  return {0,0};
577 
578  if (num_params > 3)
579  libmesh_not_implemented_msg("ERROR: Deterministic training sample generation "
580  "not implemented for more than three parameters.");
581 
582  // TODO - we don't support vector-data here yet. This would only apply in the case where
583  // min or max are vector-valued, and all the generated points need to stay within those ranges.
584  // But typically we expect that if we're calling this function, we only have 1 min and 1 max,
585  // so the generated values are single-valued as well. The .get_value() calls will throw an error
586  // if this is not the case.
587 
588  // Reinitialize training_parameters_in (but don't remove existing keys!)
589  const auto &[n_local_training_samples, first_local_index] =
590  calculate_n_local_samples_and_index(communicator, n_global_training_samples_in,
592  const auto last_local_index = first_local_index + n_local_training_samples;
593  for (const auto & pr : min_parameters)
594  local_training_parameters_in[pr.first] = std::vector<RBParameter>(n_local_training_samples);
595 
596  // n_training_samples_per_param has 3 entries, but entries after "num_params"
597  // are unused so we just set their value to 1. We need to set it to 1 (rather
598  // than 0) so that we don't skip the inner part of the triply-nested loop over
599  // n_training_samples_per_param below.
600  std::vector<unsigned int> n_training_samples_per_param(3);
601  for (unsigned int param=0; param<3; param++)
602  {
603  if (param < num_params)
604  {
605  n_training_samples_per_param[param] =
606  static_cast<unsigned int>( std::round(std::pow(static_cast<Real>(n_global_training_samples_in), 1./num_params)) );
607  }
608  else
609  {
610  n_training_samples_per_param[param] = 1;
611  }
612  }
613 
614  {
615  // The current implementation assumes that we have the same number of
616  // samples in each parameter, so we check that n_training_samples_in
617  // is consistent with this assumption.
618  unsigned int total_samples_check = 1;
619  for (unsigned int n_samples : n_training_samples_per_param)
620  {
621  total_samples_check *= n_samples;
622  }
623 
624  libmesh_error_msg_if(total_samples_check != n_global_training_samples_in,
625  "Error: Number of training samples = "
626  << n_global_training_samples_in
627  << " does not enable a uniform grid of samples with "
628  << num_params << " parameters. Try "
629  << total_samples_check << " samples instead?");
630  }
631 
632  // First we make a list of training samples associated with each parameter,
633  // then we take a tensor product to obtain the final set of training samples.
634  std::vector<std::vector<Real>> training_samples_per_param(num_params);
635  {
636  unsigned int i = 0;
637  for (const auto & pr : min_parameters)
638  {
639  const std::string & param_name = pr.first;
640  const bool use_log_scaling = libmesh_map_find(log_param_scale, param_name);
641  Real min_param = min_parameters.get_value(param_name);
642  Real max_param = max_parameters.get_value(param_name);
643 
644  training_samples_per_param[i].resize(n_training_samples_per_param[i]);
645 
646  for (unsigned int j=0; j<n_training_samples_per_param[i]; j++)
647  {
648  // Generate log10 scaled training parameters
649  if (use_log_scaling)
650  {
651  Real epsilon = 1.e-6; // Prevent rounding errors triggering asserts
652  Real log_min = std::log10(min_param + epsilon);
653  Real log_range = std::log10( (max_param-epsilon) / (min_param+epsilon) );
654  Real step_size = log_range /
655  std::max((unsigned int)1,(n_training_samples_per_param[i]-1));
656 
657  if (j<(n_training_samples_per_param[i]-1))
658  {
659  training_samples_per_param[i][j] = std::pow(10., log_min + j*step_size );
660  }
661  else
662  {
663  // due to rounding error, the last parameter can be slightly
664  // bigger than max_parameters, hence snap back to the max
665  training_samples_per_param[i][j] = max_param;
666  }
667  }
668  else
669  {
670  // Generate linearly scaled training parameters
671  Real step_size = (max_param - min_param) /
672  std::max((unsigned int)1,(n_training_samples_per_param[i]-1));
673  training_samples_per_param[i][j] = j*step_size + min_param;
674  }
675 
676  }
677  i++;
678  }
679  }
680 
681  // Now load into training_samples_in
682  {
683  std::vector<unsigned int> indices(3);
684  unsigned int index_count = 0;
685  for (indices[0]=0; indices[0]<n_training_samples_per_param[0]; indices[0]++)
686  {
687  for (indices[1]=0; indices[1]<n_training_samples_per_param[1]; indices[1]++)
688  {
689  for (indices[2]=0; indices[2]<n_training_samples_per_param[2]; indices[2]++)
690  {
691  unsigned int param_count = 0;
692  for (const auto & pr : min_parameters)
693  {
694  std::vector<RBParameter> & training_vector =
695  libmesh_map_find(local_training_parameters_in, pr.first);
696  if (first_local_index <= index_count && index_count < last_local_index)
697  training_vector[index_count - first_local_index] =
698  {training_samples_per_param[param_count][indices[param_count]]};
699 
700  param_count++;
701  }
702  index_count++;
703  }
704  }
705  }
706  }
707  return {first_local_index, first_local_index+n_local_training_samples};
708 }
T pow(const T &x)
Definition: utility.h:328
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...

◆ generate_training_parameters_random()

template<class Base >
std::pair< std::size_t, std::size_t > libMesh::RBConstructionBase< Base >::generate_training_parameters_random ( const Parallel::Communicator communicator,
const std::map< std::string, bool > &  log_param_scale,
std::map< std::string, std::vector< RBParameter >> &  local_training_parameters_in,
const unsigned int  n_global_training_samples_in,
const RBParameters min_parameters,
const RBParameters max_parameters,
const int  training_parameters_random_seed = -1,
const bool  serial_training_set = false 
)
static

Static helper function for generating a randomized set of parameters.

The parameter n_global_training_samples_in is the total number of parameters to generate, and they will be split across all the processors (unless serial_training_set=true) in the local_training_parameters_in map.

Returns
a pair of {first_local_index,last_local_index}

Definition at line 464 of file rb_construction_base.C.

472 {
473  const unsigned int num_params = min_parameters.n_parameters();
474  libmesh_error_msg_if(num_params!=max_parameters.n_parameters(),
475  "Number of parameters must be identical for min/max.");
476 
477  // Clear training_parameters_in
478  local_training_parameters_in.clear();
479 
480  if (num_params == 0)
481  return {0,0};
482 
483  if (training_parameters_random_seed < 0)
484  {
485  if (!serial_training_set)
486  {
487  // seed the random number generator with the system time
488  // and the processor ID so that the seed is different
489  // on different processors
490  std::srand( static_cast<unsigned>( std::time(0)*(1+communicator.rank()) ));
491  }
492  else
493  {
494  // seed the random number generator with the system time
495  // only so that the seed is the same on all processors
496  //
497  // Note that we broadcast the time on processor 0 to make
498  // sure all processors agree.
499  unsigned int current_time = static_cast<unsigned>( std::time(0) );
500  communicator.broadcast(current_time, 0);
501  std::srand(current_time);
502  }
503  }
504  else
505  {
506  if (!serial_training_set)
507  {
508  // seed the random number generator with the provided value
509  // and the processor ID so that the seed is different
510  // on different processors
511  std::srand( static_cast<unsigned>( training_parameters_random_seed*(1+communicator.rank()) ));
512  }
513  else
514  {
515  // seed the random number generator with the provided value
516  // so that the seed is the same on all processors
517  std::srand( static_cast<unsigned>( training_parameters_random_seed ));
518  }
519  }
520 
521  // TODO - we don't support vector-data here yet. This would only apply in the case where
522  // min or max are vector-valued, and all the generated points need to stay within those ranges.
523  // But typically we expect that if we're calling this function, we only have 1 min and 1 max,
524  // so the generated values are single-valued as well. The .get_value() calls will throw an error
525  // if this is not the case.
526 
527  // initialize training_parameters_in
528  const auto & [n_local_training_samples, first_local_index] =
529  calculate_n_local_samples_and_index(communicator, n_global_training_samples_in,
531  for (const auto & pr : min_parameters)
532  local_training_parameters_in[pr.first] = std::vector<RBParameter>(n_local_training_samples);
533 
534  // finally, set the values
535  for (auto & [param_name, sample_vector] : local_training_parameters_in)
536  {
537  for (auto i : make_range(n_local_training_samples))
538  {
539  Real random_number = static_cast<Real>(std::rand()) / RAND_MAX; // in range [0,1]
540 
541  // Generate log10 scaled training parameters
542  if (libmesh_map_find(log_param_scale, param_name))
543  {
544  Real log_min = std::log10(min_parameters.get_value(param_name));
545  Real log_range = std::log10(max_parameters.get_value(param_name) / min_parameters.get_value(param_name));
546 
547  sample_vector[i] = {std::pow(10., log_min + random_number*log_range )};
548  }
549  // Generate linearly scaled training parameters
550  else
551  {
552  sample_vector[i] = {
553  random_number * (max_parameters.get_value(param_name) -
554  min_parameters.get_value(param_name)) +
555  min_parameters.get_value(param_name)};
556  }
557  }
558  }
559  return {first_local_index, first_local_index+n_local_training_samples};
560 }
T pow(const T &x)
Definition: utility.h:328
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
DIE A HORRIBLE DEATH HERE typedef MPI_Comm communicator
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...

◆ get_closest_value()

Real libMesh::RBParametrized::get_closest_value ( Real  value,
const std::vector< Real > &  list_of_values 
)
staticinherited
Returns
The closest entry to value from list_of_values.

Definition at line 443 of file rb_parametrized.C.

References std::abs(), distance(), libMesh::Real, and value.

Referenced by libMesh::RBParametrized::is_value_in_list().

444 {
445  libmesh_error_msg_if(list_of_values.empty(), "Error: list_of_values is empty.");
446 
447  Real min_distance = std::numeric_limits<Real>::max();
448  Real closest_val = 0.;
449  for (const auto & current_value : list_of_values)
450  {
451  Real distance = std::abs(value - current_value);
452  if (distance < min_distance)
453  {
454  min_distance = distance;
455  closest_val = current_value;
456  }
457  }
458 
459  return closest_val;
460 }
Real distance(const Point &p)
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const bool value
Definition: xdr_io.C:54

◆ get_deterministic_training_parameter_name()

template<class Base>
const std::string& libMesh::RBConstructionBase< Base >::get_deterministic_training_parameter_name ( ) const

Get the name of the parameter that we will generate deterministic training parameters for.

◆ get_discrete_parameter_values()

const std::map< std::string, std::vector< Real > > & libMesh::RBParametrized::get_discrete_parameter_values ( ) const
inherited

Get a const reference to the discrete parameter values.

Definition at line 370 of file rb_parametrized.C.

References libMesh::RBParametrized::_discrete_parameter_values, and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), libMesh::RBParametrized::check_if_valid_params(), libMesh::RBParametrized::get_n_discrete_params(), libMesh::RBParametrized::initialize_parameters(), libMesh::RBParametrized::print_discrete_parameter_values(), and libMesh::RBParametrized::write_discrete_parameter_values_to_file().

371 {
372  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_discrete_parameter_values");
373 
375 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
std::map< std::string, std::vector< Real > > _discrete_parameter_values
Map that defines the allowable values of any discrete parameters.

◆ get_first_local_training_index()

template<class Base >
numeric_index_type libMesh::RBConstructionBase< Base >::get_first_local_training_index ( ) const

Get the first local index of the training parameters.

Definition at line 184 of file rb_construction_base.C.

185 {
186  libmesh_error_msg_if(!_training_parameters_initialized,
187  "Error: training parameters must first be initialized.");
188 
189  // First we check if there are no parameters here, and in that case we
190  // return 0 for a serial training set and comm().rank() for a parallel
191  // training set. This is consistent with get_n_training_samples(), and
192  // avoids accessing training_parameters.begin() when training_parameters
193  // is empty.
194  if (_training_parameters.empty())
195  {
197  return 0;
198  else
199  return this->comm().rank();
200  }
201 
202  return _first_local_index;
203 }
numeric_index_type _first_local_index
The first sample-vector index from the global vector which is stored in the _training_parameters on t...
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ get_global_max_error_pair()

template<class Base >
void libMesh::RBConstructionBase< Base >::get_global_max_error_pair ( const Parallel::Communicator communicator,
std::pair< numeric_index_type, Real > &  error_pair 
)
staticprotected

Static function to return the error pair (index,error) that is corresponds to the largest error on all processors.

Definition at line 133 of file rb_construction_base.C.

135 {
136  // Set error_pair.second to the maximum global value and also
137  // find which processor contains the maximum value
138  unsigned int proc_ID_index;
139  communicator.maxloc(error_pair.second, proc_ID_index);
140 
141  // Then broadcast error_pair.first from proc_ID_index
142  communicator.broadcast(error_pair.first, proc_ID_index);
143 }
DIE A HORRIBLE DEATH HERE typedef MPI_Comm communicator

◆ get_info()

std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (const auto & [name, cd] : _counts)
59  oss << "| " << name << " reference count information:\n"
60  << "| Creations: " << cd.first << '\n'
61  << "| Destructions: " << cd.second << '\n';
62 
63  oss << " ---------------------------------------------------------------------------- \n";
64 
65  return oss.str();
66 
67 #else
68 
69  return "";
70 
71 #endif
72 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
static Counts _counts
Actually holds the data.

◆ get_last_local_training_index()

template<class Base >
numeric_index_type libMesh::RBConstructionBase< Base >::get_last_local_training_index ( ) const

Get the last local index of the training parameters.

Definition at line 206 of file rb_construction_base.C.

207 {
208  libmesh_error_msg_if(!_training_parameters_initialized,
209  "Error: training parameters must first be initialized.");
210 
211  if (_training_parameters.empty())
212  return 0;
213 
215 }
numeric_index_type _first_local_index
The first sample-vector index from the global vector which is stored in the _training_parameters on t...
numeric_index_type _n_local_training_samples
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ get_local_n_training_samples()

template<class Base >
numeric_index_type libMesh::RBConstructionBase< Base >::get_local_n_training_samples ( ) const

Get the total number of training samples local to this processor.

Definition at line 168 of file rb_construction_base.C.

169 {
170  libmesh_error_msg_if(!_training_parameters_initialized,
171  "Error: training parameters must first be initialized.");
172 
173  // First we check if there are no parameters here, and in that case we
174  // return 1 for both serial and parallel training sets. This is consistent
175  // with get_n_training_samples(), and avoids accessing
176  // training_parameters.begin() when training_parameters is empty.
177  if (_training_parameters.empty())
178  return 1;
179 
181 }
numeric_index_type _n_local_training_samples
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ get_n_continuous_params()

unsigned int libMesh::RBParametrized::get_n_continuous_params ( ) const
inherited

Get the number of continuous parameters.

Definition at line 112 of file rb_parametrized.C.

References libMesh::RBParametrized::get_n_discrete_params(), libMesh::RBParametrized::get_n_params(), libMesh::libmesh_assert(), and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), and libMesh::RBParametrized::write_parameter_ranges_to_file().

113 {
114  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_n_continuous_params");
115 
117 
118  return static_cast<unsigned int>(get_n_params() - get_n_discrete_params());
119 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
unsigned int get_n_discrete_params() const
Get the number of discrete parameters.
libmesh_assert(ctx)
unsigned int get_n_params() const
Get the number of parameters.

◆ get_n_discrete_params()

unsigned int libMesh::RBParametrized::get_n_discrete_params ( ) const
inherited

Get the number of discrete parameters.

Definition at line 121 of file rb_parametrized.C.

References libMesh::RBParametrized::get_discrete_parameter_values(), and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), libMesh::RBParametrized::get_n_continuous_params(), and libMesh::RBParametrized::write_discrete_parameter_values_to_file().

122 {
123  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_n_discrete_params");
124 
125  return cast_int<unsigned int>
127 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values() const
Get a const reference to the discrete parameter values.

◆ get_n_params()

unsigned int libMesh::RBParametrized::get_n_params ( ) const
inherited

Get the number of parameters.

Definition at line 103 of file rb_parametrized.C.

References libMesh::RBParameters::n_parameters(), libMesh::RBParametrized::parameters_initialized, libMesh::RBParametrized::parameters_max, and libMesh::RBParametrized::parameters_min.

Referenced by libMesh::RBParametrized::check_if_valid_params(), libMesh::RBEIMConstruction::compute_max_eim_error(), libMesh::RBConstruction::compute_max_error_bound(), libMesh::RBParametrized::get_n_continuous_params(), libMesh::RBSCMConstruction::print_info(), libMesh::RBEIMConstruction::print_info(), libMesh::RBConstruction::print_info(), libMesh::RBEIMEvaluation::set_eim_error_indicator_active(), and libMesh::RBConstruction::train_reduced_basis_with_POD().

104 {
105  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_n_params");
106 
107  libmesh_assert_equal_to ( parameters_min.n_parameters(), parameters_max.n_parameters() );
108 
109  return parameters_min.n_parameters();
110 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters_min
Vectors that define the ranges (min and max) for the parameters.
unsigned int n_parameters() const
Get the number of parameters that have been added.

◆ get_n_training_samples()

template<class Base >
numeric_index_type libMesh::RBConstructionBase< Base >::get_n_training_samples ( ) const

Get the number of global training samples.

Definition at line 146 of file rb_construction_base.C.

147 {
148  libmesh_error_msg_if(!_training_parameters_initialized,
149  "Error: training parameters must first be initialized.");
150 
151  // First we check if there are no parameters here, and in that case we
152  // return 1 since a single training sample is sufficient to generate an
153  // RB approximation if there are no parameters. Note that in parallel,
154  // and when we don't have a serial training set, set return comm().size()
155  // so that each processor is assigned a single (empty) training sample.
156  if (_training_parameters.empty())
157  {
159  return 1;
160  else
161  return this->comm().size();
162  }
163 
165 }
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
numeric_index_type _n_global_training_samples
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ get_parameter_max()

Real libMesh::RBParametrized::get_parameter_max ( const std::string &  param_name) const
inherited

Get maximum allowable value of parameter param_name.

Definition at line 183 of file rb_parametrized.C.

References libMesh::RBParameters::get_value(), libMesh::RBParametrized::parameters_initialized, and libMesh::RBParametrized::parameters_max.

Referenced by libMesh::RBParametrized::check_if_valid_params(), libMesh::RBSCMConstruction::print_info(), libMesh::RBEIMConstruction::print_info(), and libMesh::RBConstruction::print_info().

184 {
185  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameter_max");
186 
187  return parameters_max.get_value(param_name);
188 }
Real get_value(const std::string &param_name) const
Get the value of the specified parameter, throw an error if it does not exist.
Definition: rb_parameters.C:64
bool parameters_initialized
Flag indicating whether the parameters have been initialized.

◆ get_parameter_min()

Real libMesh::RBParametrized::get_parameter_min ( const std::string &  param_name) const
inherited

Get minimum allowable value of parameter param_name.

Definition at line 176 of file rb_parametrized.C.

References libMesh::RBParameters::get_value(), libMesh::RBParametrized::parameters_initialized, and libMesh::RBParametrized::parameters_min.

Referenced by libMesh::RBParametrized::check_if_valid_params(), libMesh::RBSCMConstruction::print_info(), libMesh::RBEIMConstruction::print_info(), and libMesh::RBConstruction::print_info().

177 {
178  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameter_min");
179 
180  return parameters_min.get_value(param_name);
181 }
Real get_value(const std::string &param_name) const
Get the value of the specified parameter, throw an error if it does not exist.
Definition: rb_parameters.C:64
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters_min
Vectors that define the ranges (min and max) for the parameters.

◆ get_parameter_names()

std::set< std::string > libMesh::RBParametrized::get_parameter_names ( ) const
inherited

Get a set that stores the parameter names.

Definition at line 129 of file rb_parametrized.C.

References libMesh::RBParametrized::parameters_initialized, and libMesh::RBParametrized::parameters_min.

130 {
131  libmesh_deprecated();
132  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameter_names");
133 
134  std::set<std::string> parameter_names;
135  for (const auto & pr : parameters_min)
136  parameter_names.insert(pr.first);
137 
138  return parameter_names;
139 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters_min
Vectors that define the ranges (min and max) for the parameters.

◆ get_parameters()

const RBParameters & libMesh::RBParametrized::get_parameters ( ) const
inherited

Get the current parameters.

Definition at line 155 of file rb_parametrized.C.

References libMesh::RBParametrized::parameters, and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::TransientRBConstruction::add_scaled_mass_matrix(), libMesh::TransientRBEvaluation::cache_online_residual_terms(), libMesh::RBEvaluation::compute_residual_dual_norm(), libMesh::RBSCMConstruction::compute_SCM_bounds_on_training_set(), libMesh::RBSCMConstruction::enrich_C_J(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_nodes(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_sides(), libMesh::RBEvaluation::eval_output_dual_norm(), libMesh::RBSCMConstruction::evaluate_stability_constant(), libMesh::RBConstruction::get_RB_error_bound(), libMesh::RBSCMEvaluation::get_SCM_LB(), libMesh::RBSCMEvaluation::get_SCM_UB(), SimpleRBEvaluation::get_stability_lower_bound(), libMesh::RBConstruction::greedy_termination_test(), libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), libMesh::RBSCMEvaluation::legacy_read_offline_data_from_files(), libMesh::TransientRBConstruction::mass_matrix_scaled_matvec(), libMesh::RBConstruction::preevaluate_thetas(), libMesh::RBSCMConstruction::print_info(), libMesh::RBEIMConstruction::print_info(), libMesh::RBConstruction::print_info(), libMesh::RBParametrized::print_parameters(), libMesh::RBSCMConstruction::process_parameters_file(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::RBSCMEvaluation::save_current_parameters(), libMesh::RBEIMConstruction::train_eim_approximation_with_greedy(), libMesh::RBEIMConstruction::train_eim_approximation_with_POD(), libMesh::TransientRBConstruction::truth_assembly(), libMesh::RBConstruction::truth_assembly(), libMesh::TransientRBConstruction::truth_solve(), libMesh::RBConstruction::truth_solve(), libMesh::TransientRBEvaluation::uncached_compute_residual_dual_norm(), and libMesh::RBConstruction::update_greedy_param_list().

156 {
157  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameters");
158 
159  return parameters;
160 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters
Vector storing the current parameters.

◆ get_parameters_max()

const RBParameters & libMesh::RBParametrized::get_parameters_max ( ) const
inherited

Get an RBParameters object that specifies the maximum allowable value for each parameter.

Definition at line 169 of file rb_parametrized.C.

References libMesh::RBParametrized::parameters_initialized, and libMesh::RBParametrized::parameters_max.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), libMesh::RBParametrized::initialize_parameters(), libMesh::RBSCMConstruction::process_parameters_file(), libMesh::RBEIMConstruction::set_rb_construction_parameters(), libMesh::RBConstruction::set_rb_construction_parameters(), and libMesh::RBParametrized::write_parameter_ranges_to_file().

170 {
171  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameters_max");
172 
173  return parameters_max;
174 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.

◆ get_parameters_min()

const RBParameters & libMesh::RBParametrized::get_parameters_min ( ) const
inherited

Get an RBParameters object that specifies the minimum allowable value for each parameter.

Definition at line 162 of file rb_parametrized.C.

References libMesh::RBParametrized::parameters_initialized, and libMesh::RBParametrized::parameters_min.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), libMesh::RBParametrized::initialize_parameters(), libMesh::RBSCMConstruction::process_parameters_file(), libMesh::RBEIMConstruction::set_rb_construction_parameters(), libMesh::RBConstruction::set_rb_construction_parameters(), and libMesh::RBParametrized::write_parameter_ranges_to_file().

163 {
164  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::get_parameters_min");
165 
166  return parameters_min;
167 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters_min
Vectors that define the ranges (min and max) for the parameters.

◆ get_params_from_training_set()

template<class Base >
RBParameters libMesh::RBConstructionBase< Base >::get_params_from_training_set ( unsigned int  global_index)
protected

Return the RBParameters in index global_index of the global training set.

Why do we use an index here? RBParameters supports loading the full sample set. This seems probably unnecessary now to load individually. Maybe it's a memory issue?

Definition at line 224 of file rb_construction_base.C.

225 {
226  libmesh_error_msg_if(!_training_parameters_initialized,
227  "Error: training parameters must first be initialized.");
228 
229  // If the _training_parameters are empty, return an empty RBParameters.
230  // Otherwise, create a new RBParameters object from the single sample requested.
231  RBParameters params;
232  if (!_training_parameters.empty())
233  {
234  libmesh_error_msg_if((global_index < this->get_first_local_training_index()) ||
235  (global_index >= this->get_last_local_training_index()),
236  "Error: index "
237  << global_index
238  << " must be within range: "
240  << " - "
241  << this->get_last_local_training_index());
242 
243  const numeric_index_type local_index = global_index - get_first_local_training_index();
244  for (const auto & [param_name, sample_vector] : _training_parameters)
245  params.set_value(param_name, sample_vector[local_index]);
246 
247  // Copy all extra values into the new RBParameters.
248  // We assume that the samples may be indexed differently for extra parameters,
249  // so we don't just copy the local_index value.
250  const auto & mine = get_parameters();
251  for (const auto & [key, extra_sample_vector] :
252  as_range(mine.extra_begin(), mine.extra_end()))
253  {
254  for (const auto idx : index_range(extra_sample_vector))
255  params.set_extra_value(key, idx, extra_sample_vector[idx]);
256  }
257  }
258 
259  return params;
260 }
numeric_index_type get_first_local_training_index() const
Get the first local index of the training parameters.
dof_id_type numeric_index_type
Definition: id_types.h:99
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
Helper function that allows us to treat a homogenous pair as a range.
Definition: simple_range.h:57
numeric_index_type get_last_local_training_index() const
Get the last local index of the training parameters.
const RBParameters & get_parameters() const
Get the current parameters.
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)
A useful inline function which replaces the macros used previously.

◆ increment_constructor_count()

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the construction counter.

Should be called in the constructor of any derived class that will be reference counted.

Definition at line 183 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

184 {
185  libmesh_try
186  {
187  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
188  std::pair<unsigned int, unsigned int> & p = _counts[name];
189  p.first++;
190  }
191  libmesh_catch (...)
192  {
193  auto stream = libMesh::err.get();
194  stream->exceptions(stream->goodbit); // stream must not throw
195  libMesh::err << "Encountered unrecoverable error while calling "
196  << "ReferenceCounter::increment_constructor_count() "
197  << "for a(n) " << name << " object." << std::endl;
198  std::terminate();
199  }
200 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ increment_destructor_count()

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the destruction counter.

Should be called in the destructor of any derived class that will be reference counted.

Definition at line 207 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

208 {
209  libmesh_try
210  {
211  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
212  std::pair<unsigned int, unsigned int> & p = _counts[name];
213  p.second++;
214  }
215  libmesh_catch (...)
216  {
217  auto stream = libMesh::err.get();
218  stream->exceptions(stream->goodbit); // stream must not throw
219  libMesh::err << "Encountered unrecoverable error while calling "
220  << "ReferenceCounter::increment_destructor_count() "
221  << "for a(n) " << name << " object." << std::endl;
222  std::terminate();
223  }
224 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ init_data()

template<class Base >
void libMesh::RBConstructionBase< Base >::init_data ( )
protectedvirtual

Initializes the member data fields associated with the system, so that, e.g., assemble() may be used.

Reimplemented in SimpleRBConstruction, SimpleRBConstruction, SimpleRBConstruction, SimpleRBConstruction, SimpleRBConstruction, SimpleRBConstruction, ElasticityRBConstruction, SimpleEIMConstruction, and SimpleEIMConstruction.

Definition at line 122 of file rb_construction_base.C.

123 {
124  Base::init_data();
125 
126  // Initialize the inner product storage vector, which is useful for
127  // storing intermediate results when evaluating inner products
129  inner_product_storage_vector->init (this->n_dofs(), this->n_local_dofs(), false, PARALLEL);
130 }
std::unique_ptr< NumericVector< Number > > inner_product_storage_vector
We keep an extra temporary vector that is useful for performing inner products (avoids unnecessary me...
static std::unique_ptr< NumericVector< T > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...

◆ initialize_parameters() [1/2]

void libMesh::RBParametrized::initialize_parameters ( const RBParameters mu_min_in,
const RBParameters mu_max_in,
const std::map< std::string, std::vector< Real >> &  discrete_parameter_values 
)
inherited

Initialize the parameter ranges and set current_parameters.

Parameter ranges are inclusive. The input min/max RBParameters should have exactly 1 sample each. Vector-valued samples are not currently supported for the min/max parameters or for discrete parameters.

Definition at line 53 of file rb_parametrized.C.

References libMesh::RBParametrized::_discrete_parameter_values, libMesh::RBParameters::begin_serialized(), libMesh::RBParameters::end_serialized(), libMesh::RBParameters::n_parameters(), libMesh::RBParameters::n_samples(), libMesh::Quality::name(), libMesh::RBParametrized::parameters_initialized, libMesh::RBParametrized::parameters_max, libMesh::RBParametrized::parameters_min, libMesh::Real, libMesh::RBParametrized::set_parameters(), and libMesh::RBParameters::set_value().

Referenced by libMesh::RBConstruction::enrich_basis_from_rhs_terms(), libMesh::RBParametrized::initialize_parameters(), libMesh::RBDataDeserialization::load_parameter_ranges(), libMesh::RBSCMConstruction::perform_SCM_greedy(), libMesh::RBSCMConstruction::process_parameters_file(), libMesh::RBParametrized::read_parameter_data_from_files(), libMesh::RBEIMConstruction::set_rb_construction_parameters(), libMesh::RBConstruction::set_rb_construction_parameters(), RBParametersTest::testRBParametrized(), libMesh::RBEIMConstruction::train_eim_approximation_with_greedy(), libMesh::RBEIMConstruction::train_eim_approximation_with_POD(), libMesh::RBConstruction::train_reduced_basis_with_greedy(), and libMesh::RBConstruction::train_reduced_basis_with_POD().

56 {
57  // Check that the min/max vectors have the same size.
58  libmesh_error_msg_if(mu_min_in.n_parameters() != mu_max_in.n_parameters(),
59  "Error: Invalid mu_min/mu_max in initialize_parameters(), different number of parameters.");
60  libmesh_error_msg_if(mu_min_in.n_samples() != 1 ||
61  mu_max_in.n_samples() != 1,
62  "Error: Invalid mu_min/mu_max in initialize_parameters(), only 1 sample supported.");
63 
64  // Ensure all the values are valid for min and max.
65  auto pr_min = mu_min_in.begin_serialized();
66  auto pr_max = mu_max_in.begin_serialized();
67  for (; pr_min != mu_min_in.end_serialized(); ++pr_min, ++pr_max)
68  libmesh_error_msg_if((*pr_min).second > (*pr_max).second,
69  "Error: Invalid mu_min/mu_max in RBParameters constructor.");
70 
71  parameters_min = mu_min_in;
72  parameters_max = mu_max_in;
73 
74  // Add in min/max values due to the discrete parameters
75  for (const auto & [name, vals] : discrete_parameter_values)
76  {
77  libmesh_error_msg_if(vals.empty(), "Error: List of discrete parameters for " << name << " is empty.");
78 
79  Real min_val = *std::min_element(vals.begin(), vals.end());
80  Real max_val = *std::max_element(vals.begin(), vals.end());
81 
82  libmesh_assert_less_equal(min_val, max_val);
83 
84  parameters_min.set_value(name, min_val);
85  parameters_max.set_value(name, max_val);
86  }
87 
88  _discrete_parameter_values = discrete_parameter_values;
89 
91 
92  // Initialize the current parameters to parameters_min
94 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
RBParameters parameters_min
Vectors that define the ranges (min and max) for the parameters.
std::map< std::string, std::vector< Real > > _discrete_parameter_values
Map that defines the allowable values of any discrete parameters.
bool set_parameters(const RBParameters &params)
Set the current parameters to params The parameters are checked for validity; an error is thrown if t...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void set_value(const std::string &param_name, Real value)
Set the value of the specified parameter.

◆ initialize_parameters() [2/2]

void libMesh::RBParametrized::initialize_parameters ( const RBParametrized rb_parametrized)
inherited

Initialize the parameter ranges and set current_parameters.

Definition at line 96 of file rb_parametrized.C.

References libMesh::RBParametrized::get_discrete_parameter_values(), libMesh::RBParametrized::get_parameters_max(), libMesh::RBParametrized::get_parameters_min(), and libMesh::RBParametrized::initialize_parameters().

97 {
98  initialize_parameters(rb_parametrized.get_parameters_min(),
99  rb_parametrized.get_parameters_max(),
100  rb_parametrized.get_discrete_parameter_values());
101 }
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.

◆ initialize_training_parameters()

template<class Base >
void libMesh::RBConstructionBase< Base >::initialize_training_parameters ( const RBParameters mu_min,
const RBParameters mu_max,
const unsigned int  n_global_training_samples,
const std::map< std::string, bool > &  log_param_scale,
const bool  deterministic = true 
)
virtual

Initialize the parameter ranges and indicate whether deterministic or random training parameters should be used and whether or not we want the parameters to be scaled logarithmically.

n_global_training_samples is the total number of samples to generate, which will be distributed across all the processors.

Definition at line 285 of file rb_construction_base.C.

290 {
291  if (!is_quiet())
292  {
293  // Print out some info about the training set initialization
294  libMesh::out << "Initializing training parameters with "
295  << (deterministic ? "deterministic " : "random " )
296  << "training set..." << std::endl;
297 
298  for (const auto & pr : log_param_scale)
299  libMesh::out << "Parameter "
300  << pr.first
301  << ": log scaling = "
302  << pr.second
303  << std::endl;
304 
305  libMesh::out << std::endl;
306  }
307 
308  if (deterministic)
309  {
310  const auto [first_local_index, last_local_index] =
312  log_param_scale,
314  n_global_training_samples,
315  mu_min,
316  mu_max,
318  _first_local_index = first_local_index;
319  _n_local_training_samples = last_local_index-first_local_index;
320  }
321  else
322  {
323  // Generate random training samples for all parameters
324  const auto [first_local_index, last_local_index] =
326  log_param_scale,
328  n_global_training_samples,
329  mu_min,
330  mu_max,
333  _first_local_index = first_local_index;
334  _n_local_training_samples = last_local_index-first_local_index;
335  }
337 
338  if (!serial_training_set)
339  this->comm().sum(_n_global_training_samples);
340 
341  // For each parameter that only allows discrete values, we "snap" to the nearest
342  // allowable discrete value
343  if (get_n_discrete_params() > 0)
344  {
345  for (auto & [param_name, sample_vector] : _training_parameters)
346  {
347  if (is_discrete_parameter(param_name))
348  {
349  std::vector<Real> discrete_values =
350  get_discrete_parameter_values().find(param_name)->second;
351  for (const auto sample_idx : index_range(sample_vector))
352  {
353  // Round all values to the closest discrete value.
354  std::vector<Real> discretized_vector(sample_vector[sample_idx].size());
355  std::transform(sample_vector[sample_idx].cbegin(),
356  sample_vector[sample_idx].cend(),
357  discretized_vector.begin(),
358  [&discrete_values](const Real & val) {
359  return get_closest_value(val, discrete_values);
360  });
361  sample_vector[sample_idx] = discretized_vector;
362  }
363  }
364  }
365  }
366 
368 }
static std::pair< std::size_t, std::size_t > generate_training_parameters_random(const Parallel::Communicator &communicator, const std::map< std::string, bool > &log_param_scale, std::map< std::string, std::vector< RBParameter >> &local_training_parameters_in, const unsigned int n_global_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, const int training_parameters_random_seed=-1, const bool serial_training_set=false)
Static helper function for generating a randomized set of parameters.
static std::pair< std::size_t, std::size_t > generate_training_parameters_deterministic(const Parallel::Communicator &communicator, const std::map< std::string, bool > &log_param_scale, std::map< std::string, std::vector< RBParameter >> &local_training_parameters_in, const unsigned int n_global_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, const bool serial_training_set=false)
Static helper function for generating a deterministic set of parameters.
bool is_quiet() const
Is the system in quiet mode?
static Real get_closest_value(Real value, const std::vector< Real > &list_of_values)
unsigned int get_n_discrete_params() const
Get the number of discrete parameters.
numeric_index_type _first_local_index
The first sample-vector index from the global vector which is stored in the _training_parameters on t...
numeric_index_type _n_local_training_samples
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values() const
Get a const reference to the discrete parameter values.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OStreamProxy out
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...
int _training_parameters_random_seed
If < 0, use std::time() * processor_id() to seed the random number generator for the training paramet...
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
numeric_index_type _n_global_training_samples
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.
bool is_discrete_parameter(const std::string &mu_name) const
Is parameter mu_name discrete?

◆ is_discrete_parameter()

bool libMesh::RBParametrized::is_discrete_parameter ( const std::string &  mu_name) const
inherited

Is parameter mu_name discrete?

Definition at line 363 of file rb_parametrized.C.

References libMesh::RBParametrized::_discrete_parameter_values, and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::RBDataSerialization::add_parameter_ranges_to_builder(), libMesh::RBEIMConstruction::print_info(), libMesh::RBConstruction::print_info(), and libMesh::RBParametrized::write_parameter_ranges_to_file().

364 {
365  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::is_discrete_parameter");
366 
367  return (_discrete_parameter_values.find(mu_name) != _discrete_parameter_values.end());
368 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
std::map< std::string, std::vector< Real > > _discrete_parameter_values
Map that defines the allowable values of any discrete parameters.

◆ is_quiet()

template<class Base>
bool libMesh::RBConstructionBase< Base >::is_quiet ( ) const
inline

Is the system in quiet mode?

Definition at line 106 of file rb_construction_base.h.

107  { return this->quiet_mode; }
bool quiet_mode
Flag to indicate whether we print out extra information during the Offline stage. ...

◆ load_training_set()

template<class Base >
void libMesh::RBConstructionBase< Base >::load_training_set ( const std::map< std::string, std::vector< RBParameter >> &  new_training_set)
virtual

Overwrite the training parameters with new_training_set.

This training set is assumed to contain only the samples local to this processor.

Definition at line 371 of file rb_construction_base.C.

372 {
373  // Make sure we're running this on all processors at the same time
374  libmesh_parallel_only(this->comm());
375 
376  // First, make sure that an initial training set has already been generated
377  libmesh_error_msg_if(!_training_parameters_initialized,
378  "Error: load_training_set cannot be used to initialize parameters");
379 
380  // Make sure that the training set has the correct number of parameters
381  const unsigned int n_params = get_n_params();
382  libmesh_error_msg_if(new_training_set.size() > n_params,
383  "Error: new_training_set should not have more than get_n_params() parameters.");
384 
385  // Check that (new_training_set.size() == get_n_params()) is the same on all processes so that
386  // we go into the same branch of the "if" statement below on all processes.
387  const bool size_matches = (new_training_set.size() == n_params);
388  this->comm().verify(size_matches);
389 
390  if (size_matches)
391  {
392  // If new_training_set stores values for all parameters, then we overwrite
393  // _training_parameters with new_training_set.
394 
395  // Get the number of local and global training parameters
396  _first_local_index = 0;
398  cast_int<numeric_index_type>(new_training_set.begin()->second.size());
400 
401  if (!serial_training_set)
402  {
403  this->comm().sum(_n_global_training_samples);
404 
405  // Set the first/last indices.
406  std::vector<numeric_index_type> local_sizes (this->n_processors(), 0);
407  local_sizes[this->processor_id()] = _n_local_training_samples;
408  this->comm().sum(local_sizes);
409 
410  // first_local_index is the sum of local_sizes
411  // for all processor ids less than ours
412  for (auto p : make_range(this->processor_id()))
413  _first_local_index += local_sizes[p];
414  }
415 
416  // Ensure that the parameters are the same.
417  for (const auto & pr : _training_parameters)
418  libmesh_error_msg_if(!new_training_set.count(pr.first),
419  "Parameters must be identical in order to overwrite dataset.");
420 
421  // Copy the values from the new_training_set to the internal training_parameters.
422  _training_parameters = new_training_set;
423  }
424  else
425  {
426  // If new_training_set stores values for a subset of the parameters, then we keep the
427  // length of training_parameters unchanged and overwrite the entries of the specified
428  // parameters from new_training_set. Note that we repeatedly loop over new_training_set
429  // to fill up the entire length of the sample_vector.
430  for (auto & [param_name, sample_vector]: _training_parameters)
431  {
432  if (new_training_set.count(param_name))
433  {
434  for (const auto i : make_range(get_local_n_training_samples()))
435  {
436  const unsigned int num_new_samples = libmesh_map_find(new_training_set,param_name).size();
437  libmesh_error_msg_if (num_new_samples==0, "new_training_set set should not be empty");
438 
439  const unsigned int new_training_set_index = i % num_new_samples;
440  sample_vector[i] = libmesh_map_find(new_training_set,param_name)[new_training_set_index];
441  }
442  }
443  }
444  }
445 }
numeric_index_type _first_local_index
The first sample-vector index from the global vector which is stored in the _training_parameters on t...
numeric_index_type _n_local_training_samples
numeric_index_type get_local_n_training_samples() const
Get the total number of training samples local to this processor.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
bool serial_training_set
This boolean flag indicates whether or not the training set should be the same on all processors...
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
unsigned int get_n_params() const
Get the number of parameters.
numeric_index_type _n_global_training_samples
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ n_objects()

static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited

Prints the number of outstanding (created, but not yet destroyed) objects.

Definition at line 85 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

Referenced by libMesh::LibMeshInit::~LibMeshInit().

86  { return _n_objects; }
static Threads::atomic< unsigned int > _n_objects
The number of objects.

◆ operator=() [1/2]

template<class Base>
RBConstructionBase& libMesh::RBConstructionBase< Base >::operator= ( RBConstructionBase< Base > &&  )
delete

◆ operator=() [2/2]

template<class Base>
RBConstructionBase& libMesh::RBConstructionBase< Base >::operator= ( const RBConstructionBase< Base > &  )
delete

◆ print_discrete_parameter_values()

void libMesh::RBParametrized::print_discrete_parameter_values ( ) const
inherited

Print out all the discrete parameter values.

Definition at line 377 of file rb_parametrized.C.

References libMesh::RBParametrized::get_discrete_parameter_values(), libMesh::Quality::name(), libMesh::out, and value.

Referenced by libMesh::RBSCMConstruction::print_info(), libMesh::RBEIMConstruction::print_info(), and libMesh::RBConstruction::print_info().

378 {
379  for (const auto & [name, values] : get_discrete_parameter_values())
380  {
381  libMesh::out << "Discrete parameter " << name << ", values: ";
382 
383  for (const auto & value : values)
384  libMesh::out << value << " ";
385  libMesh::out << std::endl;
386  }
387 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values() const
Get a const reference to the discrete parameter values.
OStreamProxy out
static const bool value
Definition: xdr_io.C:54

◆ print_info()

void libMesh::ReferenceCounter::print_info ( std::ostream &  out_stream = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 81 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

Referenced by libMesh::LibMeshInit::~LibMeshInit().

82 {
84  out_stream << ReferenceCounter::get_info();
85 }
static std::string get_info()
Gets a string containing the reference information.
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ print_parameters()

void libMesh::RBParametrized::print_parameters ( ) const
inherited

Print the current parameters.

Definition at line 190 of file rb_parametrized.C.

References libMesh::RBParametrized::get_parameters(), libMesh::RBParametrized::parameters_initialized, and libMesh::RBParameters::print().

Referenced by libMesh::RBEIMConstruction::train_eim_approximation_with_greedy(), and libMesh::RBConstruction::train_reduced_basis_with_greedy().

191 {
192  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::print_current_parameters");
193 
194  get_parameters().print();
195 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
void print() const
Print the parameters.
const RBParameters & get_parameters() const
Get the current parameters.

◆ read_parameter_data_from_files()

void libMesh::RBParametrized::read_parameter_data_from_files ( const std::string &  continuous_param_file_name,
const std::string &  discrete_param_file_name,
const bool  read_binary_data 
)
inherited

Read in the parameter ranges from files.

Definition at line 274 of file rb_parametrized.C.

References libMesh::RBParametrized::initialize_parameters(), libMesh::RBParametrized::read_discrete_parameter_values_from_file(), and libMesh::RBParametrized::read_parameter_ranges_from_file().

Referenced by libMesh::RBSCMEvaluation::legacy_read_offline_data_from_files(), and libMesh::RBEvaluation::legacy_read_offline_data_from_files().

277 {
278  RBParameters param_min;
279  RBParameters param_max;
280  read_parameter_ranges_from_file(continuous_param_file_name,
281  read_binary_data,
282  param_min,
283  param_max);
284 
285  std::map<std::string, std::vector<Real>> discrete_parameter_values_in;
286  read_discrete_parameter_values_from_file(discrete_param_file_name,
287  read_binary_data,
288  discrete_parameter_values_in);
289 
290  initialize_parameters(param_min, param_max, discrete_parameter_values_in);
291 }
void read_parameter_ranges_from_file(const std::string &file_name, const bool read_binary, RBParameters &param_min, RBParameters &param_max)
Read in the parameter ranges from file.
void read_discrete_parameter_values_from_file(const std::string &file_name, const bool read_binary_data, std::map< std::string, std::vector< Real >> &discrete_parameter_values_in)
Read in the discrete parameter values from file, if we have any.
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.

◆ set_deterministic_training_parameter_name()

template<class Base>
void libMesh::RBConstructionBase< Base >::set_deterministic_training_parameter_name ( const std::string &  name)

In some cases we only want to allow discrete parameter values, instead of parameters that may take any value in a specified interval.

Here we provide a method to set the d Set the discrete values for parameter mu that are allowed in the training set. This must be called before the training set is generated. Set the name of the parameter that we will generate deterministic training parameters for. Defaults to "NONE".

◆ set_parameters()

bool libMesh::RBParametrized::set_parameters ( const RBParameters params)
inherited

Set the current parameters to params The parameters are checked for validity; an error is thrown if the number of parameters or samples is different than expected.

We

Returns
a boolean true if the new parameters are within the min/max range, and false otherwise (but the parameters are set regardless). Enabling the "verbose_mode" flag will also print more details.

Definition at line 141 of file rb_parametrized.C.

References libMesh::RBParametrized::check_if_valid_params(), libMesh::RBParametrized::parameters, and libMesh::RBParametrized::parameters_initialized.

Referenced by libMesh::RBSCMConstruction::compute_SCM_bounds_on_training_set(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_nodes(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_sides(), libMesh::RBConstruction::get_RB_error_bound(), SimpleRBEvaluation::get_stability_lower_bound(), libMesh::RBParametrized::initialize_parameters(), libMesh::RBSCMEvaluation::reload_current_parameters(), libMesh::RBSCMEvaluation::set_current_parameters_from_C_J(), and RBParametersTest::testRBParametrized().

142 {
143  libmesh_error_msg_if(!parameters_initialized, "Error: parameters not initialized in RBParametrized::set_parameters");
144 
145  // Terminate if params has the wrong number of parameters or samples.
146  // If the parameters are outside the min/max range, return false.
147  const bool valid_params = check_if_valid_params(params);
148 
149  // Make a copy of params (default assignment operator just does memberwise copy, which is sufficient here)
150  this->parameters = params;
151 
152  return valid_params;
153 }
bool parameters_initialized
Flag indicating whether the parameters have been initialized.
bool check_if_valid_params(const RBParameters &params) const
Helper function to check that params is valid:
RBParameters parameters
Vector storing the current parameters.

◆ set_params_from_training_set()

template<class Base >
void libMesh::RBConstructionBase< Base >::set_params_from_training_set ( unsigned int  global_index)
protected

Set parameters to the RBParameters stored in index global_index of the global training set.

Definition at line 218 of file rb_construction_base.C.

219 {
221 }
bool set_parameters(const RBParameters &params)
Set the current parameters to params The parameters are checked for validity; an error is thrown if t...
RBParameters get_params_from_training_set(unsigned int global_index)
Return the RBParameters in index global_index of the global training set.

◆ set_params_from_training_set_and_broadcast()

template<class Base >
void libMesh::RBConstructionBase< Base >::set_params_from_training_set_and_broadcast ( unsigned int  global_index)
protectedvirtual

Load the specified training parameter and then broadcast to all processors.

Definition at line 263 of file rb_construction_base.C.

264 {
265  libmesh_error_msg_if(!_training_parameters_initialized,
266  "Error: training parameters must first be initialized.");
267 
268  processor_id_type root_id = 0;
269  if ((this->get_first_local_training_index() <= global_index) &&
270  (global_index < this->get_last_local_training_index()))
271  {
272  // Set parameters on only one processor
273  set_params_from_training_set(global_index);
274 
275  // set root_id, only non-zero on one processor
276  root_id = this->processor_id();
277  }
278 
279  // broadcast
280  this->comm().max(root_id);
281  broadcast_parameters(root_id);
282 }
numeric_index_type get_first_local_training_index() const
Get the first local index of the training parameters.
void broadcast_parameters(const unsigned int proc_id)
Broadcasts parameters from processor proc_id to all processors.
uint8_t processor_id_type
numeric_index_type get_last_local_training_index() const
Get the last local index of the training parameters.
void set_params_from_training_set(unsigned int global_index)
Set parameters to the RBParameters stored in index global_index of the global training set...
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ set_quiet_mode()

template<class Base>
void libMesh::RBConstructionBase< Base >::set_quiet_mode ( bool  quiet_mode_in)
inline

Set the quiet_mode flag.

If quiet == false then we print out a lot of extra information during the Offline stage.

Definition at line 100 of file rb_construction_base.h.

101  { this->quiet_mode = quiet_mode_in; }
bool quiet_mode
Flag to indicate whether we print out extra information during the Offline stage. ...

◆ set_training_parameter_values()

template<class Base >
void libMesh::RBConstructionBase< Base >::set_training_parameter_values ( const std::string &  param_name,
const std::vector< RBParameter > &  values 
)

Overwrite the local training samples for param_name using values.

This assumes that values.size() matches get_local_n_training_samples().

Definition at line 448 of file rb_construction_base.C.

450 {
451  libmesh_error_msg_if(!_training_parameters_initialized,
452  "Training parameters must be initialized before calling set_training_parameter_values");
453  libmesh_error_msg_if(values.size() != get_local_n_training_samples(),
454  "Inconsistent sizes");
455 
456  // Copy the new data, overwriting the old data.
457  auto & training_vector = libmesh_map_find(_training_parameters, param_name);
458  training_vector = values;
459 }
numeric_index_type get_local_n_training_samples() const
Get the total number of training samples local to this processor.
std::map< std::string, std::vector< RBParameter > > _training_parameters
The training samples for each parameter.
bool _training_parameters_initialized
Boolean flag to indicate whether or not the parameter ranges have been initialized.

◆ set_training_random_seed()

template<class Base >
void libMesh::RBConstructionBase< Base >::set_training_random_seed ( int  seed)

Set the seed that is used to randomly generate training parameters.

Definition at line 772 of file rb_construction_base.C.

773 {
775 }
int _training_parameters_random_seed
If < 0, use std::time() * processor_id() to seed the random number generator for the training paramet...

◆ system()

template<class Base>
sys_type& libMesh::RBConstructionBase< Base >::system ( )
inline
Returns
A reference to *this.

Definition at line 87 of file rb_construction_base.h.

87 { return *this; }

◆ write_parameter_data_to_files()

void libMesh::RBParametrized::write_parameter_data_to_files ( const std::string &  continuous_param_file_name,
const std::string &  discrete_param_file_name,
const bool  write_binary_data 
)
inherited

Write out the parameter ranges to files.

Definition at line 197 of file rb_parametrized.C.

References libMesh::RBParametrized::write_discrete_parameter_values_to_file(), and libMesh::RBParametrized::write_parameter_ranges_to_file().

Referenced by libMesh::RBSCMEvaluation::legacy_write_offline_data_to_files(), and libMesh::RBEvaluation::legacy_write_offline_data_to_files().

200 {
201  write_parameter_ranges_to_file(continuous_param_file_name, write_binary_data);
202  write_discrete_parameter_values_to_file(discrete_param_file_name, write_binary_data);
203 }
void write_discrete_parameter_values_to_file(const std::string &file_name, const bool write_binary_data)
Write out the discrete parameter values to file.
void write_parameter_ranges_to_file(const std::string &file_name, const bool write_binary)
Write out the parameter ranges to file.

Member Data Documentation

◆ _counts

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited

Actually holds the data.

Definition at line 124 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::get_info().

◆ _enable_print_counter

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 143 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

◆ _first_local_index

template<class Base>
numeric_index_type libMesh::RBConstructionBase< Base >::_first_local_index
private

The first sample-vector index from the global vector which is stored in the _training_parameters on this processor.

_n_local_training_samples is equivalent to the .size() of any vector in _training_parameters.

Definition at line 299 of file rb_construction_base.h.

◆ _mutex

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 137 of file reference_counter.h.

◆ _n_global_training_samples

template<class Base>
numeric_index_type libMesh::RBConstructionBase< Base >::_n_global_training_samples
private

Definition at line 301 of file rb_construction_base.h.

◆ _n_local_training_samples

template<class Base>
numeric_index_type libMesh::RBConstructionBase< Base >::_n_local_training_samples
private

Definition at line 300 of file rb_construction_base.h.

◆ _n_objects

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects.

Print the reference count information when the number returns to 0.

Definition at line 132 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

◆ _training_parameters

template<class Base>
std::map<std::string, std::vector<RBParameter> > libMesh::RBConstructionBase< Base >::_training_parameters
private

The training samples for each parameter.

When serial_training_set is true, the map contains all samples of all parameters. Otherwise, the sample vectors will only contain the values for the local samples as defined by _first_local_index and _n_local_training_samples. Mapped from parameter_name -> sample_vector -> value_vector.

Definition at line 292 of file rb_construction_base.h.

◆ _training_parameters_initialized

template<class Base>
bool libMesh::RBConstructionBase< Base >::_training_parameters_initialized
private

Boolean flag to indicate whether or not the parameter ranges have been initialized.

Definition at line 283 of file rb_construction_base.h.

◆ _training_parameters_random_seed

template<class Base>
int libMesh::RBConstructionBase< Base >::_training_parameters_random_seed
private

If < 0, use std::time() * processor_id() to seed the random number generator for the training parameters (default).

If >= 0, use the provided value * processor_id() as the random number generator seed.

Definition at line 309 of file rb_construction_base.h.

◆ inner_product_storage_vector

template<class Base>
std::unique_ptr<NumericVector<Number> > libMesh::RBConstructionBase< Base >::inner_product_storage_vector
protected

We keep an extra temporary vector that is useful for performing inner products (avoids unnecessary memory allocation/deallocation).

Definition at line 274 of file rb_construction_base.h.

◆ quiet_mode

template<class Base>
bool libMesh::RBConstructionBase< Base >::quiet_mode
protected

Flag to indicate whether we print out extra information during the Offline stage.

Definition at line 259 of file rb_construction_base.h.

Referenced by libMesh::RBConstructionBase< CondensedEigenSystem >::is_quiet(), and libMesh::RBConstructionBase< CondensedEigenSystem >::set_quiet_mode().

◆ serial_training_set

template<class Base>
bool libMesh::RBConstructionBase< Base >::serial_training_set
protected

This boolean flag indicates whether or not the training set should be the same on all processors.

By default it is false, but in the case of the Empirical Interpolation Method (RBEIMConstruction), for example, we need the training set to be identical on all processors.

Definition at line 267 of file rb_construction_base.h.

◆ verbose_mode

bool libMesh::RBParametrized::verbose_mode
inherited

Public boolean to toggle verbose mode.

Definition at line 191 of file rb_parametrized.h.

Referenced by libMesh::RBParametrized::check_if_valid_params().


The documentation for this class was generated from the following files: