libMesh::RBEIMConstruction Class Reference

#include <rb_eim_construction.h>

Inheritance diagram for libMesh::RBEIMConstruction:
Inheritance graph
[legend]

List of all members.

Public Types

enum  BEST_FIT_TYPE { PROJECTION_BEST_FIT, EIM_BEST_FIT }
typedef RBEIMConstruction sys_type
typedef RBConstruction Parent
typedef std::map< std::string,
SparseMatrix< Number >
* >::iterator 		matrices_iterator
typedef std::map< std::string,
SparseMatrix< Number >
* >::const_iterator 		const_matrices_iterator
typedef std::map< std::string,
NumericVector< Number >
* >::iterator 		vectors_iterator
typedef std::map< std::string,
NumericVector< Number >
* >::const_iterator 		const_vectors_iterator

Public Member Functions

 RBEIMConstruction (EquationSystems &es, const std::string &name, const unsigned int number)
virtual ~RBEIMConstruction ()
virtual void clear ()
virtual void process_parameters_file (const std::string &parameters_filename)
virtual void print_info ()
virtual void initialize_rb_construction ()
virtual Real train_reduced_basis (const std::string &directory_name="offline_data", const bool resize_rb_eval_data=true)
virtual Real truth_solve (int plot_solution)
virtual Real compute_best_fit_error ()
virtual void init_context (FEMContext &c)
Number evaluate_mesh_function (unsigned int var_number, Point p)
virtual void initialize_eim_assembly_objects ()
std::vector< ElemAssembly * > get_eim_assembly_objects ()
virtual AutoPtr< ElemAssembly > build_eim_assembly (unsigned int bf_index)=0
void set_rb_evaluation (RBEvaluation &rb_eval_in)
RBEvaluation & get_rb_evaluation ()
bool is_rb_eval_initialized () const
RBThetaExpansion & get_rb_theta_expansion ()
void set_rb_assembly_expansion (RBAssemblyExpansion &rb_assembly_expansion_in)
RBAssemblyExpansion & get_rb_assembly_expansion ()
sys_type & system ()
virtual std::string system_type () const
virtual Real compute_max_error_bound ()
const RBParameters & get_greedy_parameter (unsigned int i)
void set_training_tolerance (Real new_training_tolerance)
Real get_training_tolerance ()
unsigned int get_Nmax () const
virtual void set_Nmax (unsigned int Nmax)
void set_quiet_mode (bool quiet_mode_in)
bool is_quiet () const
virtual void load_basis_function (unsigned int i)
virtual void load_rb_solution ()
SparseMatrix< Number > * get_inner_product_matrix ()
SparseMatrix< Number > * get_non_dirichlet_inner_product_matrix ()
SparseMatrix< Number > * get_Aq (unsigned int q)
SparseMatrix< Number > * get_non_dirichlet_Aq (unsigned int q)
NumericVector< Number > * get_Fq (unsigned int q)
NumericVector< Number > * get_non_dirichlet_Fq (unsigned int q)
NumericVector< Number > * get_output_vector (unsigned int n, unsigned int q_l)
NumericVector< Number > * get_non_dirichlet_output_vector (unsigned int n, unsigned int q_l)
void assemble_inner_product_matrix (SparseMatrix< Number > *input_matrix, bool apply_dof_constraints=true)
void assemble_constraint_matrix (SparseMatrix< Number > *input_matrix)
void assemble_and_add_constraint_matrix (SparseMatrix< Number > *input_matrix)
void assemble_Aq_matrix (unsigned int q, SparseMatrix< Number > *input_matrix, bool apply_dof_constraints=true)
void assemble_Fq_vector (unsigned int q, NumericVector< Number > *input_vector, bool apply_dof_constraints=true)
void add_scaled_Aq (Number scalar, unsigned int q_a, SparseMatrix< Number > *input_matrix, bool symmetrize)
virtual void write_riesz_representors_to_files (const std::string &riesz_representors_dir, const bool write_binary_residual_representors)
virtual void read_riesz_representors_from_files (const std::string &riesz_representors_dir, const bool write_binary_residual_representors)
virtual void recompute_all_residual_terms (const bool compute_inner_products=true)
unsigned int get_delta_N () const
void set_inner_product_assembly (ElemAssembly &inner_product_assembly_in)
ElemAssembly & get_inner_product_assembly ()
void set_constraint_assembly (ElemAssembly &constraint_assembly_in)
ElemAssembly & get_constraint_assembly ()
void zero_constrained_dofs_on_vector (NumericVector< Number > &vector)
numeric_index_type get_n_training_samples () const
numeric_index_type get_local_n_training_samples () const
numeric_index_type get_first_local_training_index () const
numeric_index_type get_last_local_training_index () const
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)
virtual void load_training_set (std::map< std::string, std::vector< Number > > &new_training_set)
std::pair< std::string,
std::string > 		set_alternative_solver (AutoPtr< LinearSolver< Number > > &ls)
void reset_alternative_solver (AutoPtr< LinearSolver< Number > > &ls, const std::pair< std::string, std::string > &orig)
void broadcast_parameters (unsigned int proc_id)
void set_training_random_seed (unsigned int seed)
void set_deterministic_training_parameter_name (const std::string name)
const std::string & get_deterministic_training_parameter_name () const
void set_deterministic_training_parameter_repeats (unsigned int repeats)
unsigned int get_deterministic_training_parameter_repeats () const
virtual void reinit ()
virtual void assemble ()
virtual void restrict_solve_to (const SystemSubset *subset, const SubsetSolveMode subset_solve_mode=SUBSET_ZERO)
virtual void solve ()
virtual LinearSolver< Number > * get_linear_solver () const
virtual void release_linear_solver (LinearSolver< Number > *) const
virtual void assembly (bool get_residual, bool get_jacobian)
unsigned int n_linear_iterations () const
Real final_linear_residual () const
void attach_shell_matrix (ShellMatrix< Number > *shell_matrix)
void detach_shell_matrix (void)
ShellMatrix< Number > * get_shell_matrix (void)
virtual void disable_cache ()
virtual std::pair< unsigned
int, Real > 		get_linear_solve_parameters () const
virtual void assemble_residual_derivatives (const ParameterVector &parameters)
virtual std::pair< unsigned
int, Real > 		sensitivity_solve (const ParameterVector &parameters)
virtual std::pair< unsigned
int, Real > 		weighted_sensitivity_solve (const ParameterVector &parameters, const ParameterVector &weights)
virtual std::pair< unsigned
int, Real > 		adjoint_solve (const QoISet &qoi_indices=QoISet())
virtual std::pair< unsigned
int, Real > 		weighted_sensitivity_adjoint_solve (const ParameterVector &parameters, const ParameterVector &weights, const QoISet &qoi_indices=QoISet())
virtual void adjoint_qoi_parameter_sensitivity (const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities)
virtual void forward_qoi_parameter_sensitivity (const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities)
virtual void qoi_parameter_hessian (const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &hessian)
virtual void qoi_parameter_hessian_vector_product (const QoISet &qoi_indices, const ParameterVector &parameters, const ParameterVector &vector, SensitivityData &product)
SparseMatrix< Number > & add_matrix (const std::string &mat_name)
bool have_matrix (const std::string &mat_name) const
const SparseMatrix< Number > * request_matrix (const std::string &mat_name) const
SparseMatrix< Number > * request_matrix (const std::string &mat_name)
const SparseMatrix< Number > & get_matrix (const std::string &mat_name) const
SparseMatrix< Number > & get_matrix (const std::string &mat_name)
virtual unsigned int n_matrices () const
virtual void assemble_qoi (const QoISet &qoi_indices=QoISet())
virtual void assemble_qoi_derivative (const QoISet &qoi_indices=QoISet())
void init ()
virtual void update ()
bool is_adjoint_already_solved () const
void set_adjoint_already_solved (bool setting)
virtual void qoi_parameter_sensitivity (const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities)
virtual bool compare (const System &other_system, const Real threshold, const bool verbose) const
const std::string & name () const
void project_solution (FunctionBase< Number > *f, FunctionBase< Gradient > *g=NULL) const
void project_solution (FEMFunctionBase< Number > *f, FEMFunctionBase< Gradient > *g=NULL) const
void project_solution (Number fptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), Gradient gptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), const Parameters &parameters) const
void project_vector (NumericVector< Number > &new_vector, FunctionBase< Number > *f, FunctionBase< Gradient > *g=NULL) const
void project_vector (NumericVector< Number > &new_vector, FEMFunctionBase< Number > *f, FEMFunctionBase< Gradient > *g=NULL) const
void project_vector (Number fptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), Gradient gptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), const Parameters &parameters, NumericVector< Number > &new_vector) const
void boundary_project_solution (const std::set< boundary_id_type > &b, const std::vector< unsigned int > &variables, FunctionBase< Number > *f, FunctionBase< Gradient > *g=NULL)
void boundary_project_solution (const std::set< boundary_id_type > &b, const std::vector< unsigned int > &variables, Number fptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), Gradient gptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), const Parameters &parameters)
void boundary_project_vector (const std::set< boundary_id_type > &b, const std::vector< unsigned int > &variables, NumericVector< Number > &new_vector, FunctionBase< Number > *f, FunctionBase< Gradient > *g=NULL) const
void boundary_project_vector (const std::set< boundary_id_type > &b, const std::vector< unsigned int > &variables, Number fptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), Gradient gptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name), const Parameters &parameters, NumericVector< Number > &new_vector) const
unsigned int number () const
void update_global_solution (std::vector< Number > &global_soln) const
void update_global_solution (std::vector< Number > &global_soln, const unsigned int dest_proc) const
const MeshBase & get_mesh () const
MeshBase & get_mesh ()
const DofMap & get_dof_map () const
DofMap & get_dof_map ()
const EquationSystems & get_equation_systems () const
EquationSystems & get_equation_systems ()
bool active () const
void activate ()
void deactivate ()
void set_basic_system_only ()
vectors_iterator vectors_begin ()
const_vectors_iterator vectors_begin () const
vectors_iterator vectors_end ()
const_vectors_iterator vectors_end () const
NumericVector< Number > & add_vector (const std::string &vec_name, const bool projections=true, const ParallelType type=PARALLEL)
void remove_vector (const std::string &vec_name)
bool & project_solution_on_reinit (void)
bool have_vector (const std::string &vec_name) const
const NumericVector< Number > * request_vector (const std::string &vec_name) const
NumericVector< Number > * request_vector (const std::string &vec_name)
const NumericVector< Number > * request_vector (const unsigned int vec_num) const
NumericVector< Number > * request_vector (const unsigned int vec_num)
const NumericVector< Number > & get_vector (const std::string &vec_name) const
NumericVector< Number > & get_vector (const std::string &vec_name)
const NumericVector< Number > & get_vector (const unsigned int vec_num) const
NumericVector< Number > & get_vector (const unsigned int vec_num)
const std::string & vector_name (const unsigned int vec_num) const
const std::string & vector_name (const NumericVector< Number > &vec_reference) const
void set_vector_preservation (const std::string &vec_name, bool preserve)
bool vector_preservation (const std::string &vec_name) const
NumericVector< Number > & add_adjoint_solution (unsigned int i=0)
NumericVector< Number > & get_adjoint_solution (unsigned int i=0)
const NumericVector< Number > & get_adjoint_solution (unsigned int i=0) const
NumericVector< Number > & add_sensitivity_solution (unsigned int i=0)
NumericVector< Number > & get_sensitivity_solution (unsigned int i=0)
const NumericVector< Number > & get_sensitivity_solution (unsigned int i=0) const
NumericVector< Number > & add_weighted_sensitivity_adjoint_solution (unsigned int i=0)
NumericVector< Number > & get_weighted_sensitivity_adjoint_solution (unsigned int i=0)
const NumericVector< Number > & get_weighted_sensitivity_adjoint_solution (unsigned int i=0) const
NumericVector< Number > & add_weighted_sensitivity_solution ()
NumericVector< Number > & get_weighted_sensitivity_solution ()
const NumericVector< Number > & get_weighted_sensitivity_solution () const
NumericVector< Number > & add_adjoint_rhs (unsigned int i=0)
NumericVector< Number > & get_adjoint_rhs (unsigned int i=0)
const NumericVector< Number > & get_adjoint_rhs (unsigned int i=0) const
NumericVector< Number > & add_sensitivity_rhs (unsigned int i=0)
NumericVector< Number > & get_sensitivity_rhs (unsigned int i=0)
const NumericVector< Number > & get_sensitivity_rhs (unsigned int i=0) const
unsigned int n_vectors () const
unsigned int n_vars () const
unsigned int n_variable_groups () const
unsigned int n_components () const
dof_id_type n_dofs () const
dof_id_type n_active_dofs () const
dof_id_type n_constrained_dofs () const
dof_id_type n_local_constrained_dofs () const
dof_id_type n_local_dofs () const
unsigned int add_variable (const std::string &var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=NULL)
unsigned int add_variable (const std::string &var, const Order order=FIRST, const FEFamily=LAGRANGE, const std::set< subdomain_id_type > *const active_subdomains=NULL)
unsigned int add_variables (const std::vector< std::string > &vars, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=NULL)
unsigned int add_variables (const std::vector< std::string > &vars, const Order order=FIRST, const FEFamily=LAGRANGE, const std::set< subdomain_id_type > *const active_subdomains=NULL)
const Variable & variable (unsigned int var) const
const VariableGroup & variable_group (unsigned int vg) const
bool has_variable (const std::string &var) const
const std::string & variable_name (const unsigned int i) const
unsigned short int variable_number (const std::string &var) const
void get_all_variable_numbers (std::vector< unsigned int > &all_variable_numbers) const
unsigned int variable_scalar_number (const std::string &var, unsigned int component) const
unsigned int variable_scalar_number (unsigned int var_num, unsigned int component) const
const FEType & variable_type (const unsigned int i) const
const FEType & variable_type (const std::string &var) const
bool identify_variable_groups () const
void identify_variable_groups (const bool)
Real calculate_norm (const NumericVector< Number > &v, unsigned int var=0, FEMNormType norm_type=L2) const
Real calculate_norm (const NumericVector< Number > &v, const SystemNorm &norm) const
void read_header (Xdr &io, const std::string &version, const bool read_header=true, const bool read_additional_data=true, const bool read_legacy_format=false)
void read_legacy_data (Xdr &io, const bool read_additional_data=true)
void read_serialized_data (Xdr &io, const bool read_additional_data=true)
dof_id_type read_serialized_vectors (Xdr &io, const std::vector< NumericVector< Number > * > &vectors) const
void read_parallel_data (Xdr &io, const bool read_additional_data)
void write_header (Xdr &io, const std::string &version, const bool write_additional_data) const
void write_serialized_data (Xdr &io, const bool write_additional_data=true) const
dof_id_type write_serialized_vectors (Xdr &io, const std::vector< const NumericVector< Number > * > &vectors) const
void write_parallel_data (Xdr &io, const bool write_additional_data) const
std::string get_info () const
void attach_init_function (void fptr(EquationSystems &es, const std::string &name))
void attach_init_object (Initialization &init)
void attach_assemble_function (void fptr(EquationSystems &es, const std::string &name))
void attach_assemble_object (Assembly &assemble)
void attach_constraint_function (void fptr(EquationSystems &es, const std::string &name))
void attach_constraint_object (Constraint &constrain)
void attach_QOI_function (void fptr(EquationSystems &es, const std::string &name, const QoISet &qoi_indices))
void attach_QOI_object (QOI &qoi)
void attach_QOI_derivative (void fptr(EquationSystems &es, const std::string &name, const QoISet &qoi_indices))
void attach_QOI_derivative_object (QOIDerivative &qoi_derivative)
virtual void user_initialization ()
virtual void user_assembly ()
virtual void user_constrain ()
virtual void user_QOI (const QoISet &qoi_indices)
virtual void user_QOI_derivative (const QoISet &qoi_indices)
virtual void re_update ()
virtual void restrict_vectors ()
virtual void prolong_vectors ()
Number current_solution (const dof_id_type global_dof_number) const
Number point_value (unsigned int var, const Point &p, const bool insist_on_success=true) const
Number point_value (unsigned int var, const Point &p, const Elem &e) const
Gradient point_gradient (unsigned int var, const Point &p, const bool insist_on_success=true) const
Gradient point_gradient (unsigned int var, const Point &p, const Elem &e) const
Tensor point_hessian (unsigned int var, const Point &p, const bool insist_on_success=true) const
Tensor point_hessian (unsigned int var, const Point &p, const Elem &e) const
void local_dof_indices (const unsigned int var, std::set< dof_id_type > &var_indices) const
void zero_variable (NumericVector< Number > &v, unsigned int var_num) const
void initialize_parameters (const RBParameters &mu_min_in, const RBParameters &mu_max_in, const RBParameters &mu_in)
void initialize_parameters (const RBParametrized &rb_parametrized)
virtual void initialize_parameters (const std::string &parameters_filename)
unsigned int get_n_params () const
const RBParameters & get_parameters () const
void set_parameters (const RBParameters &params)
const RBParameters & get_parameters_min () const
const RBParameters & get_parameters_max () const
Real get_parameter_min (const std::string &param_name) const
Real get_parameter_max (const std::string &param_name) const
void print_parameters () const
void write_parameter_ranges_to_file (const std::string &file_name, const bool write_binary)
void read_parameter_ranges_from_file (const std::string &file_name, const bool read_binary)

Static Public Member Functions

static void print_info (std::ostream &out=libMesh::out)
static void print_info (std::ostream &out=libMesh::out)
static AutoPtr< DirichletBoundary > build_zero_dirichlet_boundary_object ()
static std::string get_info ()
static std::string get_info ()
static unsigned int n_objects ()
static unsigned int n_objects ()
static void enable_print_counter_info ()
static void enable_print_counter_info ()
static void disable_print_counter_info ()
static void disable_print_counter_info ()

Public Attributes

BEST_FIT_TYPE best_fit_type_flag
CouplingMatrix _coupling_matrix
std::vector< Real > training_error_bounds
AutoPtr< SparseMatrix< Number > > inner_product_matrix
AutoPtr< SparseMatrix< Number > > non_dirichlet_inner_product_matrix
AutoPtr< SparseMatrix< Number > > constraint_matrix
std::vector< Number > truth_outputs
std::vector< std::vector
< Number > > 		output_dual_innerprods
std::vector< NumericVector
< Number > * > 		Fq_representor
std::vector< Number > Fq_representor_innerprods
bool constrained_problem
bool single_matrix_mode
bool reuse_preconditioner
bool use_relative_bound_in_greedy
bool exit_on_repeated_greedy_parameters
bool write_data_during_training
bool impose_internal_dirichlet_BCs
bool impose_internal_fluxes
bool compute_RB_inner_product
bool store_non_dirichlet_operators
bool enforce_constraints_exactly
bool use_empty_rb_solve_in_greedy
AutoPtr< LinearSolver< Number > > linear_solver
SparseMatrix< Number > * matrix
NumericVector< Number > * rhs
bool assemble_before_solve
bool use_fixed_solution
int extra_quadrature_order
AutoPtr< NumericVector< Number > > solution
AutoPtr< NumericVector< Number > > current_local_solution
Real time
std::vector< Number > qoi
bool verbose_mode

Protected Types

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

Protected Member Functions

virtual void init_data ()
virtual void enrich_RB_space ()
virtual void update_system ()
virtual void update_RB_system_matrices ()
virtual Real get_RB_error_bound ()
virtual bool greedy_termination_test (Real training_greedy_error, int count)
void initialize_parametrized_functions_in_training_set ()
virtual void allocate_data_structures ()
virtual void assemble_affine_expansion ()
virtual void truth_assembly ()
virtual AutoPtr< FEMContext > build_context ()
virtual void assemble_matrix_for_output_dual_solves ()
void update_greedy_param_list ()
void add_scaled_matrix_and_vector (Number scalar, ElemAssembly *elem_assembly, SparseMatrix< Number > *input_matrix, NumericVector< Number > *input_vector, bool symmetrize=false, bool apply_dof_constraints=true)
virtual void set_context_solution_vec (NumericVector< Number > &vec)
void assemble_scaled_matvec (Number scalar, ElemAssembly *elem_assembly, NumericVector< Number > &dest, NumericVector< Number > &arg)
virtual void assemble_misc_matrices ()
virtual void assemble_all_affine_operators ()
virtual void assemble_all_affine_vectors ()
virtual void assemble_all_output_vectors ()
virtual void compute_output_dual_innerprods ()
virtual void compute_Fq_representor_innerprods (bool compute_inner_products=true)
virtual void update_residual_terms (bool compute_inner_products=true)
RBParameters get_params_from_training_set (unsigned int index)
void set_params_from_training_set (unsigned int index)
virtual void set_params_from_training_set_and_broadcast (unsigned int index)
virtual void init_matrices ()
void project_vector (NumericVector< Number > &) const
void project_vector (const NumericVector< Number > &, NumericVector< Number > &) const
void increment_constructor_count (const std::string &name)
void increment_constructor_count (const std::string &name)
void increment_destructor_count (const std::string &name)
void increment_destructor_count (const std::string &name)

Static Protected Member Functions

static void get_global_max_error_pair (std::pair< unsigned int, Real > &error_pair)
static void generate_training_parameters_random (std::map< std::string, bool > log_param_scale, std::map< std::string, NumericVector< Number > * > &training_parameters_in, unsigned int n_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, int training_parameters_random_seed=-1, bool serial_training_set=false)
static void generate_training_parameters_partially_random (const std::string &deterministic_parameter_name, const unsigned int deterministic_parameter_repeats, std::map< std::string, bool > log_param_scale, std::map< std::string, NumericVector< Number > * > &training_parameters_in, unsigned int n_deterministic_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, int training_parameters_random_seed=-1, bool serial_training_set=false)
static void generate_training_parameters_deterministic (std::map< std::string, bool > log_param_scale, std::map< std::string, NumericVector< Number > * > &training_parameters_in, unsigned int n_training_samples_in, const RBParameters &min_parameters, const RBParameters &max_parameters, bool serial_training_set=false)

Protected Attributes

bool _parametrized_functions_in_training_set_initialized
std::vector< NumericVector
< Number > * > 		_parametrized_functions_in_training_set
unsigned int Nmax
unsigned int delta_N
bool quiet_mode
bool output_dual_innerprods_computed
bool Fq_representor_innerprods_computed
bool serial_training_set
AutoPtr< NumericVector< Number > > inner_product_storage_vector
std::string alternative_solver
unsigned int _n_linear_iterations
Real _final_linear_residual
ShellMatrix< Number > * _shell_matrix
const SystemSubset * _subset
SubsetSolveMode _subset_solve_mode

Static Protected Attributes

static Counts _counts
static Counts _counts
static Threads::atomic
< unsigned int > 		_n_objects
static Threads::atomic
< unsigned int > 		_n_objects
static Threads::spin_mutex _mutex
static Threads::spin_mutex _mutex
static bool _enable_print_counter = true
static bool _enable_print_counter = true

Private Attributes

MeshFunction * _mesh_function
bool _performing_extra_greedy_step
AutoPtr< NumericVector< Number > > _ghosted_meshfunction_vector
RBAssemblyExpansion _empty_rb_assembly_expansion
std::vector< ElemAssembly * > _rb_eim_assembly_objects

Detailed Description

This class is part of the rbOOmit framework.

RBEIMConstruction implements the Construction stage of the Empirical Interpolation Method (EIM). This can be used to generate an affine approximation to non-affine operators.

Author:
David J. Knezevic, 2010

Definition at line 51 of file rb_eim_construction.h.

Member Typedef Documentation

typedef std::map<std::string, SparseMatrix<Number>* >::const_iterator libMesh::ImplicitSystem::const_matrices_iterator [inherited]

Definition at line 277 of file implicit_system.h.

typedef std::map<std::string, NumericVector<Number>* >::const_iterator libMesh::System::const_vectors_iterator [inherited]

Definition at line 716 of file system.h.

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts [protected, inherited]

Data structure to log the information. The log is identified by the class name.

Definition at line 113 of file reference_counter.h.

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts [protected, inherited]

Data structure to log the information. The log is identified by the class name.

Definition at line 113 of file reference_counter.h.

typedef std::map<std::string, SparseMatrix<Number>* >::iterator libMesh::ImplicitSystem::matrices_iterator [inherited]

Matrix iterator typedefs.

Definition at line 276 of file implicit_system.h.

typedef RBConstruction libMesh::RBEIMConstruction::Parent

The type of the parent.

Reimplemented from libMesh::RBConstruction.

Definition at line 78 of file rb_eim_construction.h.

typedef RBEIMConstruction libMesh::RBEIMConstruction::sys_type

The type of system.

Reimplemented from libMesh::RBConstruction.

Definition at line 73 of file rb_eim_construction.h.

typedef std::map<std::string, NumericVector<Number>* >::iterator libMesh::System::vectors_iterator [inherited]

Vector iterator typedefs.

Definition at line 715 of file system.h.

Member Enumeration Documentation

enum libMesh::RBEIMConstruction::BEST_FIT_TYPE
Enumerator:
PROJECTION_BEST_FIT 
EIM_BEST_FIT 

Definition at line 55 of file rb_eim_construction.h.

00055 { PROJECTION_BEST_FIT, EIM_BEST_FIT };

Constructor & Destructor Documentation

libMesh::RBEIMConstruction::RBEIMConstruction ( EquationSystems &  es,
const std::string &  name,
const unsigned int  number 
)

Constructor. Optionally initializes required data structures.

virtual libMesh::RBEIMConstruction::~RBEIMConstruction (  )  [virtual]

Destructor.

Member Function Documentation

void libMesh::System::activate (  )  [inline, inherited]

Activates the system. Only active systems are solved.

Definition at line 1874 of file system.h.

References libMesh::System::_active.

01875 {
01876   _active = true;
01877 }

bool libMesh::System::active (  )  const [inline, inherited]
Returns:
true if the system is active, false otherwise. An active system will be solved.

Definition at line 1866 of file system.h.

References libMesh::System::_active.

01867 {
01868   return _active;
01869 }

NumericVector< Number > & libMesh::System::add_adjoint_rhs ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's adjoint rhs vectors, by default the one corresponding to the first qoi. Creates the vector if it doesn't already exist.

Definition at line 1016 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::FEMSystem::assemble_qoi_derivative(), and libMesh::ExplicitSystem::assemble_qoi_derivative().

01017 {
01018   std::ostringstream adjoint_rhs_name;
01019   adjoint_rhs_name << "adjoint_rhs" << i;
01020 
01021   return this->add_vector(adjoint_rhs_name.str(), false);
01022 }

NumericVector< Number > & libMesh::System::add_adjoint_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's adjoint solution vectors, by default the one corresponding to the first qoi. Creates the vector if it doesn't already exist.

Definition at line 956 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::ImplicitSystem::adjoint_solve().

00957 {
00958   std::ostringstream adjoint_name;
00959   adjoint_name << "adjoint_solution" << i;
00960 
00961   return this->add_vector(adjoint_name.str());
00962 }

SparseMatrix< Number > & libMesh::ImplicitSystem::add_matrix ( const std::string &  mat_name  )  [inherited]

Adds the additional matrix mat_name to this system. Only allowed prior to assemble(). All additional matrices have the same sparsity pattern as the matrix used during solution. When not System but the user wants to initialize the mayor matrix, then all the additional matrices, if existent, have to be initialized by the user, too.

Definition at line 207 of file implicit_system.C.

References libMesh::ImplicitSystem::_can_add_matrices, libMesh::ImplicitSystem::_matrices, libMesh::err, and libMesh::ImplicitSystem::have_matrix().

Referenced by libMesh::ImplicitSystem::add_system_matrix(), libMesh::EigenTimeSolver::init(), and libMesh::NewmarkSystem::NewmarkSystem().

00208 {
00209   // only add matrices before initializing...
00210   if (!_can_add_matrices)
00211     {
00212       libMesh::err << "ERROR: Too late.  Cannot add matrices to the system after initialization"
00213                     << std::endl
00214                     << " any more.  You should have done this earlier."
00215                     << std::endl;
00216       libmesh_error();
00217     }
00218 
00219   // Return the matrix if it is already there.
00220   if (this->have_matrix(mat_name))
00221     return *(_matrices[mat_name]);
00222 
00223   // Otherwise build the matrix and return it.
00224   SparseMatrix<Number>* buf = SparseMatrix<Number>::build().release();
00225   _matrices.insert (std::make_pair (mat_name, buf));
00226 
00227   return *buf;
00228 }

void libMesh::RBConstruction::add_scaled_Aq ( Number  scalar,
unsigned int  q_a,
SparseMatrix< Number > *  input_matrix,
bool  symmetrize 
) [inherited]

Add the scaled q^th affine matrix to input_matrix. If symmetrize==true, then we symmetrize Aq before adding it.

void libMesh::RBConstruction::add_scaled_matrix_and_vector ( Number  scalar,
ElemAssembly *  elem_assembly,
SparseMatrix< Number > *  input_matrix,
NumericVector< Number > *  input_vector,
bool  symmetrize = false,
bool  apply_dof_constraints = true 
) [protected, inherited]

This function loops over the mesh and applies the specified interior and/or boundary assembly routines, then adds the scaled result to input_matrix and/or input_vector. If symmetrize==true then we assemble the symmetric part of the matrix, 0.5*(A + A^T)

NumericVector< Number > & libMesh::System::add_sensitivity_rhs ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's sensitivity rhs vectors, by default the one corresponding to the first parameter. Creates the vector if it doesn't already exist.

Definition at line 1046 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::ImplicitSystem::assemble_residual_derivatives().

01047 {
01048   std::ostringstream sensitivity_rhs_name;
01049   sensitivity_rhs_name << "sensitivity_rhs" << i;
01050 
01051   return this->add_vector(sensitivity_rhs_name.str(), false);
01052 }

NumericVector< Number > & libMesh::System::add_sensitivity_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's solution sensitivity vectors, by default the one corresponding to the first parameter. Creates the vector if it doesn't already exist.

Definition at line 905 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::ImplicitSystem::sensitivity_solve().

00906 {
00907   std::ostringstream sensitivity_name;
00908   sensitivity_name << "sensitivity_solution" << i;
00909 
00910   return this->add_vector(sensitivity_name.str());
00911 }

unsigned int libMesh::System::add_variable ( const std::string &  var,
const Order  order = FIRST,
const FEFamily  family = LAGRANGE,
const std::set< subdomain_id_type > *const   active_subdomains = NULL 
) [inherited]

Adds the variable var to the list of variables for this system. Same as before, but assumes LAGRANGE as default value for FEType.family.

Definition at line 1155 of file system.C.

References libMesh::System::add_variable().

01159 {
01160   return this->add_variable(var,
01161                             FEType(order, family),
01162                             active_subdomains);
01163 }

unsigned int libMesh::System::add_variable ( const std::string &  var,
const FEType &  type,
const std::set< subdomain_id_type > *const   active_subdomains = NULL 
) [inherited]

Adds the variable var to the list of variables for this system. Returns the index number for the new variable.

Definition at line 1076 of file system.C.

References libMesh::System::_variable_groups, libMesh::System::_variable_numbers, libMesh::System::_variables, libMesh::System::add_variables(), libMesh::err, libMesh::System::identify_variable_groups(), libMesh::System::n_variable_groups(), libMesh::System::n_vars(), libMesh::System::number(), libMesh::System::variable_name(), and libMesh::System::variable_type().

Referenced by libMesh::System::add_variable(), libMesh::ErrorVector::plot_error(), and libMesh::System::read_header().

01079 {
01080   // Make sure the variable isn't there already
01081   // or if it is, that it's the type we want
01082   for (unsigned int v=0; v<this->n_vars(); v++)
01083     if (this->variable_name(v) == var)
01084       {
01085         if (this->variable_type(v) == type)
01086           return _variables[v].number();
01087 
01088         libMesh::err << "ERROR: incompatible variable "
01089                      << var
01090                      << " has already been added for this system!"
01091                      << std::endl;
01092         libmesh_error();
01093       }
01094 
01095   // Optimize for VariableGroups here - if the user is adding multiple
01096   // variables of the same FEType and subdomain restriction, catch
01097   // that here and add them as members of the same VariableGroup.
01098   //
01099   // start by setting this flag to whatever the user has requested
01100   // and then consider the conditions which should negate it.
01101   bool should_be_in_vg = this->identify_variable_groups();
01102 
01103   // No variable groups, nothing to add to
01104   if (!this->n_variable_groups())
01105     should_be_in_vg = false;
01106 
01107   else
01108     {
01109       VariableGroup &vg(_variable_groups.back());
01110 
01111       // get a pointer to their subdomain restriction, if any.
01112       const std::set<subdomain_id_type> * const
01113         their_active_subdomains (vg.implicitly_active() ?
01114                                  NULL : &vg.active_subdomains());
01115 
01116       // Different types?
01117       if (vg.type() != type)
01118         should_be_in_vg = false;
01119 
01120       // they are restricted, we aren't?
01121       if (their_active_subdomains && !active_subdomains)
01122         should_be_in_vg = false;
01123 
01124       // they aren't restriced, we are?
01125       if (!their_active_subdomains && active_subdomains)
01126         should_be_in_vg = false;
01127 
01128       if (their_active_subdomains && active_subdomains)
01129         // restricted to different sets?
01130         if (*their_active_subdomains != *active_subdomains)
01131           should_be_in_vg = false;
01132 
01133       // OK, after all that, append the variable to the vg if none of the conditions
01134       // were violated
01135       if (should_be_in_vg)
01136         {
01137           const unsigned int curr_n_vars = this->n_vars();
01138 
01139           vg.append (var);
01140 
01141           _variables.push_back(vg(vg.n_variables()-1));
01142           _variable_numbers[var] = curr_n_vars;
01143           return curr_n_vars;
01144         }
01145     }
01146 
01147   // otherwise, fall back to adding a single variable group
01148   return this->add_variables (std::vector<std::string>(1, var),
01149                               type,
01150                               active_subdomains);
01151 }

unsigned int libMesh::System::add_variables ( const std::vector< std::string > &  vars,
const Order  order = FIRST,
const FEFamily  family = LAGRANGE,
const std::set< subdomain_id_type > *const   active_subdomains = NULL 
) [inherited]

Adds the variable var to the list of variables for this system. Same as before, but assumes LAGRANGE as default value for FEType.family.

Definition at line 1220 of file system.C.

References libMesh::System::add_variables().

01224 {
01225   return this->add_variables(vars,
01226                              FEType(order, family),
01227                              active_subdomains);
01228 }

unsigned int libMesh::System::add_variables ( const std::vector< std::string > &  vars,
const FEType &  type,
const std::set< subdomain_id_type > *const   active_subdomains = NULL 
) [inherited]

Adds the variable var to the list of variables for this system. Returns the index number for the new variable.

Definition at line 1167 of file system.C.

References libMesh::System::_variable_groups, libMesh::System::_variable_numbers, libMesh::System::_variables, libMesh::err, libMesh::System::n_components(), libMesh::System::n_vars(), libMesh::System::number(), libMesh::System::variable_name(), and libMesh::System::variable_type().

Referenced by libMesh::System::add_variable(), and libMesh::System::add_variables().

01170 {
01171   // Make sure the variable isn't there already
01172   // or if it is, that it's the type we want
01173   for (unsigned int ov=0; ov<vars.size(); ov++)
01174     for (unsigned int v=0; v<this->n_vars(); v++)
01175       if (this->variable_name(v) == vars[ov])
01176         {
01177           if (this->variable_type(v) == type)
01178             return _variables[v].number();
01179 
01180           libMesh::err << "ERROR: incompatible variable "
01181                        << vars[ov]
01182                        << " has already been added for this system!"
01183                        << std::endl;
01184           libmesh_error();
01185         }
01186 
01187   const unsigned int curr_n_vars = this->n_vars();
01188 
01189   const unsigned int next_first_component = this->n_components();
01190 
01191   // Add the variable group to the list
01192   _variable_groups.push_back((active_subdomains == NULL) ?
01193                              VariableGroup(this, vars, curr_n_vars,
01194                                            next_first_component, type) :
01195                              VariableGroup(this, vars, curr_n_vars,
01196                                            next_first_component, type, *active_subdomains));
01197 
01198   const VariableGroup &vg (_variable_groups.back());
01199 
01200   // Add each component of the group individually
01201   for (unsigned int v=0; v<vars.size(); v++)
01202     {
01203       _variables.push_back (vg(v));
01204       _variable_numbers[vars[v]] = curr_n_vars+v;
01205     }
01206 
01207   libmesh_assert_equal_to ((curr_n_vars+vars.size()), this->n_vars());
01208 
01209   // BSK - Defer this now to System::init_data() so we can detect
01210   // VariableGroups 12/28/2012
01211   // // Add the variable group to the _dof_map
01212   // _dof_map->add_variable_group (vg);
01213 
01214   // Return the number of the new variable
01215   return curr_n_vars+vars.size()-1;
01216 }

NumericVector< Number > & libMesh::System::add_vector ( const std::string &  vec_name,
const bool  projections = true,
const ParallelType  type = PARALLEL 
) [inherited]

Adds the additional vector vec_name to this system. All the additional vectors are similarly distributed, like the solution, and inititialized to zero.

By default vectors added by add_vector are projected to changed grids by reinit(). To zero them instead (more efficient), pass "false" as the second argument

Definition at line 675 of file system.C.

References libMesh::System::_can_add_vectors, libMesh::System::_dof_map, libMesh::System::_vector_projections, libMesh::System::_vector_types, libMesh::System::_vectors, libMeshEnums::GHOSTED, libMesh::System::have_vector(), libMesh::NumericVector< T >::init(), libMesh::System::n_dofs(), and libMesh::System::n_local_dofs().

Referenced by libMesh::System::add_adjoint_rhs(), libMesh::System::add_adjoint_solution(), libMesh::System::add_sensitivity_rhs(), libMesh::System::add_sensitivity_solution(), libMesh::ExplicitSystem::add_system_rhs(), libMesh::System::add_weighted_sensitivity_adjoint_solution(), libMesh::System::add_weighted_sensitivity_solution(), libMesh::UnsteadySolver::init(), libMesh::ContinuationSystem::init_data(), libMesh::NewmarkSystem::NewmarkSystem(), libMesh::System::read_header(), libMesh::FrequencySystem::set_frequencies(), libMesh::FrequencySystem::set_frequencies_by_range(), and libMesh::FrequencySystem::set_frequencies_by_steps().

00678 {
00679   // Return the vector if it is already there.
00680   if (this->have_vector(vec_name))
00681     return *(_vectors[vec_name]);
00682 
00683   // Otherwise build the vector
00684   NumericVector<Number>* buf = NumericVector<Number>::build().release();
00685   _vectors.insert (std::make_pair (vec_name, buf));
00686   _vector_projections.insert (std::make_pair (vec_name, projections));
00687 
00688   _vector_types.insert (std::make_pair (vec_name, type));
00689 
00690   // Initialize it if necessary
00691   if (!_can_add_vectors)
00692   {
00693     if(type == GHOSTED)
00694     {
00695 #ifdef LIBMESH_ENABLE_GHOSTED
00696       buf->init (this->n_dofs(), this->n_local_dofs(),
00697                  _dof_map->get_send_list(), false,
00698                  GHOSTED);
00699 #else
00700       std::cerr << "Cannot initialize ghosted vectors when they are not enabled." << std::endl;
00701       libmesh_error();
00702 #endif
00703     }
00704     else
00705       buf->init (this->n_dofs(), this->n_local_dofs(), false, type);
00706   }
00707 
00708   return *buf;
00709 }

NumericVector< Number > & libMesh::System::add_weighted_sensitivity_adjoint_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's weighted sensitivity adjoint solution vectors, by default the one corresponding to the first qoi. Creates the vector if it doesn't already exist.

Definition at line 986 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

00987 {
00988   std::ostringstream adjoint_name;
00989   adjoint_name << "weighted_sensitivity_adjoint_solution" << i;
00990 
00991   return this->add_vector(adjoint_name.str());
00992 }

NumericVector< Number > & libMesh::System::add_weighted_sensitivity_solution (  )  [inherited]
Returns:
a reference to the solution of the last weighted sensitivity solve Creates the vector if it doesn't already exist.

Definition at line 935 of file system.C.

References libMesh::System::add_vector().

Referenced by libMesh::ImplicitSystem::weighted_sensitivity_solve().

00936 {
00937   return this->add_vector("weighted_sensitivity_solution");
00938 }

void libMesh::ImplicitSystem::adjoint_qoi_parameter_sensitivity ( const QoISet &  qoi_indices,
const ParameterVector &  parameters,
SensitivityData &  sensitivities 
) [virtual, inherited]

Solves for the derivative of each of the system's quantities of interest q in qoi[qoi_indices] with respect to each parameter in parameters, placing the result for qoi i and parameter j into sensitivities[i][j].

Uses adjoint_solve() and the adjoint sensitivity method.

Currently uses finite differenced derivatives (partial q / partial p) and (partial R / partial p).

Reimplemented from libMesh::System.

Definition at line 696 of file implicit_system.C.

References libMesh::SensitivityData::allocate_data(), libMesh::QoISet::has_index(), libMesh::Real, libMesh::ParameterVector::size(), and libMesh::TOLERANCE.

00699 {
00700   const unsigned int Np = libmesh_cast_int<unsigned int>
00701     (parameters.size());
00702   const unsigned int Nq = libmesh_cast_int<unsigned int>
00703     (qoi.size());
00704 
00705   // We currently get partial derivatives via central differencing
00706   const Real delta_p = TOLERANCE;
00707 
00708   // An introduction to the problem:
00709   //
00710   // Residual R(u(p),p) = 0
00711   // partial R / partial u = J = system matrix
00712   //
00713   // This implies that:
00714   // d/dp(R) = 0
00715   // (partial R / partial p) +
00716   // (partial R / partial u) * (partial u / partial p) = 0
00717 
00718   // We first do an adjoint solve:
00719   // J^T * z = (partial q / partial u)
00720   // if we havent already or dont have an initial condition for the adjoint
00721   if (!this->is_adjoint_already_solved())
00722     {
00723       this->adjoint_solve(qoi_indices);
00724     }
00725 
00726   // Get ready to fill in senstivities:
00727   sensitivities.allocate_data(qoi_indices, *this, parameters);
00728 
00729   // We use the identities:
00730   // dq/dp = (partial q / partial p) + (partial q / partial u) *
00731   //         (partial u / partial p)
00732   // dq/dp = (partial q / partial p) + (J^T * z) *
00733   //         (partial u / partial p)
00734   // dq/dp = (partial q / partial p) + z * J *
00735   //         (partial u / partial p)
00736 
00737   // Leading to our final formula:
00738   // dq/dp = (partial q / partial p) - z * (partial R / partial p)
00739 
00740   for (unsigned int j=0; j != Np; ++j)
00741     {
00742       // (partial q / partial p) ~= (q(p+dp)-q(p-dp))/(2*dp)
00743       // (partial R / partial p) ~= (rhs(p+dp) - rhs(p-dp))/(2*dp)
00744 
00745       Number old_parameter = *parameters[j];
00746       // Number old_qoi = this->qoi;
00747 
00748       *parameters[j] = old_parameter - delta_p;
00749       this->assemble_qoi(qoi_indices);
00750       std::vector<Number> qoi_minus = this->qoi;
00751 
00752       this->assembly(true, false);
00753       this->rhs->close();
00754 
00755       // FIXME - this can and should be optimized to avoid the clone()
00756       AutoPtr<NumericVector<Number> > partialR_partialp = this->rhs->clone();
00757       *partialR_partialp *= -1;
00758 
00759       *parameters[j] = old_parameter + delta_p;
00760       this->assemble_qoi(qoi_indices);
00761       std::vector<Number>& qoi_plus = this->qoi;
00762 
00763       std::vector<Number> partialq_partialp(Nq, 0);
00764       for (unsigned int i=0; i != Nq; ++i)
00765         if (qoi_indices.has_index(i))
00766           partialq_partialp[i] = (qoi_plus[i] - qoi_minus[i]) / (2.*delta_p);
00767 
00768       this->assembly(true, false);
00769       this->rhs->close();
00770       *partialR_partialp += *this->rhs;
00771       *partialR_partialp /= (2.*delta_p);
00772 
00773       // Don't leave the parameter changed
00774       *parameters[j] = old_parameter;
00775 
00776       for (unsigned int i=0; i != Nq; ++i)
00777         if (qoi_indices.has_index(i))
00778           sensitivities[i][j] = partialq_partialp[i] -
00779             partialR_partialp->dot(this->get_adjoint_solution(i));
00780     }
00781 
00782   // All parameters have been reset.
00783   // We didn't cache the original rhs or matrix for memory reasons,
00784   // but we can restore them to a state consistent solution -
00785   // principle of least surprise.
00786   this->assembly(true, true);
00787   this->rhs->close();
00788   this->matrix->close();
00789   this->assemble_qoi(qoi_indices);
00790 }

std::pair< unsigned int, Real > libMesh::ImplicitSystem::adjoint_solve ( const QoISet &  qoi_indices = QoISet()  )  [virtual, inherited]

Assembles & solves the linear system (dR/du)^T*z = dq/du, for those quantities of interest q specified by qoi_indices.

Leave qoi_indices empty to solve all adjoint problems.

Returns a pair with the total number of linear iterations performed and the (sum of the) final residual norms

Reimplemented from libMesh::System.

Reimplemented in libMesh::DifferentiableSystem.

Definition at line 382 of file implicit_system.C.

References libMesh::System::add_adjoint_solution(), libMesh::LinearSolver< T >::adjoint_solve(), libMesh::System::assemble_before_solve, libMesh::ExplicitSystem::assemble_qoi_derivative(), libMesh::ImplicitSystem::assembly(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_adjoint_rhs(), libMesh::System::get_adjoint_solution(), libMesh::System::get_dof_map(), libMesh::ImplicitSystem::get_linear_solve_parameters(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::QoISet::has_index(), libMesh::ImplicitSystem::matrix, libMesh::System::qoi, and libMesh::ImplicitSystem::release_linear_solver().

00383 {
00384   // Log how long the linear solve takes.
00385   START_LOG("adjoint_solve()", "ImplicitSystem");
00386 
00387   if (this->assemble_before_solve)
00388     // Assemble the linear system
00389     this->assembly (/* get_residual = */ false,
00390                     /* get_jacobian = */ true);
00391 
00392   // The adjoint problem is linear
00393   LinearSolver<Number> *linear_solver = this->get_linear_solver();
00394 
00395   // Reset and build the RHS from the QOI derivative
00396   this->assemble_qoi_derivative(qoi_indices);
00397 
00398   // Our iteration counts and residuals will be sums of the individual
00399   // results
00400   std::pair<unsigned int, Real> solver_params =
00401     this->get_linear_solve_parameters();
00402   std::pair<unsigned int, Real> totalrval = std::make_pair(0,0.0);
00403 
00404   for (unsigned int i=0; i != this->qoi.size(); ++i)
00405     if (qoi_indices.has_index(i))
00406       {
00407         const std::pair<unsigned int, Real> rval =
00408           linear_solver->adjoint_solve (*matrix, this->add_adjoint_solution(i),
00409                                         this->get_adjoint_rhs(i),
00410                                         solver_params.second,
00411                                         solver_params.first);
00412 
00413             totalrval.first  += rval.first;
00414             totalrval.second += rval.second;
00415       }
00416 
00417   this->release_linear_solver(linear_solver);
00418 
00419   // The linear solver may not have fit our constraints exactly
00420 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00421   for (unsigned int i=0; i != this->qoi.size(); ++i)
00422     if (qoi_indices.has_index(i))
00423       this->get_dof_map().enforce_constraints_exactly
00424         (*this, &this->get_adjoint_solution(i),
00425          /* homogeneous = */ true);
00426 #endif
00427 
00428   // Stop logging the nonlinear solve
00429   STOP_LOG("adjoint_solve()", "ImplicitSystem");
00430 
00431   return totalrval;
00432 }

virtual void libMesh::RBConstruction::allocate_data_structures (  )  [protected, virtual, inherited]

Helper function that actually allocates all the data structures required by this class.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::LinearImplicitSystem::assemble (  )  [inline, virtual, inherited]

Prepares matrix and _dof_map for matrix assembly. Does not actually assemble anything. For matrix assembly, use the assemble() in derived classes. Should be overloaded in derived classes.

Reimplemented from libMesh::ImplicitSystem.

Reimplemented in libMesh::FrequencySystem, and libMesh::NewmarkSystem.

Definition at line 104 of file linear_implicit_system.h.

Referenced by libMesh::LinearImplicitSystem::assembly(), and libMesh::LinearImplicitSystem::solve().

00104 { ImplicitSystem::assemble(); }

virtual void libMesh::RBConstruction::assemble_affine_expansion (  )  [protected, virtual, inherited]

Assemble and store the Dirichlet dof lists, the affine and output vectors. Optionally assemble and store all the affine matrices if we are not in low-memory mode.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::RBConstruction::assemble_all_affine_operators (  )  [protected, virtual, inherited]

Assemble and store all Q_a affine operators as well as the inner-product matrix.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::RBConstruction::assemble_all_affine_vectors (  )  [protected, virtual, inherited]

Assemble and store the affine RHS vectors.

virtual void libMesh::RBConstruction::assemble_all_output_vectors (  )  [protected, virtual, inherited]

Assemble and store the output vectors.

void libMesh::RBConstruction::assemble_and_add_constraint_matrix ( SparseMatrix< Number > *  input_matrix  )  [inherited]

Assemble the constraint matrix and add it to input_matrix.

void libMesh::RBConstruction::assemble_Aq_matrix ( unsigned int  q,
SparseMatrix< Number > *  input_matrix,
bool  apply_dof_constraints = true 
) [inherited]

Assemble the q^th affine matrix and store it in input_matrix.

void libMesh::RBConstruction::assemble_constraint_matrix ( SparseMatrix< Number > *  input_matrix  )  [inherited]

Assemble the constraint matrix and store it in input_matrix.

void libMesh::RBConstruction::assemble_Fq_vector ( unsigned int  q,
NumericVector< Number > *  input_vector,
bool  apply_dof_constraints = true 
) [inherited]

Assemble the q^th affine vector and store it in input_matrix.

void libMesh::RBConstruction::assemble_inner_product_matrix ( SparseMatrix< Number > *  input_matrix,
bool  apply_dof_constraints = true 
) [inherited]

Assemble the inner product matrix and store it in input_matrix.

virtual void libMesh::RBConstruction::assemble_matrix_for_output_dual_solves (  )  [protected, virtual, inherited]

Define the matrix assembly for the output residual dual norm solves. By default we use the inner product matrix for steady state problems.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::RBConstruction::assemble_misc_matrices (  )  [protected, virtual, inherited]

Assemble and store all the inner-product matrix, the constraint matrix (for constrained problems) and the mass matrix (for time-dependent problems).

Reimplemented in libMesh::TransientRBConstruction.

void libMesh::ExplicitSystem::assemble_qoi ( const QoISet &  qoi_indices = QoISet()  )  [virtual, inherited]

Prepares qoi for quantity of interest assembly, then calls user qoi function. Can be overloaded in derived classes.

Reimplemented from libMesh::System.

Reimplemented in libMesh::FEMSystem.

Definition at line 89 of file explicit_system.C.

References libMesh::System::assemble_qoi(), libMesh::QoISet::has_index(), and libMesh::System::qoi.

00090 {
00091   // The user quantity of interest assembly gets to expect to
00092   // accumulate on initially zero values
00093   for (unsigned int i=0; i != qoi.size(); ++i)
00094     if (qoi_indices.has_index(i))
00095       qoi[i] = 0;
00096 
00097   Parent::assemble_qoi (qoi_indices);
00098 }

void libMesh::ExplicitSystem::assemble_qoi_derivative ( const QoISet &  qoi_indices = QoISet()  )  [virtual, inherited]

Prepares adjoint_rhs for quantity of interest derivative assembly, then calls user qoi derivative function. Can be overloaded in derived classes.

Reimplemented from libMesh::System.

Reimplemented in libMesh::FEMSystem.

Definition at line 102 of file explicit_system.C.

References libMesh::System::add_adjoint_rhs(), libMesh::System::assemble_qoi_derivative(), libMesh::QoISet::has_index(), libMesh::System::qoi, and libMesh::NumericVector< T >::zero().

Referenced by libMesh::ImplicitSystem::adjoint_solve(), and libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

00103 {
00104   // The user quantity of interest derivative assembly gets to expect
00105   // to accumulate on initially zero vectors
00106   for (unsigned int i=0; i != qoi.size(); ++i)
00107     if (qoi_indices.has_index(i))
00108       this->add_adjoint_rhs(i).zero();
00109 
00110   Parent::assemble_qoi_derivative (qoi_indices);
00111 }

void libMesh::ImplicitSystem::assemble_residual_derivatives ( const ParameterVector &  parameters  )  [virtual, inherited]

Residual parameter derivative function.

Uses finite differences by default.

This will assemble the sensitivity rhs vectors to hold -(partial R / partial p_i), making them ready to solve the forward sensitivity equation.

Can be overloaded in derived classes.

Reimplemented from libMesh::System.

Definition at line 660 of file implicit_system.C.

References libMesh::System::add_sensitivity_rhs(), libMesh::ImplicitSystem::assembly(), libMesh::NumericVector< T >::close(), libMesh::Real, libMesh::ExplicitSystem::rhs, libMesh::ParameterVector::size(), and libMesh::TOLERANCE.

Referenced by libMesh::ImplicitSystem::sensitivity_solve().

00661 {
00662   const unsigned int Np = libmesh_cast_int<unsigned int>
00663     (parameters.size());
00664   Real deltap = TOLERANCE;
00665 
00666   for (unsigned int p=0; p != Np; ++p)
00667     {
00668       NumericVector<Number> &sensitivity_rhs = this->add_sensitivity_rhs(p);
00669 
00670       // Approximate -(partial R / partial p) by
00671       // (R(p-dp) - R(p+dp)) / (2*dp)
00672 
00673       Number old_parameter = *parameters[p];
00674       *parameters[p] -= deltap;
00675 
00676       this->assembly(true, false);
00677       this->rhs->close();
00678       sensitivity_rhs = *this->rhs;
00679 
00680       *parameters[p] = old_parameter + deltap;
00681 
00682       this->assembly(true, false);
00683       this->rhs->close();
00684 
00685       sensitivity_rhs -= *this->rhs;
00686       sensitivity_rhs /= (2*deltap);
00687       sensitivity_rhs.close();
00688 
00689       *parameters[p] = old_parameter;
00690     }
00691 }

void libMesh::RBConstruction::assemble_scaled_matvec ( Number  scalar,
ElemAssembly *  elem_assembly,
NumericVector< Number > &  dest,
NumericVector< Number > &  arg 
) [protected, inherited]

This function loops over the mesh and assembles the matrix-vector product and stores the scaled result in dest.

void libMesh::LinearImplicitSystem::assembly ( bool  get_residual,
bool  get_jacobian 
) [virtual, inherited]

Assembles a residual in rhs and/or a jacobian in matrix, as requested.

Reimplemented from libMesh::ImplicitSystem.

Definition at line 373 of file linear_implicit_system.C.

References libMesh::NumericVector< T >::add_vector(), libMesh::LinearImplicitSystem::assemble(), libMesh::SparseMatrix< T >::close(), libMesh::NumericVector< T >::close(), libMesh::ImplicitSystem::matrix, libMesh::ExplicitSystem::rhs, and libMesh::System::solution.

00375 {
00376   // Residual R(u(p),p) := A(p)*u(p) - b(p)
00377   // partial R / partial u = A
00378 
00379   this->assemble();
00380   this->rhs->close();
00381   this->matrix->close();
00382 
00383   *(this->rhs) *= -1.0;
00384   this->rhs->add_vector(*this->solution, *this->matrix);
00385 }

void libMesh::System::attach_assemble_function ( void   fptrEquationSystems &es,const std::string &name  )  [inherited]

Register a user function to use in assembling the system matrix and RHS.

Definition at line 1753 of file system.C.

References libMesh::System::_assemble_system_function, libMesh::System::_assemble_system_object, and libMesh::out.

01755 {
01756   libmesh_assert(fptr);
01757 
01758   if (_assemble_system_object != NULL)
01759     {
01760       libmesh_here();
01761       libMesh::out << "WARNING:  Cannot specify both assembly function and object!"
01762                    << std::endl;
01763 
01764       _assemble_system_object = NULL;
01765     }
01766 
01767   _assemble_system_function = fptr;
01768 }

void libMesh::System::attach_assemble_object ( System::Assembly &  assemble_in  )  [inherited]

Register a user object to use in assembling the system matrix and RHS.

Definition at line 1772 of file system.C.

References libMesh::System::_assemble_system_function, libMesh::System::_assemble_system_object, and libMesh::out.

01773 {
01774   if (_assemble_system_function != NULL)
01775     {
01776       libmesh_here();
01777       libMesh::out << "WARNING:  Cannot specify both assembly object and function!"
01778                    << std::endl;
01779 
01780       _assemble_system_function = NULL;
01781     }
01782 
01783   _assemble_system_object = &assemble_in;
01784 }

void libMesh::System::attach_constraint_function ( void   fptrEquationSystems &es,const std::string &name  )  [inherited]

Register a user function for imposing constraints.

Definition at line 1788 of file system.C.

References libMesh::System::_constrain_system_function, libMesh::System::_constrain_system_object, and libMesh::out.

01790 {
01791   libmesh_assert(fptr);
01792 
01793   if (_constrain_system_object != NULL)
01794     {
01795       libmesh_here();
01796       libMesh::out << "WARNING:  Cannot specify both constraint function and object!"
01797                    << std::endl;
01798 
01799       _constrain_system_object = NULL;
01800     }
01801 
01802   _constrain_system_function = fptr;
01803 }

void libMesh::System::attach_constraint_object ( System::Constraint &  constrain  )  [inherited]

Register a user object for imposing constraints.

Definition at line 1807 of file system.C.

References libMesh::System::_constrain_system_function, libMesh::System::_constrain_system_object, and libMesh::out.

01808 {
01809   if (_constrain_system_function != NULL)
01810     {
01811       libmesh_here();
01812       libMesh::out << "WARNING:  Cannot specify both constraint object and function!"
01813                    << std::endl;
01814 
01815       _constrain_system_function = NULL;
01816     }
01817 
01818   _constrain_system_object = &constrain;
01819 }

void libMesh::System::attach_init_function ( void   fptrEquationSystems &es,const std::string &name  )  [inherited]

Register a user function to use in initializing the system.

Definition at line 1718 of file system.C.

References libMesh::System::_init_system_function, libMesh::System::_init_system_object, and libMesh::out.

01720 {
01721   libmesh_assert(fptr);
01722 
01723   if (_init_system_object != NULL)
01724     {
01725       libmesh_here();
01726       libMesh::out << "WARNING:  Cannot specify both initialization function and object!"
01727                    << std::endl;
01728 
01729       _init_system_object = NULL;
01730     }
01731 
01732   _init_system_function = fptr;
01733 }

void libMesh::System::attach_init_object ( System::Initialization &  init_in  )  [inherited]

Register a user class to use to initialize the system. Note this is exclusive with the attach_init_function.

Definition at line 1737 of file system.C.

References libMesh::System::_init_system_function, libMesh::System::_init_system_object, and libMesh::out.

01738 {
01739   if (_init_system_function != NULL)
01740     {
01741       libmesh_here();
01742       libMesh::out << "WARNING:  Cannot specify both initialization object and function!"
01743                    << std::endl;
01744 
01745       _init_system_function = NULL;
01746     }
01747 
01748   _init_system_object = &init_in;
01749 }

void libMesh::System::attach_QOI_derivative ( void   fptrEquationSystems &es, const std::string &name, const QoISet &qoi_indices  )  [inherited]

Register a user function for evaluating derivatives of a quantity of interest with respect to test functions, whose values should be placed in System::rhs

void libMesh::System::attach_QOI_derivative_object ( QOIDerivative &  qoi_derivative  )  [inherited]

Register a user object for evaluating derivatives of a quantity of interest with respect to test functions, whose values should be placed in System::rhs

Definition at line 1879 of file system.C.

References libMesh::System::_qoi_evaluate_derivative_function, libMesh::System::_qoi_evaluate_derivative_object, and libMesh::out.

01880 {
01881   if (_qoi_evaluate_derivative_function != NULL)
01882     {
01883       libmesh_here();
01884       libMesh::out << "WARNING:  Cannot specify both QOI derivative object and function!"
01885                    << std::endl;
01886 
01887       _qoi_evaluate_derivative_function = NULL;
01888     }
01889 
01890   _qoi_evaluate_derivative_object = &qoi_derivative;
01891 }

void libMesh::System::attach_QOI_function ( void   fptrEquationSystems &es, const std::string &name, const QoISet &qoi_indices  )  [inherited]

Register a user function for evaluating the quantities of interest, whose values should be placed in System::qoi

void libMesh::System::attach_QOI_object ( QOI &  qoi  )  [inherited]

Register a user object for evaluating the quantities of interest, whose values should be placed in System::qoi

Definition at line 1843 of file system.C.

References libMesh::System::_qoi_evaluate_function, libMesh::System::_qoi_evaluate_object, and libMesh::out.

01844 {
01845   if (_qoi_evaluate_function != NULL)
01846     {
01847       libmesh_here();
01848       libMesh::out << "WARNING:  Cannot specify both QOI object and function!"
01849                    << std::endl;
01850 
01851       _qoi_evaluate_function = NULL;
01852     }
01853 
01854   _qoi_evaluate_object = &qoi_in;
01855 }

void libMesh::LinearImplicitSystem::attach_shell_matrix ( ShellMatrix< Number > *  shell_matrix  )  [inherited]

This function enables the user to provide a shell matrix, i.e. a matrix that is not stored element-wise, but as a function. When you register your shell matrix using this function, calling solve() will no longer use the matrix member but the registered shell matrix instead. You can reset this behaviour to its original state by supplying a NULL pointer to this function.

Definition at line 165 of file linear_implicit_system.C.

References libMesh::LinearImplicitSystem::_shell_matrix.

Referenced by libMesh::LinearImplicitSystem::detach_shell_matrix().

00166 {
00167   _shell_matrix = shell_matrix;
00168 }

void libMesh::System::boundary_project_solution ( const std::set< boundary_id_type > &  b,
const std::vector< unsigned int > &  variables,
Number   fptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
Gradient   gptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
const Parameters &  parameters 
) [inherited]

Projects arbitrary boundary functions onto a vector of degree of freedom values for the current system. Only degrees of freedom which affect the function's trace on a boundary in the set b are affected. Only degrees of freedom associated with the variables listed in the vector variables are projected. The function value fptr and its gradient gptr are represented by function pointers. A gradient gptr is only required/used for projecting onto finite element spaces with continuous derivatives.

This method projects components of an arbitrary boundary function onto the solution via L2 projections and nodal interpolations on each element.

Definition at line 641 of file system_projection.C.

00652 {
00653   WrappedFunction<Number> f(*this, fptr, &parameters);
00654   WrappedFunction<Gradient> g(*this, gptr, &parameters);
00655   this->boundary_project_solution(b, variables, &f, &g);
00656 }

void libMesh::System::boundary_project_solution ( const std::set< boundary_id_type > &  b,
const std::vector< unsigned int > &  variables,
FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = NULL 
) [inherited]

Projects arbitrary boundary functions onto a vector of degree of freedom values for the current system. Only degrees of freedom which affect the function's trace on a boundary in the set b are affected. Only degrees of freedom associated with the variables listed in the vector variables are projected. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitary boundary function onto the solution via L2 projections and nodal interpolations on each element.

Definition at line 665 of file system_projection.C.

00669 {
00670   this->boundary_project_vector(b, variables, *solution, f, g);
00671 
00672   solution->localize(*current_local_solution);
00673 }

void libMesh::System::boundary_project_vector ( const std::set< boundary_id_type > &  b,
const std::vector< unsigned int > &  variables,
Number   fptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
Gradient   gptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
const Parameters &  parameters,
NumericVector< Number > &  new_vector 
) const [inherited]

Projects arbitrary boundary functions onto a vector of degree of freedom values for the current system. Only degrees of freedom which affect the function's trace on a boundary in the set b are affected. Only degrees of freedom associated with the variables listed in the vector variables are projected. The function value fptr and its gradient gptr are represented by function pointers. A gradient gptr is only required/used for projecting onto finite element spaces with continuous derivatives.

This method projects an arbitrary boundary function via L2 projections and nodal interpolations on each element.

Definition at line 684 of file system_projection.C.

00696 {
00697   WrappedFunction<Number> f(*this, fptr, &parameters);
00698   WrappedFunction<Gradient> g(*this, gptr, &parameters);
00699   this->boundary_project_vector(b, variables, new_vector, &f, &g);
00700 }

void libMesh::System::boundary_project_vector ( const std::set< boundary_id_type > &  b,
const std::vector< unsigned int > &  variables,
NumericVector< Number > &  new_vector,
FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = NULL 
) const [inherited]

Projects arbitrary boundary functions onto a vector of degree of freedom values for the current system. Only degrees of freedom which affect the function's trace on a boundary in the set b are affected. Only degrees of freedom associated with the variables listed in the vector variables are projected. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitrary function via L2 projections and nodal interpolations on each element.

Definition at line 707 of file system_projection.C.

References libMesh::NumericVector< T >::close(), and libMesh::Threads::parallel_for().

00712 {
00713   START_LOG ("boundary_project_vector()", "System");
00714 
00715   Threads::parallel_for
00716     (ConstElemRange (this->get_mesh().active_local_elements_begin(),
00717                      this->get_mesh().active_local_elements_end() ),
00718      BoundaryProjectSolution(b, variables, *this, f, g,
00719                              this->get_equation_systems().parameters,
00720                              new_vector)
00721     );
00722 
00723   // We don't do SCALAR dofs when just projecting the boundary, so
00724   // we're done here.
00725 
00726   new_vector.close();
00727 
00728 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00729   this->get_dof_map().enforce_constraints_exactly(*this, &new_vector);
00730 #endif
00731 
00732   STOP_LOG("boundary_project_vector()", "System");
00733 }

void libMesh::RBConstructionBase< LinearImplicitSystem >::broadcast_parameters ( unsigned int  proc_id  )  [inherited]

Broadcasts parameters on processor proc_id to all processors.

virtual AutoPtr<FEMContext> libMesh::RBConstruction::build_context (  )  [protected, virtual, inherited]

Builds a FEMContext object with enough information to do evaluations on each element.

virtual AutoPtr<ElemAssembly> libMesh::RBEIMConstruction::build_eim_assembly ( unsigned int  bf_index  )  [pure virtual]

Build an element assembly object that will access basis function bf_index. This is pure virtual, override in subclasses to specify the appropriate ElemAssembly object.

static AutoPtr<DirichletBoundary> libMesh::RBConstruction::build_zero_dirichlet_boundary_object (  )  [static, inherited]

It's helpful to be able to generate a DirichletBoundary that stores a ZeroFunction in order to impose Dirichlet boundary conditions.

Real libMesh::System::calculate_norm ( const NumericVector< Number > &  v,
const SystemNorm &  norm 
) const [inherited]
Returns:
a norm of the vector v, using component_norm and component_scale to choose and weight the norms of each variable.

Definition at line 1400 of file system.C.

References libMesh::System::_dof_map, std::abs(), libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::TypeTensor< T >::add_scaled(), libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< Real >::build(), libMesh::CommWorld, libMesh::FEType::default_quadrature_rule(), libMeshEnums::DISCRETE_L1, libMeshEnums::DISCRETE_L2, libMeshEnums::DISCRETE_L_INF, libMesh::System::discrete_var_norm(), libMesh::DofMap::dof_indices(), libMesh::AutoPtr< Tp >::get(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMeshEnums::H1, libMeshEnums::H1_SEMINORM, libMeshEnums::H2, libMeshEnums::H2_SEMINORM, libMesh::SystemNorm::is_discrete(), libMeshEnums::L1, libMesh::NumericVector< T >::l1_norm(), libMeshEnums::L2, libMesh::NumericVector< T >::l2_norm(), libMeshEnums::L_INF, libMesh::NumericVector< T >::linfty_norm(), libMesh::NumericVector< T >::localize(), libMesh::Parallel::Communicator::max(), std::max(), libMesh::MeshBase::mesh_dimension(), libMesh::System::n_vars(), libMesh::TensorTools::norm_sq(), libMesh::Real, libMeshEnums::SERIAL, libMesh::TypeTensor< T >::size(), libMesh::TypeVector< T >::size(), libMesh::NumericVector< T >::size(), libMesh::TypeTensor< T >::size_sq(), libMesh::TypeVector< T >::size_sq(), libMesh::Parallel::Communicator::sum(), libMesh::SystemNorm::type(), libMesh::DofMap::variable_type(), libMeshEnums::W1_INF_SEMINORM, libMeshEnums::W2_INF_SEMINORM, libMesh::SystemNorm::weight(), and libMesh::SystemNorm::weight_sq().

01402 {
01403   // This function must be run on all processors at once
01404   parallel_only();
01405 
01406   START_LOG ("calculate_norm()", "System");
01407 
01408   // Zero the norm before summation
01409   Real v_norm = 0.;
01410 
01411   if (norm.is_discrete())
01412     {
01413       STOP_LOG ("calculate_norm()", "System");
01414       //Check to see if all weights are 1.0 and all types are equal
01415       FEMNormType norm_type0 = norm.type(0);
01416       unsigned int check_var = 0;
01417       for (; check_var != this->n_vars(); ++check_var)
01418         if((norm.weight(check_var) != 1.0) || (norm.type(check_var) != norm_type0))
01419           break;
01420 
01421       //All weights were 1.0 so just do the full vector discrete norm
01422       if(check_var == this->n_vars())
01423         {
01424           if(norm_type0 == DISCRETE_L1)
01425             return v.l1_norm();
01426           if(norm_type0 == DISCRETE_L2)
01427             return v.l2_norm();
01428           if(norm_type0 == DISCRETE_L_INF)
01429             return v.linfty_norm();
01430           else
01431             libmesh_error();
01432         }
01433 
01434       for (unsigned int var=0; var != this->n_vars(); ++var)
01435         {
01436           // Skip any variables we don't need to integrate
01437           if (norm.weight(var) == 0.0)
01438             continue;
01439 
01440           v_norm += norm.weight(var) * discrete_var_norm(v, var, norm.type(var));
01441         }
01442 
01443       return v_norm;
01444     }
01445 
01446   // Localize the potentially parallel vector
01447   AutoPtr<NumericVector<Number> > local_v = NumericVector<Number>::build();
01448   local_v->init(v.size(), true, SERIAL);
01449   v.localize (*local_v, _dof_map->get_send_list());
01450 
01451   unsigned int dim = this->get_mesh().mesh_dimension();
01452 
01453   // I'm not sure how best to mix Hilbert norms on some variables (for
01454   // which we'll want to square then sum then square root) with norms
01455   // like L_inf (for which we'll just want to take an absolute value
01456   // and then sum).
01457   bool using_hilbert_norm = true,
01458        using_nonhilbert_norm = true;
01459 
01460   // Loop over all variables
01461   for (unsigned int var=0; var != this->n_vars(); ++var)
01462     {
01463       // Skip any variables we don't need to integrate
01464       Real norm_weight_sq = norm.weight_sq(var);
01465       if (norm_weight_sq == 0.0)
01466         continue;
01467       Real norm_weight = norm.weight(var);
01468 
01469       // Check for unimplemented norms (rather than just returning 0).
01470       FEMNormType norm_type = norm.type(var);
01471       if((norm_type==H1) ||
01472          (norm_type==H2) ||
01473          (norm_type==L2) ||
01474          (norm_type==H1_SEMINORM) ||
01475          (norm_type==H2_SEMINORM))
01476         {
01477           if (!using_hilbert_norm)
01478             libmesh_not_implemented();
01479           using_nonhilbert_norm = false;
01480         }
01481       else if ((norm_type==L1) ||
01482                (norm_type==L_INF) ||
01483                (norm_type==W1_INF_SEMINORM) ||
01484                (norm_type==W2_INF_SEMINORM))
01485         {
01486           if (!using_nonhilbert_norm)
01487             libmesh_not_implemented();
01488           using_hilbert_norm = false;
01489         }
01490       else
01491         libmesh_not_implemented();
01492 
01493       const FEType& fe_type = this->get_dof_map().variable_type(var);
01494       AutoPtr<QBase> qrule =
01495         fe_type.default_quadrature_rule (dim);
01496       AutoPtr<FEBase> fe
01497         (FEBase::build(dim, fe_type));
01498       fe->attach_quadrature_rule (qrule.get());
01499 
01500       const std::vector<Real>&               JxW = fe->get_JxW();
01501       const std::vector<std::vector<Real> >* phi = NULL;
01502       if (norm_type == H1 ||
01503           norm_type == H2 ||
01504           norm_type == L2 ||
01505           norm_type == L1 ||
01506           norm_type == L_INF)
01507         phi = &(fe->get_phi());
01508 
01509       const std::vector<std::vector<RealGradient> >* dphi = NULL;
01510       if (norm_type == H1 ||
01511           norm_type == H2 ||
01512           norm_type == H1_SEMINORM ||
01513           norm_type == W1_INF_SEMINORM)
01514         dphi = &(fe->get_dphi());
01515 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
01516       const std::vector<std::vector<RealTensor> >*   d2phi = NULL;
01517       if (norm_type == H2 ||
01518           norm_type == H2_SEMINORM ||
01519           norm_type == W2_INF_SEMINORM)
01520         d2phi = &(fe->get_d2phi());
01521 #endif
01522 
01523       std::vector<dof_id_type> dof_indices;
01524 
01525       // Begin the loop over the elements
01526       MeshBase::const_element_iterator       el     =
01527         this->get_mesh().active_local_elements_begin();
01528       const MeshBase::const_element_iterator end_el =
01529         this->get_mesh().active_local_elements_end();
01530 
01531       for ( ; el != end_el; ++el)
01532         {
01533           const Elem* elem = *el;
01534 
01535           fe->reinit (elem);
01536 
01537           this->get_dof_map().dof_indices (elem, dof_indices, var);
01538 
01539           const unsigned int n_qp = qrule->n_points();
01540 
01541           const unsigned int n_sf = libmesh_cast_int<unsigned int>
01542             (dof_indices.size());
01543 
01544           // Begin the loop over the Quadrature points.
01545           for (unsigned int qp=0; qp<n_qp; qp++)
01546             {
01547               if (norm_type == L1)
01548                 {
01549                   Number u_h = 0.;
01550                   for (unsigned int i=0; i != n_sf; ++i)
01551                     u_h += (*phi)[i][qp] * (*local_v)(dof_indices[i]);
01552                   v_norm += norm_weight *
01553                             JxW[qp] * std::abs(u_h);
01554                 }
01555 
01556               if (norm_type == L_INF)
01557                 {
01558                   Number u_h = 0.;
01559                   for (unsigned int i=0; i != n_sf; ++i)
01560                     u_h += (*phi)[i][qp] * (*local_v)(dof_indices[i]);
01561                   v_norm = std::max(v_norm, norm_weight * std::abs(u_h));
01562                 }
01563 
01564               if (norm_type == H1 ||
01565                   norm_type == H2 ||
01566                   norm_type == L2)
01567                 {
01568                   Number u_h = 0.;
01569                   for (unsigned int i=0; i != n_sf; ++i)
01570                     u_h += (*phi)[i][qp] * (*local_v)(dof_indices[i]);
01571                   v_norm += norm_weight_sq *
01572                             JxW[qp] * TensorTools::norm_sq(u_h);
01573                 }
01574 
01575               if (norm_type == H1 ||
01576                   norm_type == H2 ||
01577                   norm_type == H1_SEMINORM)
01578                 {
01579                   Gradient grad_u_h;
01580                   for (unsigned int i=0; i != n_sf; ++i)
01581                     grad_u_h.add_scaled((*dphi)[i][qp], (*local_v)(dof_indices[i]));
01582                   v_norm += norm_weight_sq *
01583                             JxW[qp] * grad_u_h.size_sq();
01584                 }
01585 
01586               if (norm_type == W1_INF_SEMINORM)
01587                 {
01588                   Gradient grad_u_h;
01589                   for (unsigned int i=0; i != n_sf; ++i)
01590                     grad_u_h.add_scaled((*dphi)[i][qp], (*local_v)(dof_indices[i]));
01591                   v_norm = std::max(v_norm, norm_weight * grad_u_h.size());
01592                 }
01593 
01594 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
01595               if (norm_type == H2 ||
01596                   norm_type == H2_SEMINORM)
01597                 {
01598                   Tensor hess_u_h;
01599                   for (unsigned int i=0; i != n_sf; ++i)
01600                     hess_u_h.add_scaled((*d2phi)[i][qp], (*local_v)(dof_indices[i]));
01601                   v_norm += norm_weight_sq *
01602                             JxW[qp] * hess_u_h.size_sq();
01603                 }
01604 
01605               if (norm_type == W2_INF_SEMINORM)
01606                 {
01607                   Tensor hess_u_h;
01608                   for (unsigned int i=0; i != n_sf; ++i)
01609                     hess_u_h.add_scaled((*d2phi)[i][qp], (*local_v)(dof_indices[i]));
01610                   v_norm = std::max(v_norm, norm_weight * hess_u_h.size());
01611                 }
01612 #endif
01613             }
01614         }
01615     }
01616 
01617   if (using_hilbert_norm)
01618     {
01619       CommWorld.sum(v_norm);
01620       v_norm = std::sqrt(v_norm);
01621     }
01622   else
01623     {
01624       CommWorld.max(v_norm);
01625     }
01626 
01627   STOP_LOG ("calculate_norm()", "System");
01628 
01629   return v_norm;
01630 }

Real libMesh::System::calculate_norm ( const NumericVector< Number > &  v,
unsigned int  var = 0,
FEMNormType  norm_type = L2 
) const [inherited]
Returns:
a norm of variable var in the vector v, in the specified norm (e.g. L2, L_INF, H1)

Definition at line 1379 of file system.C.

References libMeshEnums::DISCRETE_L1, libMeshEnums::DISCRETE_L2, libMeshEnums::DISCRETE_L_INF, libMesh::System::discrete_var_norm(), libMeshEnums::L2, libMesh::System::n_vars(), and libMesh::Real.

Referenced by libMesh::AdaptiveTimeSolver::calculate_norm(), and libMesh::UnsteadySolver::du().

01382 {
01383   //short circuit to save time
01384   if(norm_type == DISCRETE_L1 ||
01385      norm_type == DISCRETE_L2 ||
01386      norm_type == DISCRETE_L_INF)
01387     return discrete_var_norm(v,var,norm_type);
01388 
01389   // Not a discrete norm
01390   std::vector<FEMNormType> norms(this->n_vars(), L2);
01391   std::vector<Real> weights(this->n_vars(), 0.0);
01392   norms[var] = norm_type;
01393   weights[var] = 1.0;
01394   Real val = this->calculate_norm(v, SystemNorm(norms, weights));
01395   return val;
01396 }

virtual void libMesh::RBEIMConstruction::clear (  )  [virtual]

Clear this object.

Reimplemented from libMesh::RBConstruction.

bool libMesh::System::compare ( const System &  other_system,
const Real  threshold,
const bool  verbose 
) const [virtual, inherited]
Returns:
true when the other system contains identical data, up to the given threshold. Outputs some diagnostic info when verbose is set.

Definition at line 525 of file system.C.

References libMesh::System::_can_add_vectors, libMesh::System::_sys_name, libMesh::System::_vectors, libMesh::AutoPtr< Tp >::get(), libMesh::System::get_vector(), libMesh::System::n_vectors(), libMesh::System::name(), libMesh::out, and libMesh::System::solution.

Referenced by libMesh::EquationSystems::compare().

00528 {
00529   // we do not care for matrices, but for vectors
00530   libmesh_assert (!_can_add_vectors);
00531   libmesh_assert (!other_system._can_add_vectors);
00532 
00533   if (verbose)
00534     {
00535       libMesh::out << "  Systems \"" << _sys_name << "\"" << std::endl;
00536       libMesh::out << "   comparing matrices not supported." << std::endl;
00537       libMesh::out << "   comparing names...";
00538     }
00539 
00540   // compare the name: 0 means identical
00541   const int name_result = _sys_name.compare(other_system.name());
00542   if (verbose)
00543     {
00544       if (name_result == 0)
00545         libMesh::out << " identical." << std::endl;
00546       else
00547         libMesh::out << "  names not identical." << std::endl;
00548       libMesh::out << "   comparing solution vector...";
00549     }
00550 
00551 
00552   // compare the solution: -1 means identical
00553   const int solu_result = solution->compare (*other_system.solution.get(),
00554                                              threshold);
00555 
00556   if (verbose)
00557     {
00558       if (solu_result == -1)
00559         libMesh::out << " identical up to threshold." << std::endl;
00560       else
00561         libMesh::out << "  first difference occured at index = "
00562                   << solu_result << "." << std::endl;
00563     }
00564 
00565 
00566   // safety check, whether we handle at least the same number
00567   // of vectors
00568   std::vector<int> ov_result;
00569 
00570   if (this->n_vectors() != other_system.n_vectors())
00571     {
00572       if (verbose)
00573         {
00574           libMesh::out << "   Fatal difference. This system handles "
00575                     << this->n_vectors() << " add'l vectors," << std::endl
00576                     << "   while the other system handles "
00577                     << other_system.n_vectors()
00578                     << " add'l vectors." << std::endl
00579                     << "   Aborting comparison." << std::endl;
00580         }
00581       return false;
00582     }
00583   else if (this->n_vectors() == 0)
00584     {
00585       // there are no additional vectors...
00586       ov_result.clear ();
00587     }
00588   else
00589     {
00590       // compare other vectors
00591       for (const_vectors_iterator pos = _vectors.begin();
00592            pos != _vectors.end(); ++pos)
00593         {
00594           if (verbose)
00595               libMesh::out << "   comparing vector \""
00596                         << pos->first << "\" ...";
00597 
00598           // assume they have the same name
00599           const NumericVector<Number>& other_system_vector =
00600               other_system.get_vector(pos->first);
00601 
00602           ov_result.push_back(pos->second->compare (other_system_vector,
00603                                                     threshold));
00604 
00605           if (verbose)
00606             {
00607               if (ov_result[ov_result.size()-1] == -1)
00608                 libMesh::out << " identical up to threshold." << std::endl;
00609               else
00610                 libMesh::out << " first difference occured at" << std::endl
00611                           << "   index = " << ov_result[ov_result.size()-1] << "." << std::endl;
00612             }
00613 
00614         }
00615 
00616     } // finished comparing additional vectors
00617 
00618 
00619   bool overall_result;
00620 
00621   // sum up the results
00622   if ((name_result==0) && (solu_result==-1))
00623     {
00624       if (ov_result.size()==0)
00625         overall_result = true;
00626       else
00627         {
00628           bool ov_identical;
00629           unsigned int n    = 0;
00630           do
00631             {
00632               ov_identical = (ov_result[n]==-1);
00633               n++;
00634             }
00635           while (ov_identical && n<ov_result.size());
00636           overall_result = ov_identical;
00637         }
00638     }
00639   else
00640     overall_result = false;
00641 
00642   if (verbose)
00643     {
00644       libMesh::out << "   finished comparisons, ";
00645       if (overall_result)
00646         libMesh::out << "found no differences." << std::endl << std::endl;
00647       else
00648         libMesh::out << "found differences." << std::endl << std::endl;
00649     }
00650 
00651   return overall_result;
00652 }

virtual Real libMesh::RBEIMConstruction::compute_best_fit_error (  )  [virtual]

We compute the best fit of parametrized_function into the EIM space and then evaluate the error in the norm defined by inner_product_matrix.

Returns:
the error in the best fit
virtual void libMesh::RBConstruction::compute_Fq_representor_innerprods ( bool  compute_inner_products = true  )  [protected, virtual, inherited]

Compute the terms that are combined `online' to determine the dual norm of the residual. Here we compute the terms associated with the right-hand side. These terms are basis independent, hence we separate them from the rest of the calculations that are done in update_residual_terms. By default, inner product terms are also computed, but you can turn this feature off e.g. if you are already reading in that data from files.

virtual Real libMesh::RBConstruction::compute_max_error_bound (  )  [virtual, inherited]

(i) Compute the a posteriori error bound for each set of parameters in the training set, (ii) set current_parameters to the parameters that maximize the error bound, and (iii) return the maximum error bound.

virtual void libMesh::RBConstruction::compute_output_dual_innerprods (  )  [protected, virtual, inherited]

Compute and store the dual norm of each output functional.

Number libMesh::System::current_solution ( const dof_id_type  global_dof_number  )  const [inherited]
Returns:
the current solution for the specified global DOF.

Definition at line 181 of file system.C.

References libMesh::System::_dof_map, and libMesh::System::current_local_solution.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::HPCoarsenTest::add_projection(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::FEMContext::pre_fe_reinit(), libMesh::HPCoarsenTest::select_refinement(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

00182 {
00183   // Check the sizes
00184   libmesh_assert_less (global_dof_number, _dof_map->n_dofs());
00185   libmesh_assert_less (global_dof_number, current_local_solution->size());
00186 
00187   return (*current_local_solution)(global_dof_number);
00188 }

void libMesh::System::deactivate (  )  [inline, inherited]

Deactivates the system. Only active systems are solved.

Definition at line 1882 of file system.h.

References libMesh::System::_active.

01883 {
01884   _active = false;
01885 }

void libMesh::LinearImplicitSystem::detach_shell_matrix ( void   )  [inline, inherited]

Detaches a shell matrix. Same as attach_shell_matrix(NULL).

Definition at line 176 of file linear_implicit_system.h.

References libMesh::LinearImplicitSystem::attach_shell_matrix().

00176 { attach_shell_matrix(NULL); }

void libMesh::ImplicitSystem::disable_cache (  )  [virtual, inherited]

Assembles & solves the linear system Ax=b. Avoids use of any cached data that might affect any solve result. Should be overloaded in derived systems.

Reimplemented from libMesh::System.

Definition at line 313 of file implicit_system.C.

References libMesh::System::assemble_before_solve, libMesh::ImplicitSystem::get_linear_solver(), and libMesh::LinearSolver< T >::reuse_preconditioner().

00313                                     {
00314   this->assemble_before_solve = true;
00315   this->get_linear_solver()->reuse_preconditioner(false);
00316 }

void libMesh::ReferenceCounter::disable_print_counter_info (  )  [static, inherited]

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

00107 {
00108   _enable_print_counter = false;
00109   return;
00110 }

void libMesh::ReferenceCounter::disable_print_counter_info (  )  [static, inherited]

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

00107 {
00108   _enable_print_counter = false;
00109   return;
00110 }

void libMesh::ReferenceCounter::enable_print_counter_info (  )  [static, inherited]

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

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

00101 {
00102   _enable_print_counter = true;
00103   return;
00104 }

void libMesh::ReferenceCounter::enable_print_counter_info (  )  [static, inherited]

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

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

00101 {
00102   _enable_print_counter = true;
00103   return;
00104 }

virtual void libMesh::RBEIMConstruction::enrich_RB_space (  )  [protected, virtual]

Add a new basis function to the RB space. Overload to enrich with the EIM basis functions.

Reimplemented from libMesh::RBConstruction.

Number libMesh::RBEIMConstruction::evaluate_mesh_function ( unsigned int  var_number,
Point  p 
)

Evaluate the mesh function at the specified point and for the specified variable.

Real libMesh::LinearImplicitSystem::final_linear_residual (  )  const [inline, inherited]

Returns the final residual for the linear system solve.

Definition at line 160 of file linear_implicit_system.h.

References libMesh::LinearImplicitSystem::_final_linear_residual.

00160 { return _final_linear_residual; }

void libMesh::ImplicitSystem::forward_qoi_parameter_sensitivity ( const QoISet &  qoi_indices,
const ParameterVector &  parameters,
SensitivityData &  sensitivities 
) [virtual, inherited]

Solves for the derivative of each of the system's quantities of interest q in qoi[qoi_indices] with respect to each parameter in parameters, placing the result for qoi i and parameter j into sensitivities[i][j].

Uses the forward sensitivity method.

Currently uses finite differenced derivatives (partial q / partial p) and (partial R / partial p).

Reimplemented from libMesh::System.

Definition at line 795 of file implicit_system.C.

References libMesh::SensitivityData::allocate_data(), libMesh::QoISet::has_index(), libMesh::Real, libMesh::ParameterVector::size(), and libMesh::TOLERANCE.

00798 {
00799   const unsigned int Np = libmesh_cast_int<unsigned int>
00800     (parameters.size());
00801   const unsigned int Nq = libmesh_cast_int<unsigned int>
00802     (qoi.size());
00803 
00804   // We currently get partial derivatives via central differencing
00805   const Real delta_p = TOLERANCE;
00806 
00807   // An introduction to the problem:
00808   //
00809   // Residual R(u(p),p) = 0
00810   // partial R / partial u = J = system matrix
00811   //
00812   // This implies that:
00813   // d/dp(R) = 0
00814   // (partial R / partial p) +
00815   // (partial R / partial u) * (partial u / partial p) = 0
00816 
00817   // We first solve for (partial u / partial p) for each parameter:
00818   // J * (partial u / partial p) = - (partial R / partial p)
00819 
00820   this->sensitivity_solve(parameters);
00821 
00822   // Get ready to fill in senstivities:
00823   sensitivities.allocate_data(qoi_indices, *this, parameters);
00824 
00825   // We use the identity:
00826   // dq/dp = (partial q / partial p) + (partial q / partial u) *
00827   //         (partial u / partial p)
00828 
00829   // We get (partial q / partial u) from the user
00830   this->assemble_qoi_derivative(qoi_indices);
00831 
00832   // We don't need these to be closed() in this function, but libMesh
00833   // standard practice is to have them closed() by the time the
00834   // function exits
00835   for (unsigned int i=0; i != this->qoi.size(); ++i)
00836     if (qoi_indices.has_index(i))
00837       this->get_adjoint_rhs(i).close();
00838 
00839   for (unsigned int j=0; j != Np; ++j)
00840     {
00841       // (partial q / partial p) ~= (q(p+dp)-q(p-dp))/(2*dp)
00842 
00843       Number old_parameter = *parameters[j];
00844 
00845       *parameters[j] = old_parameter - delta_p;
00846       this->assemble_qoi();
00847       std::vector<Number> qoi_minus = this->qoi;
00848 
00849       *parameters[j] = old_parameter + delta_p;
00850       this->assemble_qoi();
00851       std::vector<Number>& qoi_plus = this->qoi;
00852 
00853       std::vector<Number> partialq_partialp(Nq, 0);
00854       for (unsigned int i=0; i != Nq; ++i)
00855         if (qoi_indices.has_index(i))
00856           partialq_partialp[i] = (qoi_plus[i] - qoi_minus[i]) / (2.*delta_p);
00857 
00858       // Don't leave the parameter changed
00859       *parameters[j] = old_parameter;
00860 
00861       for (unsigned int i=0; i != Nq; ++i)
00862         if (qoi_indices.has_index(i))
00863           sensitivities[i][j] = partialq_partialp[i] +
00864             this->get_adjoint_rhs(i).dot(this->get_sensitivity_solution(j));
00865     }
00866 
00867   // All parameters have been reset.
00868   // We didn't cache the original rhs or matrix for memory reasons,
00869   // but we can restore them to a state consistent solution -
00870   // principle of least surprise.
00871   this->assembly(true, true);
00872   this->rhs->close();
00873   this->matrix->close();
00874   this->assemble_qoi(qoi_indices);
00875 }

static void libMesh::RBConstructionBase< LinearImplicitSystem >::generate_training_parameters_deterministic ( std::map< std::string, bool >  log_param_scale,
std::map< std::string, NumericVector< Number > * > &  training_parameters_in,
unsigned int  n_training_samples_in,
const RBParameters &  min_parameters,
const RBParameters &  max_parameters,
bool  serial_training_set = false 
) [static, protected, inherited]

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.

static void libMesh::RBConstructionBase< LinearImplicitSystem >::generate_training_parameters_partially_random ( const std::string &  deterministic_parameter_name,
const unsigned int  deterministic_parameter_repeats,
std::map< std::string, bool >  log_param_scale,
std::map< std::string, NumericVector< Number > * > &  training_parameters_in,
unsigned int  n_deterministic_training_samples_in,
const RBParameters &  min_parameters,
const RBParameters &  max_parameters,
int  training_parameters_random_seed = -1,
bool  serial_training_set = false 
) [static, protected, inherited]

Static helper function for generating a "partially" random set of parameters, that is the parameter indicated by this->get_deterministic_training_parameter() will be deterministic.

static void libMesh::RBConstructionBase< LinearImplicitSystem >::generate_training_parameters_random ( std::map< std::string, bool >  log_param_scale,
std::map< std::string, NumericVector< Number > * > &  training_parameters_in,
unsigned int  n_training_samples_in,
const RBParameters &  min_parameters,
const RBParameters &  max_parameters,
int  training_parameters_random_seed = -1,
bool  serial_training_set = false 
) [static, protected, inherited]

Static helper function for generating a randomized set of parameters.

const NumericVector< Number > & libMesh::System::get_adjoint_rhs ( unsigned int  i = 0  )  const [inherited]
Returns:
a reference to one of the system's adjoint rhs vectors, by default the one corresponding to the first qoi.

Definition at line 1036 of file system.C.

References libMesh::System::get_vector().

01037 {
01038   std::ostringstream adjoint_rhs_name;
01039   adjoint_rhs_name << "adjoint_rhs" << i;
01040 
01041   return this->get_vector(adjoint_rhs_name.str());
01042 }

NumericVector< Number > & libMesh::System::get_adjoint_rhs ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's adjoint rhs vectors, by default the one corresponding to the first qoi. This what the user's QoI derivative code should assemble when setting up an adjoint problem

Definition at line 1026 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::adjoint_solve(), and libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

01027 {
01028   std::ostringstream adjoint_rhs_name;
01029   adjoint_rhs_name << "adjoint_rhs" << i;
01030 
01031   return this->get_vector(adjoint_rhs_name.str());
01032 }

const NumericVector< Number > & libMesh::System::get_adjoint_solution ( unsigned int  i = 0  )  const [inherited]
Returns:
a reference to one of the system's adjoint solution vectors, by default the one corresponding to the first qoi.

Definition at line 976 of file system.C.

References libMesh::System::get_vector().

00977 {
00978   std::ostringstream adjoint_name;
00979   adjoint_name << "adjoint_solution" << i;
00980 
00981   return this->get_vector(adjoint_name.str());
00982 }

NumericVector< Number > & libMesh::System::get_adjoint_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's adjoint solution vectors, by default the one corresponding to the first qoi.

Definition at line 966 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::adjoint_solve(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), and libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

00967 {
00968   std::ostringstream adjoint_name;
00969   adjoint_name << "adjoint_solution" << i;
00970 
00971   return this->get_vector(adjoint_name.str());
00972 }

void libMesh::System::get_all_variable_numbers ( std::vector< unsigned int > &  all_variable_numbers  )  const [inherited]

Fills all_variable_numbers with all the variable numbers for the variables that have been added to this system.

Definition at line 1259 of file system.C.

References libMesh::System::_variable_numbers, and libMesh::System::n_vars().

01260 {
01261   all_variable_numbers.resize(n_vars());
01262 
01263   // Make sure the variable exists
01264   std::map<std::string, unsigned short int>::const_iterator
01265     it = _variable_numbers.begin();
01266   std::map<std::string, unsigned short int>::const_iterator
01267     it_end = _variable_numbers.end();
01268 
01269   unsigned int count = 0;
01270   for( ; it != it_end; ++it)
01271   {
01272     all_variable_numbers[count] = it->second;
01273     count++;
01274   }
01275 }

SparseMatrix<Number>* libMesh::RBConstruction::get_Aq ( unsigned int  q  )  [inherited]

Get a pointer to Aq.

ElemAssembly& libMesh::RBConstruction::get_constraint_assembly (  )  [inherited]
Returns:
a reference to the constraint assembly object
unsigned int libMesh::RBConstruction::get_delta_N (  )  const [inline, inherited]

Get delta_N, the number of basis functions we add to the RB space per iteration of the greedy algorithm. For steady-state systems, this should be 1, but can be more than 1 for time-dependent systems.

Definition at line 344 of file rb_construction.h.

References libMesh::RBConstruction::delta_N.

00344 { return delta_N; }

const std::string& libMesh::RBConstructionBase< LinearImplicitSystem >::get_deterministic_training_parameter_name (  )  const [inherited]

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

unsigned int libMesh::RBConstructionBase< LinearImplicitSystem >::get_deterministic_training_parameter_repeats (  )  const [inherited]

Get the number of times each sample of the deterministic training parameter is repeated.

DofMap & libMesh::System::get_dof_map (  )  [inline, inherited]
Returns:
a writeable reference to this system's _dof_map.

Definition at line 1858 of file system.h.

References libMesh::System::_dof_map.

01859 {
01860   return *_dof_map;
01861 }

const DofMap & libMesh::System::get_dof_map (  )  const [inline, inherited]
Returns:
a constant reference to this system's _dof_map.

Definition at line 1850 of file system.h.

References libMesh::System::_dof_map.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::ExactSolution::_compute_error(), libMesh::HPCoarsenTest::add_projection(), libMesh::UnsteadySolver::adjoint_advance_timestep(), libMesh::ImplicitSystem::adjoint_solve(), libMesh::UnsteadySolver::advance_timestep(), libMesh::EquationSystems::allgather(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_solution_vector(), libMesh::System::calculate_norm(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), DMCreateDomainDecomposition_libMesh(), DMCreateFieldDecomposition_libMesh(), DMFunction_libMesh(), DMJacobian_libMesh(), DMLibMeshSetSystem(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_info(), libMesh::EquationSystems::get_solution(), libMesh::SystemSubsetBySubdomain::init(), libMesh::UnsteadySolver::init_data(), libMesh::ImplicitSystem::init_matrices(), libMesh::EigenSystem::init_matrices(), libMesh::CondensedEigenSystem::initialize_condensed_dofs(), libMesh::System::local_dof_indices(), libMesh::DofMap::max_constraint_error(), libMesh::UnsteadySolver::old_nonlinear_solution(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), libMesh::ProjectSolution::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::ErrorVector::plot_error(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::FEMContext::pre_fe_reinit(), libMesh::System::project_vector(), libMesh::System::re_update(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::UnsteadySolver::reinit(), libMesh::ImplicitSystem::reinit(), libMesh::EquationSystems::reinit(), libMesh::EigenSystem::reinit(), libMesh::UnsteadySolver::retrieve_timestep(), libMesh::HPCoarsenTest::select_refinement(), libMesh::ImplicitSystem::sensitivity_solve(), libMesh::PetscDiffSolver::solve(), libMesh::NewtonSolver::solve(), libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve(), libMesh::ImplicitSystem::weighted_sensitivity_solve(), libMesh::System::write_parallel_data(), libMesh::EnsightIO::write_scalar_ascii(), libMesh::System::write_SCALAR_dofs(), and libMesh::EnsightIO::write_vector_ascii().

01851 {
01852   return *_dof_map;
01853 }

std::vector<ElemAssembly*> libMesh::RBEIMConstruction::get_eim_assembly_objects (  ) 
Returns:
the vector of assembly objects that point to this RBEIMConstruction.
EquationSystems& libMesh::System::get_equation_systems (  )  [inline, inherited]
Returns:
a reference to this system's parent EquationSystems object.

Definition at line 684 of file system.h.

References libMesh::System::_equation_systems.

00684 { return _equation_systems; }

const EquationSystems& libMesh::System::get_equation_systems (  )  const [inline, inherited]
Returns:
a constant reference to this system's parent EquationSystems object.

Definition at line 679 of file system.h.

References libMesh::System::_equation_systems.

Referenced by libMesh::NewmarkSystem::clear(), libMesh::FrequencySystem::clear_all(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::ImplicitSystem::get_linear_solve_parameters(), libMesh::FrequencySystem::init_data(), libMesh::FrequencySystem::n_frequencies(), libMesh::System::project_vector(), libMesh::FrequencySystem::set_current_frequency(), libMesh::FrequencySystem::set_frequencies(), libMesh::FrequencySystem::set_frequencies_by_range(), libMesh::FrequencySystem::set_frequencies_by_steps(), libMesh::NewmarkSystem::set_newmark_parameters(), libMesh::NonlinearImplicitSystem::set_solver_parameters(), libMesh::LinearImplicitSystem::solve(), libMesh::FrequencySystem::solve(), libMesh::EigenSystem::solve(), libMesh::CondensedEigenSystem::solve(), and libMesh::WrappedFunction< Output >::WrappedFunction().

00679 { return _equation_systems; }

numeric_index_type libMesh::RBConstructionBase< LinearImplicitSystem >::get_first_local_training_index (  )  const [inherited]

Get the first local index of the training parameters.

NumericVector<Number>* libMesh::RBConstruction::get_Fq ( unsigned int  q  )  [inherited]

Get a pointer to Fq.

static void libMesh::RBConstructionBase< LinearImplicitSystem >::get_global_max_error_pair ( std::pair< unsigned int, Real > &  error_pair  )  [static, protected, inherited]

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

const RBParameters& libMesh::RBConstruction::get_greedy_parameter ( unsigned int  i  )  [inherited]

Return the parameters chosen during the i^th step of the Greedy algorithm.

std::string libMesh::ReferenceCounter::get_info (  )  [static, inherited]

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().

00048 {
00049 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
00050 
00051   std::ostringstream oss;
00052 
00053   oss << '\n'
00054       << " ---------------------------------------------------------------------------- \n"
00055       << "| Reference count information                                                |\n"
00056       << " ---------------------------------------------------------------------------- \n";
00057 
00058   for (Counts::iterator it = _counts.begin();
00059        it != _counts.end(); ++it)
00060     {
00061       const std::string name(it->first);
00062       const unsigned int creations    = it->second.first;
00063       const unsigned int destructions = it->second.second;
00064 
00065       oss << "| " << name << " reference count information:\n"
00066           << "|  Creations:    " << creations    << '\n'
00067           << "|  Destructions: " << destructions << '\n';
00068     }
00069 
00070   oss << " ---------------------------------------------------------------------------- \n";
00071 
00072   return oss.str();
00073 
00074 #else
00075 
00076   return "";
00077 
00078 #endif
00079 }

std::string libMesh::ReferenceCounter::get_info (  )  [static, inherited]

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().

00048 {
00049 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
00050 
00051   std::ostringstream oss;
00052 
00053   oss << '\n'
00054       << " ---------------------------------------------------------------------------- \n"
00055       << "| Reference count information                                                |\n"
00056       << " ---------------------------------------------------------------------------- \n";
00057 
00058   for (Counts::iterator it = _counts.begin();
00059        it != _counts.end(); ++it)
00060     {
00061       const std::string name(it->first);
00062       const unsigned int creations    = it->second.first;
00063       const unsigned int destructions = it->second.second;
00064 
00065       oss << "| " << name << " reference count information:\n"
00066           << "|  Creations:    " << creations    << '\n'
00067           << "|  Destructions: " << destructions << '\n';
00068     }
00069 
00070   oss << " ---------------------------------------------------------------------------- \n";
00071 
00072   return oss.str();
00073 
00074 #else
00075 
00076   return "";
00077 
00078 #endif
00079 }

std::string libMesh::System::get_info (  )  const [inherited]
Returns:
a string containing information about the system.

Definition at line 1634 of file system.C.

References libMesh::Utility::enum_to_string< FEFamily >(), libMesh::Utility::enum_to_string< InfMapType >(), libMesh::Utility::enum_to_string< Order >(), libMesh::FEType::family, libMesh::System::get_dof_map(), libMesh::DofMap::get_info(), libMesh::FEType::inf_map, libMesh::System::n_constrained_dofs(), libMesh::System::n_dofs(), libMesh::System::n_local_constrained_dofs(), libMesh::System::n_local_dofs(), libMesh::System::n_matrices(), libMesh::System::n_variable_groups(), libMesh::VariableGroup::n_variables(), libMesh::System::n_vectors(), libMesh::VariableGroup::name(), libMesh::System::name(), libMesh::System::number(), libMesh::FEType::order, libMesh::FEType::radial_family, libMesh::FEType::radial_order, libMesh::System::system_type(), libMesh::Variable::type(), libMesh::DofMap::variable_group(), and libMesh::System::variable_group().

01635 {
01636   std::ostringstream oss;
01637 
01638 
01639   const std::string& sys_name = this->name();
01640 
01641   oss << "   System #"  << this->number() << ", \"" << sys_name << "\"\n"
01642       << "    Type \""  << this->system_type() << "\"\n"
01643       << "    Variables=";
01644 
01645   for (unsigned int vg=0; vg<this->n_variable_groups(); vg++)
01646     {
01647       const VariableGroup &vg_description (this->variable_group(vg));
01648 
01649       if (vg_description.n_variables() > 1) oss << "{ ";
01650       for (unsigned int vn=0; vn<vg_description.n_variables(); vn++)
01651         oss << "\"" << vg_description.name(vn) << "\" ";
01652       if (vg_description.n_variables() > 1) oss << "} ";
01653     }
01654 
01655   oss << '\n';
01656 
01657   oss << "    Finite Element Types=";
01658 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
01659   for (unsigned int vg=0; vg<this->n_variable_groups(); vg++)
01660     oss << "\""
01661         << Utility::enum_to_string<FEFamily>(this->get_dof_map().variable_group(vg).type().family)
01662         << "\" ";
01663 #else
01664   for (unsigned int vg=0; vg<this->n_variable_groups(); vg++)
01665     {
01666       oss << "\""
01667           << Utility::enum_to_string<FEFamily>(this->get_dof_map().variable_group(vg).type().family)
01668           << "\", \""
01669           << Utility::enum_to_string<FEFamily>(this->get_dof_map().variable_group(vg).type().radial_family)
01670           << "\" ";
01671     }
01672 
01673   oss << '\n' << "    Infinite Element Mapping=";
01674   for (unsigned int vg=0; vg<this->n_variable_groups(); vg++)
01675     oss << "\""
01676         << Utility::enum_to_string<InfMapType>(this->get_dof_map().variable_group(vg).type().inf_map)
01677         << "\" ";
01678 #endif
01679 
01680   oss << '\n';
01681 
01682   oss << "    Approximation Orders=";
01683   for (unsigned int vg=0; vg<this->n_variable_groups(); vg++)
01684     {
01685 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
01686       oss << "\""
01687           << Utility::enum_to_string<Order>(this->get_dof_map().variable_group(vg).type().order)
01688           << "\" ";
01689 #else
01690       oss << "\""
01691           << Utility::enum_to_string<Order>(this->get_dof_map().variable_group(vg).type().order)
01692           << "\", \""
01693           << Utility::enum_to_string<Order>(this->get_dof_map().variable_group(vg).type().radial_order)
01694           << "\" ";
01695 #endif
01696     }
01697 
01698   oss << '\n';
01699 
01700   oss << "    n_dofs()="             << this->n_dofs()             << '\n';
01701   oss << "    n_local_dofs()="       << this->n_local_dofs()       << '\n';
01702 #ifdef LIBMESH_ENABLE_CONSTRAINTS
01703   oss << "    n_constrained_dofs()=" << this->n_constrained_dofs() << '\n';
01704   oss << "    n_local_constrained_dofs()=" << this->n_local_constrained_dofs() << '\n';
01705 #endif
01706 
01707   oss << "    " << "n_vectors()="  << this->n_vectors()  << '\n';
01708   oss << "    " << "n_matrices()="  << this->n_matrices()  << '\n';
01709 //   oss << "    " << "n_additional_matrices()=" << this->n_additional_matrices() << '\n';
01710 
01711   oss << this->get_dof_map().get_info();
01712 
01713   return oss.str();
01714 }

ElemAssembly& libMesh::RBConstruction::get_inner_product_assembly (  )  [inherited]
Returns:
a reference to the inner product assembly object
SparseMatrix<Number>* libMesh::RBConstruction::get_inner_product_matrix (  )  [inherited]

Get a pointer to inner_product_matrix. Accessing via this function, rather than directly through the class member allows us to do error checking (e.g. inner_product_matrix is not defined in low-memory mode).

numeric_index_type libMesh::RBConstructionBase< LinearImplicitSystem >::get_last_local_training_index (  )  const [inherited]

Get the last local index of the training parameters.

std::pair< unsigned int, Real > libMesh::ImplicitSystem::get_linear_solve_parameters (  )  const [virtual, inherited]

Returns an integer corresponding to the upper iteration count limit and a Real corresponding to the convergence tolerance to be used in linear adjoint and/or sensitivity solves

Reimplemented in libMesh::DifferentiableSystem, and libMesh::NonlinearImplicitSystem.

Definition at line 1359 of file implicit_system.C.

References libMesh::System::get_equation_systems(), and libMesh::Real.

Referenced by libMesh::ImplicitSystem::adjoint_solve(), libMesh::ImplicitSystem::sensitivity_solve(), libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve(), and libMesh::ImplicitSystem::weighted_sensitivity_solve().

01360 {
01361   return std::make_pair(this->get_equation_systems().parameters.get<unsigned int>("linear solver maximum iterations"),
01362                         this->get_equation_systems().parameters.get<Real>("linear solver tolerance"));
01363 }

LinearSolver< Number > * libMesh::LinearImplicitSystem::get_linear_solver (  )  const [virtual, inherited]

Returns a pointer to a linear solver appropriate for use in adjoint and/or sensitivity solves

Reimplemented from libMesh::ImplicitSystem.

Definition at line 360 of file linear_implicit_system.C.

References libMesh::AutoPtr< Tp >::get(), and libMesh::LinearImplicitSystem::linear_solver.

00361 {
00362   return linear_solver.get();
00363 }

numeric_index_type libMesh::RBConstructionBase< LinearImplicitSystem >::get_local_n_training_samples (  )  const [inherited]

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

SparseMatrix< Number > & libMesh::ImplicitSystem::get_matrix ( const std::string &  mat_name  )  [inherited]
Returns:
a writeable reference to this system's additional matrix named mat_name. None of these matrices is involved in the solution process. Access is only granted when the matrix is already properly initialized.

Definition at line 277 of file implicit_system.C.

References libMesh::ImplicitSystem::_matrices, and libMesh::err.

00278 {
00279   // Make sure the matrix exists
00280   matrices_iterator pos = _matrices.find (mat_name);
00281 
00282   if (pos == _matrices.end())
00283     {
00284       libMesh::err << "ERROR: matrix "
00285                     << mat_name
00286                     << " does not exist in this system!"
00287                     << std::endl;
00288       libmesh_error();
00289     }
00290 
00291   return *(pos->second);
00292 }

const SparseMatrix< Number > & libMesh::ImplicitSystem::get_matrix ( const std::string &  mat_name  )  const [inherited]
Returns:
a const reference to this system's additional matrix named mat_name. None of these matrices is involved in the solution process. Access is only granted when the matrix is already properly initialized.

Definition at line 258 of file implicit_system.C.

References libMesh::ImplicitSystem::_matrices, and libMesh::err.

Referenced by libMesh::NewmarkSystem::compute_matrix(), libMesh::EigenTimeSolver::solve(), and libMesh::NewmarkSystem::update_rhs().

00259 {
00260   // Make sure the matrix exists
00261   const_matrices_iterator pos = _matrices.find (mat_name);
00262 
00263   if (pos == _matrices.end())
00264     {
00265       libMesh::err << "ERROR: matrix "
00266                     << mat_name
00267                     << " does not exist in this system!"
00268                     << std::endl;
00269       libmesh_error();
00270     }
00271 
00272   return *(pos->second);
00273 }

MeshBase & libMesh::System::get_mesh (  )  [inline, inherited]
Returns:
a reference to this systems's _mesh.

Definition at line 1842 of file system.h.

References libMesh::System::_mesh.

01843 {
01844   return _mesh;
01845 }

const MeshBase & libMesh::System::get_mesh (  )  const [inline, inherited]
Returns:
a constant reference to this systems's _mesh.

Definition at line 1834 of file system.h.

References libMesh::System::_mesh.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::HPCoarsenTest::add_projection(), libMesh::FEMSystem::assemble_qoi(), libMesh::FEMSystem::assemble_qoi_derivative(), libMesh::FEMSystem::assembly(), libMesh::System::calculate_norm(), DMCreateDomainDecomposition_libMesh(), DMCreateFieldDecomposition_libMesh(), DMLibMeshSetSystem(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::SystemSubsetBySubdomain::init(), libMesh::System::init_data(), libMesh::ImplicitSystem::init_matrices(), libMesh::EigenSystem::init_matrices(), libMesh::System::local_dof_indices(), libMesh::DofMap::max_constraint_error(), libMesh::FEMSystem::mesh_position_get(), libMesh::FEMSystem::mesh_position_set(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::ProjectSolution::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::FEMSystem::postprocess(), libMesh::System::project_vector(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::System::read_parallel_data(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::ImplicitSystem::reinit(), libMesh::EigenSystem::reinit(), libMesh::HPSingularity::select_refinement(), libMesh::HPCoarsenTest::select_refinement(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), and libMesh::System::zero_variable().

01835 {
01836   return _mesh;
01837 }

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

Get the number of parameters.

numeric_index_type libMesh::RBConstructionBase< LinearImplicitSystem >::get_n_training_samples (  )  const [inherited]

Get the total number of training samples.

unsigned int libMesh::RBConstruction::get_Nmax (  )  const [inline, inherited]

Get/set Nmax, the maximum number of RB functions we are willing to compute.

Definition at line 182 of file rb_construction.h.

References libMesh::RBConstruction::Nmax.

00182 { return Nmax; }

SparseMatrix<Number>* libMesh::RBConstruction::get_non_dirichlet_Aq ( unsigned int  q  )  [inherited]

Get a pointer to non_dirichlet_Aq.

NumericVector<Number>* libMesh::RBConstruction::get_non_dirichlet_Fq ( unsigned int  q  )  [inherited]

Get a pointer to non-Dirichlet Fq.

SparseMatrix<Number>* libMesh::RBConstruction::get_non_dirichlet_inner_product_matrix (  )  [inherited]

Get a pointer to non_dirichlet_inner_product_matrix. Accessing via this function, rather than directly through the class member allows us to do error checking (e.g. non_dirichlet_inner_product_matrix is not defined in low-memory mode, and we need store_non_dirichlet_operators==true).

NumericVector<Number>* libMesh::RBConstruction::get_non_dirichlet_output_vector ( unsigned int  n,
unsigned int  q_l 
) [inherited]

Get a pointer to non-Dirichlet output vector.

NumericVector<Number>* libMesh::RBConstruction::get_output_vector ( unsigned int  n,
unsigned int  q_l 
) [inherited]

Get a pointer to the n^th output.

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

Get maximum allowable value of parameter param_name.

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

Get minimum allowable value of parameter param_name.

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

Get the current parameters.

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

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

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

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

RBParameters libMesh::RBConstructionBase< LinearImplicitSystem >::get_params_from_training_set ( unsigned int  index  )  [protected, inherited]

Return the RBParameters in index index of training set.

RBAssemblyExpansion& libMesh::RBConstruction::get_rb_assembly_expansion (  )  [inherited]
Returns:
a reference to the rb_assembly_expansion object
virtual Real libMesh::RBEIMConstruction::get_RB_error_bound (  )  [protected, virtual]

Overload to return the best fit error. This function is used in the Greedy algorithm to select the next parameter.

Reimplemented from libMesh::RBConstruction.

RBEvaluation& libMesh::RBConstruction::get_rb_evaluation (  )  [inherited]

Get a reference to the RBEvaluation object.

RBThetaExpansion& libMesh::RBConstruction::get_rb_theta_expansion (  )  [inherited]

Get a reference to the RBThetaExpansion object that that belongs to rb_eval.

const NumericVector< Number > & libMesh::System::get_sensitivity_rhs ( unsigned int  i = 0  )  const [inherited]
Returns:
a reference to one of the system's sensitivity rhs vectors, by default the one corresponding to the first parameter.

Definition at line 1066 of file system.C.

References libMesh::System::get_vector().

01067 {
01068   std::ostringstream sensitivity_rhs_name;
01069   sensitivity_rhs_name << "sensitivity_rhs" << i;
01070 
01071   return this->get_vector(sensitivity_rhs_name.str());
01072 }

NumericVector< Number > & libMesh::System::get_sensitivity_rhs ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's sensitivity rhs vectors, by default the one corresponding to the first parameter. By default these vectors are built by the library, using finite differences, when assemble_residual_derivatives() is called.

When assembled, this vector should hold -(partial R / partial p_i)

Definition at line 1056 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::sensitivity_solve().

01057 {
01058   std::ostringstream sensitivity_rhs_name;
01059   sensitivity_rhs_name << "sensitivity_rhs" << i;
01060 
01061   return this->get_vector(sensitivity_rhs_name.str());
01062 }

const NumericVector< Number > & libMesh::System::get_sensitivity_solution ( unsigned int  i = 0  )  const [inherited]
Returns:
a reference to one of the system's solution sensitivity vectors, by default the one corresponding to the first parameter.

Definition at line 925 of file system.C.

References libMesh::System::get_vector().

00926 {
00927   std::ostringstream sensitivity_name;
00928   sensitivity_name << "sensitivity_solution" << i;
00929 
00930   return this->get_vector(sensitivity_name.str());
00931 }

NumericVector< Number > & libMesh::System::get_sensitivity_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's solution sensitivity vectors, by default the one corresponding to the first parameter.

Definition at line 915 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::sensitivity_solve().

00916 {
00917   std::ostringstream sensitivity_name;
00918   sensitivity_name << "sensitivity_solution" << i;
00919 
00920   return this->get_vector(sensitivity_name.str());
00921 }

ShellMatrix<Number>* libMesh::LinearImplicitSystem::get_shell_matrix ( void   )  [inline, inherited]

Returns a pointer to the currently attached shell matrix, if any, or NULL else.

Definition at line 182 of file linear_implicit_system.h.

References libMesh::LinearImplicitSystem::_shell_matrix.

00182 { return _shell_matrix; }

Real libMesh::RBConstruction::get_training_tolerance (  )  [inline, inherited]

Definition at line 176 of file rb_construction.h.

References libMesh::RBConstruction::training_tolerance.

00176 { return training_tolerance; }

NumericVector< Number > & libMesh::System::get_vector ( const unsigned int  vec_num  )  [inherited]
Returns:
a writeable reference to this system's additional vector number vec_num (where the vectors are counted starting with 0).

Definition at line 838 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00839 {
00840   vectors_iterator v = vectors_begin();
00841   vectors_iterator v_end = vectors_end();
00842   unsigned int num = 0;
00843   while((num<vec_num) && (v!=v_end))
00844     {
00845       num++;
00846       ++v;
00847     }
00848   libmesh_assert (v != v_end);
00849   return *(v->second);
00850 }

const NumericVector< Number > & libMesh::System::get_vector ( const unsigned int  vec_num  )  const [inherited]
Returns:
a const reference to this system's additional vector number vec_num (where the vectors are counted starting with 0).

Definition at line 822 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00823 {
00824   const_vectors_iterator v = vectors_begin();
00825   const_vectors_iterator v_end = vectors_end();
00826   unsigned int num = 0;
00827   while((num<vec_num) && (v!=v_end))
00828     {
00829       num++;
00830       ++v;
00831     }
00832   libmesh_assert (v != v_end);
00833   return *(v->second);
00834 }

NumericVector< Number > & libMesh::System::get_vector ( const std::string &  vec_name  )  [inherited]
Returns:
a writeable reference to this system's additional vector named vec_name. Access is only granted when the vector is already properly initialized.

Definition at line 803 of file system.C.

References libMesh::System::_vectors, and libMesh::err.

00804 {
00805   // Make sure the vector exists
00806   vectors_iterator pos = _vectors.find(vec_name);
00807 
00808   if (pos == _vectors.end())
00809     {
00810       libMesh::err << "ERROR: vector "
00811                     << vec_name
00812                     << " does not exist in this system!"
00813                     << std::endl;
00814       libmesh_error();
00815     }
00816 
00817   return *(pos->second);
00818 }

const NumericVector< Number > & libMesh::System::get_vector ( const std::string &  vec_name  )  const [inherited]
Returns:
a const reference to this system's additional vector named vec_name. Access is only granted when the vector is already properly initialized.

Definition at line 784 of file system.C.

References libMesh::System::_vectors, and libMesh::err.

Referenced by libMesh::UnsteadySolver::adjoint_advance_timestep(), libMesh::UnsteadySolver::advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), libMesh::System::compare(), libMesh::UnsteadySolver::du(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_adjoint_rhs(), libMesh::System::get_adjoint_solution(), libMesh::System::get_sensitivity_rhs(), libMesh::System::get_sensitivity_solution(), libMesh::System::get_weighted_sensitivity_adjoint_solution(), libMesh::System::get_weighted_sensitivity_solution(), libMesh::NewmarkSystem::initial_conditions(), libMesh::MemorySolutionHistory::retrieve(), libMesh::UnsteadySolver::retrieve_timestep(), libMesh::TwostepTimeSolver::solve(), libMesh::FrequencySystem::solve(), libMesh::NewmarkSystem::update_rhs(), and libMesh::NewmarkSystem::update_u_v_a().

00785 {
00786   // Make sure the vector exists
00787   const_vectors_iterator pos = _vectors.find(vec_name);
00788 
00789   if (pos == _vectors.end())
00790     {
00791       libMesh::err << "ERROR: vector "
00792                     << vec_name
00793                     << " does not exist in this system!"
00794                     << std::endl;
00795       libmesh_error();
00796     }
00797 
00798   return *(pos->second);
00799 }

const NumericVector< Number > & libMesh::System::get_weighted_sensitivity_adjoint_solution ( unsigned int  i = 0  )  const [inherited]
Returns:
a reference to one of the system's weighted sensitivity adjoint solution vectors, by default the one corresponding to the first qoi.

Definition at line 1006 of file system.C.

References libMesh::System::get_vector().

01007 {
01008   std::ostringstream adjoint_name;
01009   adjoint_name << "weighted_sensitivity_adjoint_solution" << i;
01010 
01011   return this->get_vector(adjoint_name.str());
01012 }

NumericVector< Number > & libMesh::System::get_weighted_sensitivity_adjoint_solution ( unsigned int  i = 0  )  [inherited]
Returns:
a reference to one of the system's weighted sensitivity adjoint solution vectors, by default the one corresponding to the first qoi.

Definition at line 996 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

00997 {
00998   std::ostringstream adjoint_name;
00999   adjoint_name << "weighted_sensitivity_adjoint_solution" << i;
01000 
01001   return this->get_vector(adjoint_name.str());
01002 }

const NumericVector< Number > & libMesh::System::get_weighted_sensitivity_solution (  )  const [inherited]
Returns:
a reference to the solution of the last weighted sensitivity solve

Definition at line 949 of file system.C.

References libMesh::System::get_vector().

00950 {
00951   return this->get_vector("weighted_sensitivity_solution");
00952 }

NumericVector< Number > & libMesh::System::get_weighted_sensitivity_solution (  )  [inherited]
Returns:
a reference to the solution of the last weighted sensitivity solve

Definition at line 942 of file system.C.

References libMesh::System::get_vector().

Referenced by libMesh::ImplicitSystem::weighted_sensitivity_solve().

00943 {
00944   return this->get_vector("weighted_sensitivity_solution");
00945 }

virtual bool libMesh::RBEIMConstruction::greedy_termination_test ( Real  training_greedy_error,
int  count 
) [protected, virtual]

Function that indicates when to terminate the Greedy basis training. Overload in subclasses to specialize.

Reimplemented from libMesh::RBConstruction.

bool libMesh::System::has_variable ( const std::string &  var  )  const [inherited]
Returns:
true if a variable named var exists in this System

Definition at line 1232 of file system.C.

References libMesh::System::_variable_numbers.

Referenced by libMesh::GMVIO::copy_nodal_solution().

01233 {
01234   return _variable_numbers.count(var);
01235 }

bool libMesh::ImplicitSystem::have_matrix ( const std::string &  mat_name  )  const [inline, inherited]
Returns:
true if this System has a matrix associated with the given name, false otherwise.

Definition at line 378 of file implicit_system.h.

References libMesh::ImplicitSystem::_matrices.

Referenced by libMesh::ImplicitSystem::add_matrix(), and libMesh::EigenTimeSolver::init().

00379 {
00380   return (_matrices.count(mat_name));
00381 }

bool libMesh::System::have_vector ( const std::string &  vec_name  )  const [inline, inherited]
Returns:
true if this System has a vector associated with the given name, false otherwise.

Definition at line 2018 of file system.h.

References libMesh::System::_vectors.

Referenced by libMesh::System::add_vector(), and libMesh::System::remove_vector().

02019 {
02020   return (_vectors.count(vec_name));
02021 }

void libMesh::System::identify_variable_groups ( const bool  ivg  )  [inline, inherited]

Toggle automatic VariableGroup identification.

Definition at line 2002 of file system.h.

References libMesh::System::_identify_variable_groups.

02003 {
02004   _identify_variable_groups = ivg;
02005 }

bool libMesh::System::identify_variable_groups (  )  const [inline, inherited]
Returns:
true when VariableGroup structures should be automatically identified, false otherwise.

Definition at line 1994 of file system.h.

References libMesh::System::_identify_variable_groups.

Referenced by libMesh::System::add_variable().

01995 {
01996   return _identify_variable_groups;
01997 }

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name  )  [inline, protected, inherited]

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

Definition at line 163 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, and libMesh::Threads::spin_mtx.

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

00164 {
00165   Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
00166   std::pair<unsigned int, unsigned int>& p = _counts[name];
00167 
00168   p.first++;
00169 }

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name  )  [inline, protected, inherited]

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

Definition at line 163 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, and libMesh::Threads::spin_mtx.

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

00164 {
00165   Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
00166   std::pair<unsigned int, unsigned int>& p = _counts[name];
00167 
00168   p.first++;
00169 }

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name  )  [inline, protected, inherited]

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

Definition at line 176 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, and libMesh::Threads::spin_mtx.

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

00177 {
00178   Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
00179   std::pair<unsigned int, unsigned int>& p = _counts[name];
00180 
00181   p.second++;
00182 }

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name  )  [inline, protected, inherited]

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

Definition at line 176 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, and libMesh::Threads::spin_mtx.

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

00177 {
00178   Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
00179   std::pair<unsigned int, unsigned int>& p = _counts[name];
00180 
00181   p.second++;
00182 }

void libMesh::System::init (  )  [inherited]

Initializes degrees of freedom on the current mesh. Sets the

Definition at line 223 of file system.C.

References libMesh::System::_basic_system_only, libMesh::System::init_data(), libMesh::System::n_vars(), and libMesh::System::user_initialization().

00224 {
00225   // First initialize any required data:
00226   // either only the basic System data
00227   if (_basic_system_only)
00228     System::init_data();
00229   // or all the derived class' data too
00230   else
00231     this->init_data();
00232 
00233   // If no variables have been added to this system
00234   // don't do anything
00235   if(!this->n_vars())
00236     return;
00237 
00238   // Then call the user-provided intialization function
00239   this->user_initialization();
00240 }

virtual void libMesh::RBEIMConstruction::init_context ( FEMContext &  c  )  [virtual]

Provide an implementation of init_context that is relevant to the projection calculations in load_calN_parametrized_function.

Reimplemented from libMesh::RBConstruction.

virtual void libMesh::RBEIMConstruction::init_data (  )  [protected, virtual]

Override to initialize the coupling matrix to decouple variables in this system.

Reimplemented from libMesh::RBConstructionBase< LinearImplicitSystem >.

void libMesh::ImplicitSystem::init_matrices (  )  [protected, virtual, inherited]

Initializes the matrices associated with this system.

Definition at line 105 of file implicit_system.C.

References libMesh::ImplicitSystem::_can_add_matrices, libMesh::ImplicitSystem::_matrices, libMesh::DofMap::attach_matrix(), libMesh::DofMap::compute_sparsity(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::SparseMatrix< T >::initialized(), libMesh::DofMap::is_attached(), and libMesh::ImplicitSystem::matrix.

Referenced by libMesh::ImplicitSystem::init_data().

00106 {
00107   libmesh_assert(matrix);
00108 
00109   // Check for quick return in case the system matrix
00110   // (and by extension all the matrices) has already
00111   // been initialized
00112   if (matrix->initialized())
00113     return;
00114 
00115   // Get a reference to the DofMap
00116   DofMap& dof_map = this->get_dof_map();
00117 
00118   // no chance to add other matrices
00119   _can_add_matrices = false;
00120 
00121   // Tell the matrices about the dof map, and vice versa
00122   for (matrices_iterator pos = _matrices.begin();
00123        pos != _matrices.end(); ++pos)
00124     {
00125       SparseMatrix<Number> &m = *(pos->second);
00126       libmesh_assert (!m.initialized());
00127 
00128       // We want to allow repeated init() on systems, but we don't
00129       // want to attach the same matrix to the DofMap twice
00130       if (!dof_map.is_attached(m))
00131         dof_map.attach_matrix (m);
00132     }
00133 
00134   // Compute the sparsity pattern for the current
00135   // mesh and DOF distribution.  This also updates
00136   // additional matrices, \p DofMap now knows them
00137   dof_map.compute_sparsity (this->get_mesh());
00138 
00139   // Initialize matrices
00140   for (matrices_iterator pos = _matrices.begin();
00141        pos != _matrices.end(); ++pos)
00142     pos->second->init ();
00143 
00144   // Set the additional matrices to 0.
00145   for (matrices_iterator pos = _matrices.begin();
00146        pos != _matrices.end(); ++pos)
00147     pos->second->zero ();
00148 }

virtual void libMesh::RBEIMConstruction::initialize_eim_assembly_objects (  )  [virtual]

Build a vector of ElemAssembly objects that accesses the basis functions stored in this RBEIMConstruction object. This is useful for performing the Offline stage of the Reduced Basis method where we want to use assembly functions based on this EIM approximation.

virtual void libMesh::RBParametrized::initialize_parameters ( const std::string &  parameters_filename  )  [virtual, inherited]

Initialize the parameter ranges and set current_parameters by reading in data from the file input_filename

void libMesh::RBParametrized::initialize_parameters ( const RBParametrized &  rb_parametrized  )  [inherited]

Initialize the parameter ranges and set current_parameters.

void libMesh::RBParametrized::initialize_parameters ( const RBParameters &  mu_min_in,
const RBParameters &  mu_max_in,
const RBParameters &  mu_in 
) [inherited]

Initialize the parameter ranges and set current_parameters.

void libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set (  )  [protected]

Loop over the training set and compute the parametrized function for each training index.

virtual void libMesh::RBEIMConstruction::initialize_rb_construction (  )  [virtual]

Initialize this system so that we can perform the Construction stage of the RB method.

Reimplemented from libMesh::RBConstruction.

virtual void libMesh::RBConstructionBase< LinearImplicitSystem >::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 
) [virtual, inherited]

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.

bool libMesh::System::is_adjoint_already_solved (  )  const [inline, inherited]

Accessor for the adjoint_already_solved boolean

Definition at line 359 of file system.h.

References libMesh::System::adjoint_already_solved.

Referenced by libMesh::AdjointResidualErrorEstimator::estimate_error().

00360   { return adjoint_already_solved;}

bool libMesh::RBConstruction::is_quiet (  )  const [inline, inherited]

Is the system in quiet mode?

Definition at line 196 of file rb_construction.h.

References libMesh::RBConstruction::quiet_mode.

00197     { return this->quiet_mode; }

bool libMesh::RBConstruction::is_rb_eval_initialized (  )  const [inherited]
Returns:
true if rb_eval is initialized. False, otherwise.
virtual void libMesh::RBConstruction::load_basis_function ( unsigned int  i  )  [virtual, inherited]

Load the i^th RB function into the RBConstruction solution vector.

virtual void libMesh::RBConstruction::load_rb_solution (  )  [virtual, inherited]

Load the RB solution from the most recent solve with rb_eval into this system's solution vector.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::RBConstructionBase< LinearImplicitSystem >::load_training_set ( std::map< std::string, std::vector< Number > > &  new_training_set  )  [virtual, inherited]

Overwrite the training parameters with new_training_set.

void libMesh::System::local_dof_indices ( const unsigned int  var,
std::set< dof_id_type > &  var_indices 
) const [inherited]

Fills the std::set with the degrees of freedom on the local processor corresponding the the variable number passed in.

Definition at line 1278 of file system.C.

References libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::DofMap::dof_indices(), libMesh::DofMap::end_dof(), libMesh::DofMap::first_dof(), libMesh::System::get_dof_map(), and libMesh::System::get_mesh().

Referenced by libMesh::System::discrete_var_norm().

01280 {
01281   // Make sure the set is clear
01282   var_indices.clear();
01283 
01284   std::vector<dof_id_type> dof_indices;
01285 
01286   // Begin the loop over the elements
01287   MeshBase::const_element_iterator       el     =
01288     this->get_mesh().active_local_elements_begin();
01289   const MeshBase::const_element_iterator end_el =
01290     this->get_mesh().active_local_elements_end();
01291 
01292   const dof_id_type
01293     first_local = this->get_dof_map().first_dof(),
01294     end_local   = this->get_dof_map().end_dof();
01295 
01296   for ( ; el != end_el; ++el)
01297     {
01298       const Elem* elem = *el;
01299       this->get_dof_map().dof_indices (elem, dof_indices, var);
01300 
01301       for(unsigned int i=0; i<dof_indices.size(); i++)
01302         {
01303           dof_id_type dof = dof_indices[i];
01304 
01305           //If the dof is owned by the local processor
01306           if(first_local <= dof && dof < end_local)
01307             var_indices.insert(dof_indices[i]);
01308         }
01309     }
01310 }

dof_id_type libMesh::System::n_active_dofs (  )  const [inline, inherited]

Returns the number of active degrees of freedom for this System.

Definition at line 2010 of file system.h.

References libMesh::System::n_constrained_dofs(), and libMesh::System::n_dofs().

02011 {
02012   return this->n_dofs() - this->n_constrained_dofs();
02013 }

unsigned int libMesh::System::n_components (  )  const [inline, inherited]
Returns:
the total number of scalar components in the system's variables. This will equal n_vars() in the case of all scalar-valued variables.

Definition at line 1914 of file system.h.

References libMesh::System::_variables, libMesh::Variable::first_scalar_number(), and libMesh::Variable::n_components().

Referenced by libMesh::System::add_variables(), libMesh::WrappedFunction< Output >::operator()(), and libMesh::System::project_vector().

01915 {
01916   if (_variables.empty())
01917     return 0;
01918 
01919   const Variable& last = _variables.back();
01920   return last.first_scalar_number() + last.n_components();
01921 }

dof_id_type libMesh::System::n_constrained_dofs (  )  const [inherited]
Returns:
the total number of constrained degrees of freedom in the system.

Definition at line 144 of file system.C.

References libMesh::System::_dof_map.

Referenced by libMesh::System::get_info(), and libMesh::System::n_active_dofs().

00145 {
00146 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00147 
00148   return _dof_map->n_constrained_dofs();
00149 
00150 #else
00151 
00152   return 0;
00153 
00154 #endif
00155 }

dof_id_type libMesh::System::n_dofs (  )  const [inherited]
Returns:
the number of degrees of freedom in the system

Definition at line 137 of file system.C.

References libMesh::System::_dof_map.

Referenced by libMesh::System::add_vector(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_info(), libMesh::UnsteadySolver::init_data(), libMesh::System::init_data(), libMesh::System::n_active_dofs(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::System::project_vector(), libMesh::System::read_legacy_data(), libMesh::UnsteadySolver::reinit(), libMesh::System::reinit(), and libMesh::System::restrict_vectors().

00138 {
00139   return _dof_map->n_dofs();
00140 }

unsigned int libMesh::LinearImplicitSystem::n_linear_iterations (  )  const [inline, inherited]

Returns the number of iterations taken for the most recent linear solve.

Definition at line 155 of file linear_implicit_system.h.

References libMesh::LinearImplicitSystem::_n_linear_iterations.

00155 { return _n_linear_iterations; }

dof_id_type libMesh::System::n_local_constrained_dofs (  )  const [inherited]
Returns:
the number of constrained degrees of freedom on this processor.

Definition at line 159 of file system.C.

References libMesh::System::_dof_map.

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

00160 {
00161 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00162 
00163   return _dof_map->n_local_constrained_dofs();
00164 
00165 #else
00166 
00167   return 0;
00168 
00169 #endif
00170 }

dof_id_type libMesh::System::n_local_dofs (  )  const [inherited]
Returns:
the number of degrees of freedom local to this processor

Definition at line 174 of file system.C.

References libMesh::System::_dof_map, and libMesh::processor_id().

Referenced by libMesh::System::add_vector(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_info(), libMesh::UnsteadySolver::init_data(), libMesh::System::init_data(), libMesh::System::project_vector(), libMesh::UnsteadySolver::reinit(), libMesh::System::reinit(), and libMesh::System::restrict_vectors().

00175 {
00176   return _dof_map->n_dofs_on_processor (libMesh::processor_id());
00177 }

unsigned int libMesh::ImplicitSystem::n_matrices (  )  const [inline, virtual, inherited]
Returns:
the number of matrices handled by this system

Reimplemented from libMesh::System.

Definition at line 385 of file implicit_system.h.

References libMesh::ImplicitSystem::_matrices.

00386 {
00387  return libmesh_cast_int<unsigned int>(_matrices.size());
00388 }

static unsigned int libMesh::ReferenceCounter::n_objects (  )  [inline, static, inherited]

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

Definition at line 79 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

00080   { return _n_objects; }

static unsigned int libMesh::ReferenceCounter::n_objects (  )  [inline, static, inherited]

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

Definition at line 79 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

00080   { return _n_objects; }

unsigned int libMesh::System::n_variable_groups (  )  const [inline, inherited]
Returns:
the number of VariableGroup variable groups in the system

Definition at line 1906 of file system.h.

References libMesh::System::_variable_groups.

Referenced by libMesh::System::add_variable(), libMesh::System::get_info(), and libMesh::System::init_data().

01907 {
01908   return libmesh_cast_int<unsigned int>(_variable_groups.size());
01909 }

unsigned int libMesh::System::n_vars (  )  const [inline, inherited]
Returns:
the number of variables in the system

Definition at line 1898 of file system.h.

References libMesh::System::_variables.

Referenced by libMesh::DiffContext::add_localized_vector(), libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_solution_vector(), libMesh::System::calculate_norm(), libMesh::WrappedFunction< Output >::component(), libMesh::DiffContext::DiffContext(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactSolution::ExactSolution(), libMesh::FEMContext::FEMContext(), libMesh::System::get_all_variable_numbers(), libMesh::EquationSystems::get_solution(), libMesh::System::init(), libMesh::FEMSystem::init_context(), libMesh::WrappedFunction< Output >::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::ProjectSolution::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::FEMContext::pre_fe_reinit(), libMesh::System::project_vector(), libMesh::System::re_update(), libMesh::System::read_legacy_data(), libMesh::System::read_parallel_data(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::System::reinit(), libMesh::EquationSystems::reinit(), libMesh::HPCoarsenTest::select_refinement(), libMesh::SystemSubsetBySubdomain::set_var_nums(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), and libMesh::System::zero_variable().

01899 {
01900   return libmesh_cast_int<unsigned int>(_variables.size());
01901 }

unsigned int libMesh::System::n_vectors (  )  const [inline, inherited]
Returns:
the number of vectors (in addition to the solution) handled by this system This is the size of the _vectors map

Definition at line 2026 of file system.h.

References libMesh::System::_vectors.

Referenced by libMesh::ExplicitSystem::add_system_rhs(), libMesh::System::compare(), libMesh::System::get_info(), and libMesh::System::write_header().

02027 {
02028   return libmesh_cast_int<unsigned int>(_vectors.size());
02029 }

const std::string & libMesh::System::name (  )  const [inline, inherited]
Returns:
the system name.

Definition at line 1818 of file system.h.

References libMesh::System::_sys_name.

Referenced by libMesh::System::compare(), libMesh::WrappedFunction< Output >::component(), DMLibMeshSetUpName_Private(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ExactSolution::ExactSolution(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::System::get_info(), libMesh::WrappedFunction< Output >::operator()(), libMesh::FrequencySystem::solve(), libMesh::System::user_assembly(), libMesh::System::user_constrain(), libMesh::System::user_initialization(), libMesh::System::user_QOI(), libMesh::System::user_QOI_derivative(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), and libMesh::System::write_serialized_data().

01819 {
01820   return _sys_name;
01821 }

unsigned int libMesh::System::number (  )  const [inline, inherited]
Returns:
the system number.

Definition at line 1826 of file system.h.

References libMesh::System::_sys_number.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::System::get_info(), libMesh::System::read_legacy_data(), libMesh::System::read_parallel_data(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::HPCoarsenTest::select_refinement(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_blocked_dof_objects(), and libMesh::System::zero_variable().

01827 {
01828   return _sys_number;
01829 }

Gradient libMesh::System::point_gradient ( unsigned int  var,
const Point &  p,
const Elem &  e 
) const [inherited]

Returns the gradient of the solution variable var at the physical point p in local Elem e in the mesh, similarly to point_value.

Definition at line 2110 of file system.C.

References libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< Real >::build(), libMesh::Elem::contains_point(), libMesh::System::current_solution(), libMesh::Elem::dim(), libMesh::System::get_dof_map(), libMesh::FEInterface::inverse_map(), libMesh::processor_id(), and libMesh::DofObject::processor_id().

02111 {
02112   libmesh_assert_equal_to (e.processor_id(), libMesh::processor_id());
02113 
02114   // Ensuring that the given point is really in the element is an
02115   // expensive assert, but as long as debugging is turned on we might
02116   // as well try to catch a particularly nasty potential error
02117   libmesh_assert (e.contains_point(p));
02118 
02119   // Get the dof map to get the proper indices for our computation
02120   const DofMap& dof_map = this->get_dof_map();
02121 
02122   // Need dof_indices for phi[i][j]
02123   std::vector<dof_id_type> dof_indices;
02124 
02125   // Fill in the dof_indices for our element
02126   dof_map.dof_indices (&e, dof_indices, var);
02127 
02128   // Get the no of dofs assciated with this point
02129   const unsigned int num_dofs = libmesh_cast_int<unsigned int>
02130     (dof_indices.size());
02131 
02132   FEType fe_type = dof_map.variable_type(0);
02133 
02134   // Build a FE again so we can calculate u(p)
02135   AutoPtr<FEBase> fe (FEBase::build(e.dim(), fe_type));
02136 
02137   // Map the physical co-ordinates to the master co-ordinates using the inverse_map from fe_interface.h
02138   // Build a vector of point co-ordinates to send to reinit
02139   std::vector<Point> coor(1, FEInterface::inverse_map(e.dim(), fe_type, &e, p));
02140 
02141   // Get the values of the shape function derivatives
02142   const std::vector<std::vector<RealGradient> >&  dphi = fe->get_dphi();
02143 
02144   // Reinitialize the element and compute the shape function values at coor
02145   fe->reinit (&e, &coor);
02146 
02147   // Get ready to accumulate a gradient
02148   Gradient grad_u;
02149 
02150   for (unsigned int l=0; l<num_dofs; l++)
02151     {
02152       grad_u.add_scaled (dphi[l][0], this->current_solution (dof_indices[l]));
02153     }
02154 
02155   return grad_u;
02156 }

Gradient libMesh::System::point_gradient ( unsigned int  var,
const Point &  p,
const bool  insist_on_success = true 
) const [inherited]

Returns the gradient of the solution variable var at the physical point p in the mesh, similarly to point_value.

Definition at line 2062 of file system.C.

References libMesh::Parallel::Communicator::broadcast(), libMesh::CommWorld, libMesh::PointLocatorBase::enable_out_of_mesh_mode(), libMesh::System::get_mesh(), mesh, libMesh::Parallel::Communicator::min(), libMesh::n_processors(), libMesh::processor_id(), libMesh::DofObject::processor_id(), libMesh::MeshBase::sub_point_locator(), and libMesh::Parallel::Communicator::verify().

02063 {
02064   // This function must be called on every processor; there's no
02065   // telling where in the partition p falls.
02066   parallel_only();
02067 
02068   // And every processor had better agree about which point we're
02069   // looking for
02070 #ifndef NDEBUG
02071   CommWorld.verify(p);
02072 #endif // NDEBUG
02073 
02074   // Get a reference to the mesh object associated with the system object that calls this function
02075   const MeshBase &mesh = this->get_mesh();
02076 
02077   // Use an existing PointLocator or create a new one
02078   AutoPtr<PointLocatorBase> locator_ptr = mesh.sub_point_locator();
02079   PointLocatorBase& locator = *locator_ptr;
02080 
02081   if (!insist_on_success)
02082     locator.enable_out_of_mesh_mode();
02083 
02084   // Get a pointer to the element that contains P
02085   const Elem *e = locator(p);
02086 
02087   Gradient grad_u;
02088 
02089   if (e && e->processor_id() == libMesh::processor_id())
02090     grad_u = point_gradient(var, p, *e);
02091 
02092   // If I have an element containing p, then let's let everyone know
02093   processor_id_type lowest_owner =
02094     (e && (e->processor_id() == libMesh::processor_id())) ?
02095     libMesh::processor_id() : libMesh::n_processors();
02096   CommWorld.min(lowest_owner);
02097 
02098   // Everybody should get their value from a processor that was able
02099   // to compute it.
02100   // If nobody admits owning the point, we may have a problem.
02101   if (lowest_owner != libMesh::n_processors())
02102     CommWorld.broadcast(grad_u, lowest_owner);
02103   else
02104     libmesh_assert(!insist_on_success);
02105 
02106   return grad_u;
02107 }

Tensor libMesh::System::point_hessian ( unsigned int  var,
const Point &  p,
const Elem &  e 
) const [inherited]

Returns the second derivative tensor of the solution variable var at the physical point p in local Elem e in the mesh, similarly to point_value.

Definition at line 2208 of file system.C.

References libMesh::TypeTensor< T >::add_scaled(), libMesh::FEGenericBase< Real >::build(), libMesh::Elem::contains_point(), libMesh::System::current_solution(), libMesh::Elem::dim(), libMesh::System::get_dof_map(), libMesh::FEInterface::inverse_map(), libMesh::processor_id(), and libMesh::DofObject::processor_id().

02209 {
02210   libmesh_assert_equal_to (e.processor_id(), libMesh::processor_id());
02211 
02212   // Ensuring that the given point is really in the element is an
02213   // expensive assert, but as long as debugging is turned on we might
02214   // as well try to catch a particularly nasty potential error
02215   libmesh_assert (e.contains_point(p));
02216 
02217   // Get the dof map to get the proper indices for our computation
02218   const DofMap& dof_map = this->get_dof_map();
02219 
02220   // Need dof_indices for phi[i][j]
02221   std::vector<dof_id_type> dof_indices;
02222 
02223   // Fill in the dof_indices for our element
02224   dof_map.dof_indices (&e, dof_indices, var);
02225 
02226   // Get the no of dofs assciated with this point
02227   const unsigned int num_dofs = libmesh_cast_int<unsigned int>
02228     (dof_indices.size());
02229 
02230   FEType fe_type = dof_map.variable_type(0);
02231 
02232   // Build a FE again so we can calculate u(p)
02233   AutoPtr<FEBase> fe (FEBase::build(e.dim(), fe_type));
02234 
02235   // Map the physical co-ordinates to the master co-ordinates using the inverse_map from fe_interface.h
02236   // Build a vector of point co-ordinates to send to reinit
02237   std::vector<Point> coor(1, FEInterface::inverse_map(e.dim(), fe_type, &e, p));
02238 
02239   // Get the values of the shape function derivatives
02240   const std::vector<std::vector<RealTensor> >&  d2phi = fe->get_d2phi();
02241 
02242   // Reinitialize the element and compute the shape function values at coor
02243   fe->reinit (&e, &coor);
02244 
02245   // Get ready to accumulate a hessian
02246   Tensor hess_u;
02247 
02248   for (unsigned int l=0; l<num_dofs; l++)
02249     {
02250       hess_u.add_scaled (d2phi[l][0], this->current_solution (dof_indices[l]));
02251     }
02252 
02253   return hess_u;
02254 }

Tensor libMesh::System::point_hessian ( unsigned int  var,
const Point &  p,
const bool  insist_on_success = true 
) const [inherited]

Returns the second derivative tensor of the solution variable var at the physical point p in the mesh, similarly to point_value.

Definition at line 2161 of file system.C.

References libMesh::Parallel::Communicator::broadcast(), libMesh::CommWorld, libMesh::PointLocatorBase::enable_out_of_mesh_mode(), libMesh::System::get_mesh(), mesh, libMesh::Parallel::Communicator::min(), libMesh::n_processors(), libMesh::processor_id(), libMesh::DofObject::processor_id(), libMesh::MeshBase::sub_point_locator(), and libMesh::Parallel::Communicator::verify().

02162 {
02163   // This function must be called on every processor; there's no
02164   // telling where in the partition p falls.
02165   parallel_only();
02166 
02167   // And every processor had better agree about which point we're
02168   // looking for
02169 #ifndef NDEBUG
02170   CommWorld.verify(p);
02171 #endif // NDEBUG
02172 
02173   // Get a reference to the mesh object associated with the system object that calls this function
02174   const MeshBase &mesh = this->get_mesh();
02175 
02176   // Use an existing PointLocator or create a new one
02177   AutoPtr<PointLocatorBase> locator_ptr = mesh.sub_point_locator();
02178   PointLocatorBase& locator = *locator_ptr;
02179 
02180   if (!insist_on_success)
02181     locator.enable_out_of_mesh_mode();
02182 
02183   // Get a pointer to the element that contains P
02184   const Elem *e = locator(p);
02185 
02186   Tensor hess_u;
02187 
02188   if (e && e->processor_id() == libMesh::processor_id())
02189     hess_u = point_hessian(var, p, *e);
02190 
02191   // If I have an element containing p, then let's let everyone know
02192   processor_id_type lowest_owner =
02193     (e && (e->processor_id() == libMesh::processor_id())) ?
02194     libMesh::processor_id() : libMesh::n_processors();
02195   CommWorld.min(lowest_owner);
02196 
02197   // Everybody should get their value from a processor that was able
02198   // to compute it.
02199   // If nobody admits owning the point, we may have a problem.
02200   if (lowest_owner != libMesh::n_processors())
02201     CommWorld.broadcast(hess_u, lowest_owner);
02202   else
02203     libmesh_assert(!insist_on_success);
02204 
02205   return hess_u;
02206 }

Number libMesh::System::point_value ( unsigned int  var,
const Point &  p,
const Elem &  e 
) const [inherited]

Returns the value of the solution variable var at the physical point p contained in local Elem e

This version of point_value can be run in serial, but assumes e is in the local mesh partition.

Definition at line 2012 of file system.C.

References libMesh::FEGenericBase< Real >::build(), libMesh::Elem::contains_point(), libMesh::System::current_solution(), libMesh::Elem::dim(), libMesh::System::get_dof_map(), libMesh::FEInterface::inverse_map(), libMesh::processor_id(), and libMesh::DofObject::processor_id().

02013 {
02014   libmesh_assert_equal_to (e.processor_id(), libMesh::processor_id());
02015 
02016   // Ensuring that the given point is really in the element is an
02017   // expensive assert, but as long as debugging is turned on we might
02018   // as well try to catch a particularly nasty potential error
02019   libmesh_assert (e.contains_point(p));
02020 
02021   // Get the dof map to get the proper indices for our computation
02022   const DofMap& dof_map = this->get_dof_map();
02023 
02024   // Need dof_indices for phi[i][j]
02025   std::vector<dof_id_type> dof_indices;
02026 
02027   // Fill in the dof_indices for our element
02028   dof_map.dof_indices (&e, dof_indices, var);
02029 
02030   // Get the no of dofs assciated with this point
02031   const unsigned int num_dofs = libmesh_cast_int<unsigned int>
02032     (dof_indices.size());
02033 
02034   FEType fe_type = dof_map.variable_type(0);
02035 
02036   // Build a FE so we can calculate u(p)
02037   AutoPtr<FEBase> fe (FEBase::build(e.dim(), fe_type));
02038 
02039   // Map the physical co-ordinates to the master co-ordinates using the inverse_map from fe_interface.h
02040   // Build a vector of point co-ordinates to send to reinit
02041   std::vector<Point> coor(1, FEInterface::inverse_map(e.dim(), fe_type, &e, p));
02042 
02043   // Get the shape function values
02044   const std::vector<std::vector<Real> >& phi = fe->get_phi();
02045 
02046   // Reinitialize the element and compute the shape function values at coor
02047   fe->reinit (&e, &coor);
02048 
02049   // Get ready to accumulate a value
02050   Number u = 0;
02051 
02052   for (unsigned int l=0; l<num_dofs; l++)
02053     {
02054       u += phi[l][0]*this->current_solution (dof_indices[l]);
02055     }
02056 
02057   return u;
02058 }

Number libMesh::System::point_value ( unsigned int  var,
const Point &  p,
const bool  insist_on_success = true 
) const [inherited]

Returns the value of the solution variable var at the physical point p in the mesh, without knowing a priori which element contains p.

Note that this function uses MeshBase::sub_point_locator(); users may or may not want to call MeshBase::clear_point_locator() afterward. Also, point_locator() is expensive (N log N for initial construction, log N for evaluations). Avoid using this function in any context where you are already looping over elements.

Because the element containing p may lie on any processor, this function is parallel-only.

By default this method expects the point to reside inside the domain and will abort if no element can be found which contains . The optional parameter insist_on_success can be set to false to allow the method to return 0 when the point is not located.

Definition at line 1965 of file system.C.

References libMesh::Parallel::Communicator::broadcast(), libMesh::CommWorld, libMesh::PointLocatorBase::enable_out_of_mesh_mode(), libMesh::System::get_mesh(), mesh, libMesh::Parallel::Communicator::min(), libMesh::n_processors(), libMesh::processor_id(), libMesh::DofObject::processor_id(), libMesh::MeshBase::sub_point_locator(), and libMesh::Parallel::Communicator::verify().

01966 {
01967   // This function must be called on every processor; there's no
01968   // telling where in the partition p falls.
01969   parallel_only();
01970 
01971   // And every processor had better agree about which point we're
01972   // looking for
01973 #ifndef NDEBUG
01974   CommWorld.verify(p);
01975 #endif // NDEBUG
01976 
01977   // Get a reference to the mesh object associated with the system object that calls this function
01978   const MeshBase &mesh = this->get_mesh();
01979 
01980   // Use an existing PointLocator or create a new one
01981   AutoPtr<PointLocatorBase> locator_ptr = mesh.sub_point_locator();
01982   PointLocatorBase& locator = *locator_ptr;
01983 
01984   if (!insist_on_success)
01985     locator.enable_out_of_mesh_mode();
01986 
01987   // Get a pointer to the element that contains P
01988   const Elem *e = locator(p);
01989 
01990   Number u = 0;
01991 
01992   if (e && e->processor_id() == libMesh::processor_id())
01993     u = point_value(var, p, *e);
01994 
01995   // If I have an element containing p, then let's let everyone know
01996   processor_id_type lowest_owner =
01997     (e && (e->processor_id() == libMesh::processor_id())) ?
01998     libMesh::processor_id() : libMesh::n_processors();
01999   CommWorld.min(lowest_owner);
02000 
02001   // Everybody should get their value from a processor that was able
02002   // to compute it.
02003   // If nobody admits owning the point, we have a problem.
02004   if (lowest_owner != libMesh::n_processors())
02005     CommWorld.broadcast(u, lowest_owner);
02006   else
02007     libmesh_assert(!insist_on_success);
02008 
02009   return u;
02010 }

void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out  )  [static, inherited]

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

Definition at line 88 of file reference_counter.C.

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

00089 {
00090   if( _enable_print_counter ) out_stream << ReferenceCounter::get_info();
00091 }

void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out  )  [static, inherited]

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

Definition at line 88 of file reference_counter.C.

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

00089 {
00090   if( _enable_print_counter ) out_stream << ReferenceCounter::get_info();
00091 }

virtual void libMesh::RBEIMConstruction::print_info (  )  [virtual]

Print out info that describes the current setup of this RBConstruction.

Reimplemented from libMesh::RBConstruction.

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

Print the current parameters.

virtual void libMesh::RBEIMConstruction::process_parameters_file ( const std::string &  parameters_filename  )  [virtual]

Read parameters in from file and set up this system accordingly.

Reimplemented from libMesh::RBConstruction.

void libMesh::System::project_solution ( Number   fptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
Gradient   gptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
const Parameters &  parameters 
) const [inherited]

Projects arbitrary functions onto the current solution. The function value fptr and its gradient gptr are represented by function pointers. A gradient gptr is only required/used for projecting onto finite element spaces with continuous derivatives.

This method projects an arbitrary function onto the solution via L2 projections and nodal interpolations on each element.

Definition at line 450 of file system_projection.C.

References libMesh::System::project_solution().

00459 {
00460   WrappedFunction<Number> f(*this, fptr, &parameters);
00461   WrappedFunction<Gradient> g(*this, gptr, &parameters);
00462   this->project_solution(&f, &g);
00463 }

void libMesh::System::project_solution ( FEMFunctionBase< Number > *  f,
FEMFunctionBase< Gradient > *  g = NULL 
) const [inherited]

Projects arbitrary functions onto the current solution. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitrary function onto the solution via L2 projections and nodal interpolations on each element.

Definition at line 483 of file system_projection.C.

References libMesh::System::_dof_map, libMesh::System::current_local_solution, libMesh::System::project_vector(), and libMesh::System::solution.

00485 {
00486   this->project_vector(*solution, f, g);
00487 
00488   solution->localize(*current_local_solution, _dof_map->get_send_list());
00489 }

void libMesh::System::project_solution ( FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = NULL 
) const [inherited]

Projects arbitrary functions onto the current solution. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitrary function onto the solution via L2 projections and nodal interpolations on each element.

Definition at line 470 of file system_projection.C.

References libMesh::System::_dof_map, libMesh::System::current_local_solution, libMesh::System::project_vector(), and libMesh::System::solution.

Referenced by libMesh::System::project_solution().

00472 {
00473   this->project_vector(*solution, f, g);
00474 
00475   solution->localize(*current_local_solution, _dof_map->get_send_list());
00476 }

bool& libMesh::System::project_solution_on_reinit ( void   )  [inline, inherited]

Tells the System whether or not to project the solution vector onto new grids when the system is reinitialized. The solution will be projected unless project_solution_on_reinit() = false is called.

Definition at line 761 of file system.h.

References libMesh::System::_solution_projection.

Referenced by libMesh::AdjointRefinementEstimator::estimate_error(), and libMesh::MemorySolutionHistory::store().

00762     { return _solution_projection; }

void libMesh::System::project_vector ( const NumericVector< Number > &  old_v,
NumericVector< Number > &  new_v 
) const [protected, inherited]

Projects the vector defined on the old mesh onto the new mesh. The original vector is unchanged and the new vector is passed through the second argument.

This method projects the vector via L2 projections or nodal interpolations on each element.

This method projects a solution from an old mesh to a current, refined mesh. The input vector old_v gives the solution on the old mesh, while the new_v gives the solution (to be computed) on the new mesh.

Definition at line 272 of file system_projection.C.

References libMesh::NumericVector< T >::clear(), libMesh::NumericVector< T >::close(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::err, libMesh::FEType::family, libMesh::AutoPtr< Tp >::get(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMeshEnums::GHOSTED, libMesh::NumericVector< T >::init(), libMesh::NumericVector< T >::local_size(), libMesh::NumericVector< T >::localize(), libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), libMesh::n_processors(), libMesh::System::n_vars(), libMeshEnums::PARALLEL, libMesh::Threads::parallel_for(), libMesh::Threads::parallel_reduce(), libMesh::processor_id(), libMeshEnums::SCALAR, libMesh::DofMap::SCALAR_dof_indices(), libMesh::BuildProjectionList::send_list, libMeshEnums::SERIAL, libMesh::NumericVector< T >::set(), libMesh::NumericVector< T >::size(), libMesh::Variable::type(), libMesh::NumericVector< T >::type(), libMesh::BuildProjectionList::unique(), and libMesh::System::variable().

00274 {
00275   START_LOG ("project_vector()", "System");
00276 
00283   new_v.clear();
00284 
00285 #ifdef LIBMESH_ENABLE_AMR
00286 
00287   // Resize the new vector and get a serial version.
00288   NumericVector<Number> *new_vector_ptr = NULL;
00289   AutoPtr<NumericVector<Number> > new_vector_built;
00290   NumericVector<Number> *local_old_vector;
00291   AutoPtr<NumericVector<Number> > local_old_vector_built;
00292   const NumericVector<Number> *old_vector_ptr = NULL;
00293 
00294   ConstElemRange active_local_elem_range
00295     (this->get_mesh().active_local_elements_begin(),
00296      this->get_mesh().active_local_elements_end());
00297 
00298   // If the old vector was uniprocessor, make the new
00299   // vector uniprocessor
00300   if (old_v.type() == SERIAL)
00301     {
00302       new_v.init (this->n_dofs(), false, SERIAL);
00303       new_vector_ptr = &new_v;
00304       old_vector_ptr = &old_v;
00305     }
00306 
00307   // Otherwise it is a parallel, distributed vector, which
00308   // we need to localize.
00309   else if (old_v.type() == PARALLEL)
00310     {
00311       // Build a send list for efficient localization
00312       BuildProjectionList projection_list(*this);
00313       Threads::parallel_reduce (active_local_elem_range,
00314                                 projection_list);
00315 
00316       // Create a sorted, unique send_list
00317       projection_list.unique();
00318 
00319       new_v.init (this->n_dofs(), this->n_local_dofs(), false, PARALLEL);
00320       new_vector_built = NumericVector<Number>::build();
00321       local_old_vector_built = NumericVector<Number>::build();
00322       new_vector_ptr = new_vector_built.get();
00323       local_old_vector = local_old_vector_built.get();
00324       new_vector_ptr->init(this->n_dofs(), false, SERIAL);
00325       local_old_vector->init(old_v.size(), false, SERIAL);
00326       old_v.localize(*local_old_vector, projection_list.send_list);
00327       local_old_vector->close();
00328       old_vector_ptr = local_old_vector;
00329     }
00330   else if (old_v.type() == GHOSTED)
00331     {
00332       // Build a send list for efficient localization
00333       BuildProjectionList projection_list(*this);
00334       Threads::parallel_reduce (active_local_elem_range,
00335                                 projection_list);
00336 
00337       // Create a sorted, unique send_list
00338       projection_list.unique();
00339 
00340       new_v.init (this->n_dofs(), this->n_local_dofs(),
00341                   this->get_dof_map().get_send_list(), false, GHOSTED);
00342 
00343       local_old_vector_built = NumericVector<Number>::build();
00344       new_vector_ptr = &new_v;
00345       local_old_vector = local_old_vector_built.get();
00346       local_old_vector->init(old_v.size(), old_v.local_size(),
00347                              projection_list.send_list, false, GHOSTED);
00348       old_v.localize(*local_old_vector, projection_list.send_list);
00349       local_old_vector->close();
00350       old_vector_ptr = local_old_vector;
00351     }
00352   else // unknown old_v.type()
00353     {
00354       libMesh::err << "ERROR: Unknown old_v.type() == " << old_v.type()
00355                     << std::endl;
00356       libmesh_error();
00357     }
00358 
00359   // Note that the above will have zeroed the new_vector.
00360   // Just to be sure, assert that new_vector_ptr and old_vector_ptr
00361   // were successfully set before trying to deref them.
00362   libmesh_assert(new_vector_ptr);
00363   libmesh_assert(old_vector_ptr);
00364 
00365   NumericVector<Number> &new_vector = *new_vector_ptr;
00366   const NumericVector<Number> &old_vector = *old_vector_ptr;
00367 
00368   Threads::parallel_for (active_local_elem_range,
00369                          ProjectVector(*this,
00370                                        old_vector,
00371                                        new_vector)
00372                          );
00373 
00374   // Copy the SCALAR dofs from old_vector to new_vector
00375   // Note: We assume that all SCALAR dofs are on the
00376   // processor with highest ID
00377   if(libMesh::processor_id() == (libMesh::n_processors()-1))
00378   {
00379     const DofMap& dof_map = this->get_dof_map();
00380     for (unsigned int var=0; var<this->n_vars(); var++)
00381       if(this->variable(var).type().family == SCALAR)
00382         {
00383           // We can just map SCALAR dofs directly across
00384           std::vector<dof_id_type> new_SCALAR_indices, old_SCALAR_indices;
00385           dof_map.SCALAR_dof_indices (new_SCALAR_indices, var, false);
00386           dof_map.SCALAR_dof_indices (old_SCALAR_indices, var, true);
00387           const unsigned int new_n_dofs = 
00388             libmesh_cast_int<unsigned int>(new_SCALAR_indices.size());
00389 
00390           for (unsigned int i=0; i<new_n_dofs; i++)
00391           {
00392             new_vector.set( new_SCALAR_indices[i], old_vector(old_SCALAR_indices[i]) );
00393           }
00394         }
00395   }
00396 
00397   new_vector.close();
00398 
00399   // If the old vector was serial, we probably need to send our values
00400   // to other processors
00401   //
00402   // FIXME: I'm not sure how to make a NumericVector do that without
00403   // creating a temporary parallel vector to use localize! - RHS
00404   if (old_v.type() == SERIAL)
00405     {
00406       AutoPtr<NumericVector<Number> > dist_v = NumericVector<Number>::build();
00407       dist_v->init(this->n_dofs(), this->n_local_dofs(), false, PARALLEL);
00408       dist_v->close();
00409 
00410       for (dof_id_type i=0; i!=dist_v->size(); i++)
00411         if (new_vector(i) != 0.0)
00412           dist_v->set(i, new_vector(i));
00413 
00414       dist_v->close();
00415 
00416       dist_v->localize (new_v, this->get_dof_map().get_send_list());
00417       new_v.close();
00418     }
00419   // If the old vector was parallel, we need to update it
00420   // and free the localized copies
00421   else if (old_v.type() == PARALLEL)
00422     {
00423       // We may have to set dof values that this processor doesn't
00424       // own in certain special cases, like LAGRANGE FIRST or
00425       // HERMITE THIRD elements on second-order meshes
00426       for (dof_id_type i=0; i!=new_v.size(); i++)
00427         if (new_vector(i) != 0.0)
00428           new_v.set(i, new_vector(i));
00429       new_v.close();
00430     }
00431 
00432   this->get_dof_map().enforce_constraints_exactly(*this, &new_v);
00433 
00434 #else
00435 
00436   // AMR is disabled: simply copy the vector
00437   new_v = old_v;
00438 
00439 #endif // #ifdef LIBMESH_ENABLE_AMR
00440 
00441   STOP_LOG("project_vector()", "System");
00442 }

void libMesh::System::project_vector ( NumericVector< Number > &  vector  )  const [protected, inherited]

Projects the vector defined on the old mesh onto the new mesh.

Definition at line 255 of file system_projection.C.

References libMesh::NumericVector< T >::clone(), and libMesh::System::project_vector().

00256 {
00257   // Create a copy of the vector, which currently
00258   // contains the old data.
00259   AutoPtr<NumericVector<Number> >
00260     old_vector (vector.clone());
00261 
00262   // Project the old vector to the new vector
00263   this->project_vector (*old_vector, vector);
00264 }

void libMesh::System::project_vector ( Number   fptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
Gradient   gptrconst Point &p,const Parameters &parameters,const std::string &sys_name,const std::string &unknown_name,
const Parameters &  parameters,
NumericVector< Number > &  new_vector 
) const [inherited]

Projects arbitrary functions onto a vector of degree of freedom values for the current system. The function value fptr and its gradient gptr are represented by function pointers. A gradient gptr is only required/used for projecting onto finite element spaces with continuous derivatives.

This method projects an arbitrary function via L2 projections and nodal interpolations on each element.

Definition at line 496 of file system_projection.C.

References libMesh::System::project_vector().

00506 {
00507   WrappedFunction<Number> f(*this, fptr, &parameters);
00508   WrappedFunction<Gradient> g(*this, gptr, &parameters);
00509   this->project_vector(new_vector, &f, &g);
00510 }

void libMesh::System::project_vector ( NumericVector< Number > &  new_vector,
FEMFunctionBase< Number > *  f,
FEMFunctionBase< Gradient > *  g = NULL 
) const [inherited]

Projects arbitrary functions onto a vector of degree of freedom values for the current system. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitrary function via L2 projections and nodal interpolations on each element.

Definition at line 579 of file system_projection.C.

References libMesh::MeshBase::active_local_elements_begin(), libMesh::NumericVector< T >::close(), libMesh::FEMFunctionBase< Output >::component(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::FEType::family, libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::n_processors(), libMesh::System::n_vars(), libMesh::Threads::parallel_for(), libMesh::FEMContext::pre_fe_reinit(), libMesh::processor_id(), libMeshEnums::SCALAR, libMesh::DofMap::SCALAR_dof_indices(), libMesh::NumericVector< T >::set(), libMesh::System::time, libMesh::Variable::type(), libMesh::System::variable(), and libMesh::System::variable_scalar_number().

00582 {
00583   START_LOG ("project_fem_vector()", "System");
00584 
00585   Threads::parallel_for
00586     (ConstElemRange (this->get_mesh().active_local_elements_begin(),
00587                      this->get_mesh().active_local_elements_end() ),
00588      ProjectFEMSolution(*this, f, g, new_vector)
00589     );
00590 
00591   // Also, load values into the SCALAR dofs
00592   // Note: We assume that all SCALAR dofs are on the
00593   // processor with highest ID
00594   if(libMesh::processor_id() == (libMesh::n_processors()-1))
00595   {
00596     // FIXME: Do we want to first check for SCALAR vars before building this? [PB]
00597     FEMContext context( *this );
00598 
00599     const DofMap& dof_map = this->get_dof_map();
00600     for (unsigned int var=0; var<this->n_vars(); var++)
00601       if(this->variable(var).type().family == SCALAR)
00602         {
00603           // FIXME: We reinit with an arbitrary element in case the user
00604           //        doesn't override FEMFunctionBase::component. Is there
00605           //        any use case we're missing? [PB]
00606           Elem *el = const_cast<Elem *>(*(this->get_mesh().active_local_elements_begin()));
00607           context.pre_fe_reinit( *this, el );
00608 
00609           std::vector<dof_id_type> SCALAR_indices;
00610           dof_map.SCALAR_dof_indices (SCALAR_indices, var);
00611           const unsigned int n_SCALAR_dofs =
00612             libmesh_cast_int<unsigned int>(SCALAR_indices.size());
00613 
00614           for (unsigned int i=0; i<n_SCALAR_dofs; i++)
00615           {
00616             const dof_id_type global_index = SCALAR_indices[i];
00617             const unsigned int component_index =
00618               this->variable_scalar_number(var,i);
00619             
00620             new_vector.set(global_index, f->component(context, component_index, Point(), this->time));
00621           }
00622         }
00623   }
00624 
00625   new_vector.close();
00626 
00627 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00628   this->get_dof_map().enforce_constraints_exactly(*this, &new_vector);
00629 #endif
00630 
00631   STOP_LOG("project_fem_vector()", "System");
00632 }

void libMesh::System::project_vector ( NumericVector< Number > &  new_vector,
FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = NULL 
) const [inherited]

Projects arbitrary functions onto a vector of degree of freedom values for the current system. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives. If non-default Parameters are to be used, they can be provided in the parameters argument.

This method projects an arbitrary function via L2 projections and nodal interpolations on each element.

Definition at line 516 of file system_projection.C.

References libMesh::NumericVector< T >::close(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::FEType::family, libMesh::System::get_dof_map(), libMesh::System::get_equation_systems(), libMesh::System::get_mesh(), libMesh::System::n_components(), libMesh::n_processors(), libMesh::System::n_vars(), libMesh::Threads::parallel_for(), libMesh::processor_id(), libMeshEnums::SCALAR, libMesh::DofMap::SCALAR_dof_indices(), libMesh::NumericVector< T >::set(), libMesh::System::time, libMesh::Variable::type(), libMesh::System::variable(), and libMesh::System::variable_scalar_number().

Referenced by libMesh::System::project_solution(), libMesh::System::project_vector(), and libMesh::System::restrict_vectors().

00519 {
00520   START_LOG ("project_vector()", "System");
00521 
00522   Threads::parallel_for
00523     (ConstElemRange (this->get_mesh().active_local_elements_begin(),
00524                      this->get_mesh().active_local_elements_end() ),
00525      ProjectSolution(*this, f, g,
00526                      this->get_equation_systems().parameters,
00527                      new_vector)
00528     );
00529 
00530   // Also, load values into the SCALAR dofs
00531   // Note: We assume that all SCALAR dofs are on the
00532   // processor with highest ID
00533   if(libMesh::processor_id() == (libMesh::n_processors()-1))
00534   {
00535     // We get different scalars as different
00536     // components from a new-style f functor.
00537     DenseVector<Number> fout(this->n_components());
00538     bool filled_fout = false;
00539 
00540     const DofMap& dof_map = this->get_dof_map();
00541     for (unsigned int var=0; var<this->n_vars(); var++)
00542       if(this->variable(var).type().family == SCALAR)
00543         {
00544           if (!filled_fout)
00545             {
00546               (*f) (Point(), this->time, fout);
00547               filled_fout = true;
00548             }
00549 
00550           std::vector<dof_id_type> SCALAR_indices;
00551           dof_map.SCALAR_dof_indices (SCALAR_indices, var);
00552           const unsigned int n_SCALAR_dofs = 
00553             libmesh_cast_int<unsigned int>(SCALAR_indices.size());
00554 
00555           for (unsigned int i=0; i<n_SCALAR_dofs; i++)
00556           {
00557             const dof_id_type global_index = SCALAR_indices[i];
00558             const unsigned int component_index =
00559               this->variable_scalar_number(var,i);
00560             new_vector.set(global_index, fout(component_index));
00561           }
00562         }
00563   }
00564 
00565   new_vector.close();
00566 
00567 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00568   this->get_dof_map().enforce_constraints_exactly(*this, &new_vector);
00569 #endif
00570 
00571   STOP_LOG("project_vector()", "System");
00572 }

void libMesh::System::prolong_vectors (  )  [virtual, inherited]

Prolong vectors after the mesh has refined

Definition at line 368 of file system.C.

References libMesh::System::restrict_vectors().

Referenced by libMesh::EquationSystems::reinit().

00369 {
00370 #ifdef LIBMESH_ENABLE_AMR
00371   // Currently project_vector handles both restriction and prolongation
00372   this->restrict_vectors();
00373 #endif
00374 }

void libMesh::ImplicitSystem::qoi_parameter_hessian ( const QoISet &  qoi_indices,
const ParameterVector &  parameters,
SensitivityData &  hessian 
) [virtual, inherited]

For each of the system's quantities of interest q in qoi[qoi_indices], and for a vector of parameters p, the parameter sensitivity Hessian H_ij is defined as H_ij = (d^2 q)/(d p_i d p_j) This Hessian is the output of this method, where for each q_i, H_jk is stored in hessian.second_derivative(i,j,k).

Reimplemented from libMesh::System.

Definition at line 1080 of file implicit_system.C.

References libMesh::SensitivityData::allocate_hessian_data(), libMesh::QoISet::has_index(), libMesh::Real, libMesh::SensitivityData::second_derivative(), libMesh::ParameterVector::size(), and libMesh::TOLERANCE.

01083 {
01084   // We currently get partial derivatives via finite differencing
01085   const Real delta_p = TOLERANCE;
01086 
01087   // We'll use one temporary vector for matrix-vector-vector products
01088   AutoPtr<NumericVector<Number> > tempvec = this->solution->zero_clone();
01089 
01090   // And another temporary vector to hold a copy of the true solution
01091   // so we can safely perturb this->solution.
01092   AutoPtr<NumericVector<Number> > oldsolution = this->solution->clone();
01093 
01094   const unsigned int Np = libmesh_cast_int<unsigned int>
01095     (parameters.size());
01096   const unsigned int Nq = libmesh_cast_int<unsigned int>
01097     (qoi.size());
01098 
01099   // For each quantity of interest q, the parameter sensitivity
01100   // Hessian is defined as q''_{kl} = {d^2 q}/{d p_k d p_l}.
01101   //
01102   // We calculate it from values and partial derivatives of the
01103   // quantity of interest function Q, solution u, adjoint solution z,
01104   // and residual R, as:
01105   //
01106   // q''_{kl} =
01107   // Q''_{kl} + Q''_{uk}(u)*u'_l + Q''_{ul}(u) * u'_k +
01108   // Q''_{uu}(u)*u'_k*u'_l -
01109   // R''_{kl}(u,z) -
01110   // R''_{uk}(u,z)*u'_l - R''_{ul}(u,z)*u'_k -
01111   // R''_{uu}(u,z)*u'_k*u'_l
01112   //
01113   // See the adjoints model document for more details.
01114 
01115   // We first do an adjoint solve to get z for each quantity of
01116   // interest
01117   // if we havent already or dont have an initial condition for the adjoint
01118   if (!this->is_adjoint_already_solved())
01119     {
01120       this->adjoint_solve(qoi_indices);
01121     }
01122 
01123   // And a sensitivity solve to get u_k for each parameter
01124   this->sensitivity_solve(parameters);
01125 
01126   // Get ready to fill in second derivatives:
01127   sensitivities.allocate_hessian_data(qoi_indices, *this, parameters);
01128 
01129   for (unsigned int k=0; k != Np; ++k)
01130     {
01131       Number old_parameterk = *parameters[k];
01132 
01133       // The Hessian is symmetric, so we just calculate the lower
01134       // triangle and the diagonal, and we get the upper triangle from
01135       // the transpose of the lower
01136 
01137       for (unsigned int l=0; l != k+1; ++l)
01138         {
01139           // The second partial derivatives with respect to parameters
01140           // are all calculated via a central finite difference
01141           // stencil:
01142           // F''_{kl} ~= (F(p+dp*e_k+dp*e_l) - F(p+dp*e_k-dp*e_l) -
01143           //              F(p-dp*e_k+dp*e_l) + F(p-dp*e_k-dp*e_l))/(4*dp^2)
01144           // We will add Q''_{kl}(u) and subtract R''_{kl}(u,z) at the
01145           // same time.
01146           //
01147           // We have to be careful with the perturbations to handle
01148           // the k=l case
01149 
01150           Number old_parameterl = *parameters[l];
01151 
01152           *parameters[k] += delta_p;
01153           *parameters[l] += delta_p;
01154           this->assemble_qoi(qoi_indices);
01155           this->assembly(true, false);
01156           this->rhs->close();
01157           std::vector<Number> partial2q_term = this->qoi;
01158           std::vector<Number> partial2R_term(this->qoi.size());
01159           for (unsigned int i=0; i != Nq; ++i)
01160             if (qoi_indices.has_index(i))
01161               partial2R_term[i] = this->rhs->dot(this->get_adjoint_solution(i));
01162 
01163           *parameters[l] -= 2.*delta_p;
01164           this->assemble_qoi(qoi_indices);
01165           this->assembly(true, false);
01166           this->rhs->close();
01167           for (unsigned int i=0; i != Nq; ++i)
01168             if (qoi_indices.has_index(i))
01169               {
01170                 partial2q_term[i] -= this->qoi[i];
01171                 partial2R_term[i] -= this->rhs->dot(this->get_adjoint_solution(i));
01172               }
01173 
01174           *parameters[k] -= 2.*delta_p;
01175           this->assemble_qoi(qoi_indices);
01176           this->assembly(true, false);
01177           this->rhs->close();
01178           for (unsigned int i=0; i != Nq; ++i)
01179             if (qoi_indices.has_index(i))
01180               {
01181                 partial2q_term[i] += this->qoi[i];
01182                 partial2R_term[i] += this->rhs->dot(this->get_adjoint_solution(i));
01183               }
01184 
01185           *parameters[l] += 2.*delta_p;
01186           this->assemble_qoi(qoi_indices);
01187           this->assembly(true, false);
01188           this->rhs->close();
01189           for (unsigned int i=0; i != Nq; ++i)
01190             if (qoi_indices.has_index(i))
01191               {
01192                 partial2q_term[i] -= this->qoi[i];
01193                 partial2R_term[i] -= this->rhs->dot(this->get_adjoint_solution(i));
01194                 partial2q_term[i] /= (4. * delta_p * delta_p);
01195                 partial2R_term[i] /= (4. * delta_p * delta_p);
01196               }
01197 
01198           for (unsigned int i=0; i != Nq; ++i)
01199             if (qoi_indices.has_index(i))
01200               {
01201                 Number current_terms = partial2q_term[i] - partial2R_term[i];
01202                 sensitivities.second_derivative(i,k,l) += current_terms;
01203                 if (k != l)
01204                   sensitivities.second_derivative(i,l,k) += current_terms;
01205               }
01206 
01207           // Don't leave the parameters perturbed
01208           *parameters[l] = old_parameterl;
01209           *parameters[k] = old_parameterk;
01210         }
01211 
01212       // We get (partial q / partial u) and
01213       // (partial R / partial u) from the user, but centrally
01214       // difference to get q_uk and R_uk terms:
01215       // (partial^2 q / partial u partial k)
01216       // q_uk*u'_l = (q_u(p+dp*e_k)*u'_l - q_u(p-dp*e_k)*u'_l)/(2*dp)
01217       // R_uk*z*u'_l = (R_u(p+dp*e_k)*z*u'_l - R_u(p-dp*e_k)*z*u'_l)/(2*dp)
01218       //
01219       // To avoid creating Nq temporary vectors, we add these
01220       // subterms to the sensitivities output one by one.
01221       //
01222       // FIXME: this is probably a bad order of operations for
01223       // controlling floating point error.
01224 
01225       *parameters[k] = old_parameterk + delta_p;
01226       this->assembly(false, true);
01227       this->matrix->close();
01228       this->assemble_qoi_derivative(qoi_indices);
01229 
01230       for (unsigned int l=0; l != Np; ++l)
01231         {
01232           this->matrix->vector_mult(*tempvec, this->get_sensitivity_solution(l));
01233           for (unsigned int i=0; i != Nq; ++i)
01234             if (qoi_indices.has_index(i))
01235               {
01236                 this->get_adjoint_rhs(i).close();
01237                 Number current_terms =
01238                   (this->get_adjoint_rhs(i).dot(this->get_sensitivity_solution(l)) -
01239                    tempvec->dot(this->get_adjoint_solution(i))) / (2.*delta_p);
01240                 sensitivities.second_derivative(i,k,l) += current_terms;
01241 
01242                 // We use the _uk terms twice; symmetry lets us reuse
01243                 // these calculations for the _ul terms.
01244 
01245                 sensitivities.second_derivative(i,l,k) += current_terms;
01246               }
01247         }
01248 
01249       *parameters[k] = old_parameterk - delta_p;
01250       this->assembly(false, true);
01251       this->matrix->close();
01252       this->assemble_qoi_derivative(qoi_indices);
01253 
01254       for (unsigned int l=0; l != Np; ++l)
01255         {
01256           this->matrix->vector_mult(*tempvec, this->get_sensitivity_solution(l));
01257           for (unsigned int i=0; i != Nq; ++i)
01258             if (qoi_indices.has_index(i))
01259               {
01260                 this->get_adjoint_rhs(i).close();
01261                 Number current_terms =
01262                   (-this->get_adjoint_rhs(i).dot(this->get_sensitivity_solution(l)) +
01263                    tempvec->dot(this->get_adjoint_solution(i))) / (2.*delta_p);
01264                 sensitivities.second_derivative(i,k,l) += current_terms;
01265 
01266                 // We use the _uk terms twice; symmetry lets us reuse
01267                 // these calculations for the _ul terms.
01268 
01269                 sensitivities.second_derivative(i,l,k) += current_terms;
01270               }
01271         }
01272 
01273       // Don't leave the parameter perturbed
01274       *parameters[k] = old_parameterk;
01275 
01276       // Our last remaining terms are -R_uu(u,z)*u_k*u_l and
01277       // Q_uu(u)*u_k*u_l
01278       //
01279       // We take directional central finite differences of R_u and Q_u
01280       // to approximate these terms, e.g.:
01281       //
01282       // Q_uu(u)*u_k ~= (Q_u(u+dp*u_k) - Q_u(u-dp*u_k))/(2*dp)
01283 
01284       *this->solution = this->get_sensitivity_solution(k);
01285       *this->solution *= delta_p;
01286       *this->solution += *oldsolution;
01287       this->assembly(false, true);
01288       this->matrix->close();
01289       this->assemble_qoi_derivative(qoi_indices);
01290 
01291       // The Hessian is symmetric, so we just calculate the lower
01292       // triangle and the diagonal, and we get the upper triangle from
01293       // the transpose of the lower
01294       //
01295       // Note that, because we took the directional finite difference
01296       // with respect to k and not l, we've added an O(delta_p^2)
01297       // error to any permutational symmetry in the Hessian...
01298       for (unsigned int l=0; l != k+1; ++l)
01299         {
01300           this->matrix->vector_mult(*tempvec, this->get_sensitivity_solution(l));
01301           for (unsigned int i=0; i != Nq; ++i)
01302             if (qoi_indices.has_index(i))
01303               {
01304                 this->get_adjoint_rhs(i).close();
01305                 Number current_terms =
01306                   (this->get_adjoint_rhs(i).dot(this->get_sensitivity_solution(l)) -
01307                    tempvec->dot(this->get_adjoint_solution(i))) / (2.*delta_p);
01308                 sensitivities.second_derivative(i,k,l) += current_terms;
01309                 if (k != l)
01310                   sensitivities.second_derivative(i,l,k) += current_terms;
01311               }
01312         }
01313 
01314       *this->solution = this->get_sensitivity_solution(k);
01315       *this->solution *= -delta_p;
01316       *this->solution += *oldsolution;
01317       this->assembly(false, true);
01318       this->matrix->close();
01319       this->assemble_qoi_derivative(qoi_indices);
01320 
01321       for (unsigned int l=0; l != k+1; ++l)
01322         {
01323           this->matrix->vector_mult(*tempvec, this->get_sensitivity_solution(l));
01324           for (unsigned int i=0; i != Nq; ++i)
01325             if (qoi_indices.has_index(i))
01326               {
01327                 this->get_adjoint_rhs(i).close();
01328                 Number current_terms =
01329                   (-this->get_adjoint_rhs(i).dot(this->get_sensitivity_solution(l)) +
01330                    tempvec->dot(this->get_adjoint_solution(i))) / (2.*delta_p);
01331                 sensitivities.second_derivative(i,k,l) += current_terms;
01332                 if (k != l)
01333                   sensitivities.second_derivative(i,l,k) += current_terms;
01334               }
01335         }
01336 
01337       // Don't leave the solution perturbed
01338       *this->solution = *oldsolution;
01339     }
01340 
01341   // All parameters have been reset.
01342   // Don't leave the qoi or system changed - principle of least
01343   // surprise.
01344   this->assembly(true, true);
01345   this->rhs->close();
01346   this->matrix->close();
01347   this->assemble_qoi(qoi_indices);
01348 }

void libMesh::ImplicitSystem::qoi_parameter_hessian_vector_product ( const QoISet &  qoi_indices,
const ParameterVector &  parameters,
const ParameterVector &  vector,
SensitivityData &  product 
) [virtual, inherited]

For each of the system's quantities of interest q in qoi[qoi_indices], and for a vector of parameters p, the parameter sensitivity Hessian H_ij is defined as H_ij = (d^2 q)/(d p_i d p_j) The Hessian-vector product, for a vector v_k in parameter space, is S_j = H_jk v_k This product is the output of this method, where for each q_i, S_j is stored in sensitivities[i][j].

Reimplemented from libMesh::System.

Definition at line 880 of file implicit_system.C.

References libMesh::SensitivityData::allocate_data(), libMesh::ParameterVector::deep_copy(), libMesh::QoISet::has_index(), libMesh::Real, libMesh::ParameterVector::size(), libMesh::TOLERANCE, and libMesh::ParameterVector::value_copy().

00884 {
00885   // We currently get partial derivatives via finite differencing
00886   const Real delta_p = TOLERANCE;
00887 
00888   // We'll use a single temporary vector for matrix-vector-vector products
00889   AutoPtr<NumericVector<Number> > tempvec = this->solution->zero_clone();
00890 
00891   const unsigned int Np = libmesh_cast_int<unsigned int>
00892     (parameters.size());
00893   const unsigned int Nq = libmesh_cast_int<unsigned int>
00894     (qoi.size());
00895 
00896   // For each quantity of interest q, the parameter sensitivity
00897   // Hessian is defined as q''_{kl} = {d^2 q}/{d p_k d p_l}.
00898   // Given a vector of parameter perturbation weights w_l, this
00899   // function evaluates the hessian-vector product sum_l(q''_{kl}*w_l)
00900   //
00901   // We calculate it from values and partial derivatives of the
00902   // quantity of interest function Q, solution u, adjoint solution z,
00903   // parameter sensitivity adjoint solutions z^l, and residual R, as:
00904   //
00905   // sum_l(q''_{kl}*w_l) =
00906   // sum_l(w_l * Q''_{kl}) + Q''_{uk}(u)*(sum_l(w_l u'_l)) -
00907   // R'_k(u, sum_l(w_l*z^l)) - R'_{uk}(u,z)*(sum_l(w_l u'_l) -
00908   // sum_l(w_l*R''_{kl}(u,z))
00909   //
00910   // See the adjoints model document for more details.
00911 
00912   // We first do an adjoint solve to get z for each quantity of
00913   // interest
00914   // if we havent already or dont have an initial condition for the adjoint
00915   if (!this->is_adjoint_already_solved())
00916     {
00917       this->adjoint_solve(qoi_indices);
00918     }
00919 
00920   // Get ready to fill in senstivities:
00921   sensitivities.allocate_data(qoi_indices, *this, parameters);
00922 
00923   // We can't solve for all the solution sensitivities u'_l or for all
00924   // of the parameter sensitivity adjoint solutions z^l without
00925   // requiring O(Nq*Np) linear solves.  So we'll solve directly for their
00926   // weighted sum - this is just O(Nq) solves.
00927 
00928   // First solve for sum_l(w_l u'_l).
00929   this->weighted_sensitivity_solve(parameters, vector);
00930 
00931   // Then solve for sum_l(w_l z^l).
00932   this->weighted_sensitivity_adjoint_solve(parameters, vector, qoi_indices);
00933 
00934   for (unsigned int k=0; k != Np; ++k)
00935     {
00936       // We approximate sum_l(w_l * Q''_{kl}) with a central
00937       // differencing pertubation:
00938       // sum_l(w_l * Q''_{kl}) ~=
00939       // (Q(p + dp*w_l*e_l + dp*e_k) - Q(p - dp*w_l*e_l + dp*e_k) -
00940       // Q(p + dp*w_l*e_l - dp*e_k) + Q(p - dp*w_l*e_l - dp*e_k))/(4*dp^2)
00941 
00942       // The sum(w_l*R''_kl) term requires the same sort of pertubation,
00943       // and so we subtract it in at the same time:
00944       // sum_l(w_l * R''_{kl}) ~=
00945       // (R(p + dp*w_l*e_l + dp*e_k) - R(p - dp*w_l*e_l + dp*e_k) -
00946       // R(p + dp*w_l*e_l - dp*e_k) + R(p - dp*w_l*e_l - dp*e_k))/(4*dp^2)
00947 
00948       ParameterVector oldparameters, parameterperturbation;
00949       parameters.deep_copy(oldparameters);
00950       vector.deep_copy(parameterperturbation);
00951       parameterperturbation *= delta_p;
00952       parameters += parameterperturbation;
00953 
00954       Number old_parameter = *parameters[k];
00955 
00956       *parameters[k] = old_parameter + delta_p;
00957       this->assemble_qoi(qoi_indices);
00958       this->assembly(true, false);
00959       this->rhs->close();
00960       std::vector<Number> partial2q_term = this->qoi;
00961       std::vector<Number> partial2R_term(this->qoi.size());
00962       for (unsigned int i=0; i != Nq; ++i)
00963         if (qoi_indices.has_index(i))
00964           partial2R_term[i] = this->rhs->dot(this->get_adjoint_solution(i));
00965 
00966       *parameters[k] = old_parameter - delta_p;
00967       this->assemble_qoi(qoi_indices);
00968       this->assembly(true, false);
00969       this->rhs->close();
00970       for (unsigned int i=0; i != Nq; ++i)
00971         if (qoi_indices.has_index(i))
00972           {
00973             partial2q_term[i] -= this->qoi[i];
00974             partial2R_term[i] -= this->rhs->dot(this->get_adjoint_solution(i));
00975           }
00976 
00977       oldparameters.value_copy(parameters);
00978       parameterperturbation *= -1.0;
00979       parameters += parameterperturbation;
00980 
00981       // Re-center old_parameter, which may be affected by vector
00982       old_parameter = *parameters[k];
00983 
00984       *parameters[k] = old_parameter + delta_p;
00985       this->assemble_qoi(qoi_indices);
00986       this->assembly(true, false);
00987       this->rhs->close();
00988       for (unsigned int i=0; i != Nq; ++i)
00989         if (qoi_indices.has_index(i))
00990           {
00991             partial2q_term[i] -= this->qoi[i];
00992             partial2R_term[i] -= this->rhs->dot(this->get_adjoint_solution(i));
00993           }
00994 
00995       *parameters[k] = old_parameter - delta_p;
00996       this->assemble_qoi(qoi_indices);
00997       this->assembly(true, false);
00998       this->rhs->close();
00999       for (unsigned int i=0; i != Nq; ++i)
01000         if (qoi_indices.has_index(i))
01001           {
01002             partial2q_term[i] += this->qoi[i];
01003             partial2R_term[i] += this->rhs->dot(this->get_adjoint_solution(i));
01004           }
01005 
01006       for (unsigned int i=0; i != Nq; ++i)
01007         if (qoi_indices.has_index(i))
01008           {
01009             partial2q_term[i] /= (4. * delta_p * delta_p);
01010             partial2R_term[i] /= (4. * delta_p * delta_p);
01011           }
01012 
01013       for (unsigned int i=0; i != Nq; ++i)
01014         if (qoi_indices.has_index(i))
01015           sensitivities[i][k] = partial2q_term[i] - partial2R_term[i];
01016 
01017       // We get (partial q / partial u), R, and
01018       // (partial R / partial u) from the user, but centrally
01019       // difference to get q_uk, R_k, and R_uk terms:
01020       // (partial R / partial k)
01021       // R_k*sum(w_l*z^l) = (R(p+dp*e_k)*sum(w_l*z^l) - R(p-dp*e_k)*sum(w_l*z^l))/(2*dp)
01022       // (partial^2 q / partial u partial k)
01023       // q_uk = (q_u(p+dp*e_k) - q_u(p-dp*e_k))/(2*dp)
01024       // (partial^2 R / partial u partial k)
01025       // R_uk*z*sum(w_l*u'_l) = (R_u(p+dp*e_k)*z*sum(w_l*u'_l) - R_u(p-dp*e_k)*z*sum(w_l*u'_l))/(2*dp)
01026 
01027       // To avoid creating Nq temporary vectors for q_uk or R_uk, we add
01028       // subterms to the sensitivities output one by one.
01029       //
01030       // FIXME: this is probably a bad order of operations for
01031       // controlling floating point error.
01032 
01033       *parameters[k] = old_parameter + delta_p;
01034       this->assembly(true, true);
01035       this->rhs->close();
01036       this->matrix->close();
01037       this->assemble_qoi_derivative(qoi_indices);
01038 
01039       this->matrix->vector_mult(*tempvec, this->get_weighted_sensitivity_solution());
01040 
01041       for (unsigned int i=0; i != Nq; ++i)
01042         if (qoi_indices.has_index(i))
01043           {
01044             this->get_adjoint_rhs(i).close();
01045             sensitivities[i][k] += (this->get_adjoint_rhs(i).dot(this->get_weighted_sensitivity_solution()) -
01046                                     this->rhs->dot(this->get_weighted_sensitivity_adjoint_solution(i)) -
01047                                     this->get_adjoint_solution(i).dot(*tempvec)) / (2.*delta_p);
01048           }
01049 
01050       *parameters[k] = old_parameter - delta_p;
01051       this->assembly(true, true);
01052       this->rhs->close();
01053       this->matrix->close();
01054       this->assemble_qoi_derivative(qoi_indices);
01055 
01056       this->matrix->vector_mult(*tempvec, this->get_weighted_sensitivity_solution());
01057 
01058       for (unsigned int i=0; i != Nq; ++i)
01059         if (qoi_indices.has_index(i))
01060           {
01061             this->get_adjoint_rhs(i).close();
01062             sensitivities[i][k] += (-this->get_adjoint_rhs(i).dot(this->get_weighted_sensitivity_solution()) +
01063                                     this->rhs->dot(this->get_weighted_sensitivity_adjoint_solution(i)) +
01064                                     this->get_adjoint_solution(i).dot(*tempvec)) / (2.*delta_p);
01065           }
01066     }
01067 
01068   // All parameters have been reset.
01069   // Don't leave the qoi or system changed - principle of least
01070   // surprise.
01071   this->assembly(true, true);
01072   this->rhs->close();
01073   this->matrix->close();
01074   this->assemble_qoi(qoi_indices);
01075 }

void libMesh::System::qoi_parameter_sensitivity ( const QoISet &  qoi_indices,
const ParameterVector &  parameters,
SensitivityData &  sensitivities 
) [virtual, inherited]

Solves for the derivative of each of the system's quantities of interest q in qoi[qoi_indices] with respect to each parameter in parameters, placing the result for qoi i and parameter j into sensitivities[i][j].

Note that parameters is a const vector, not a vector-of-const; parameter values in this vector need to be mutable for finite differencing to work.

Automatically chooses the forward method for problems with more quantities of interest than parameters, or the adjoint method otherwise.

This method is only usable in derived classes which overload an implementation.

Definition at line 509 of file system.C.

References libMesh::ParameterVector::size(), and libMesh::QoISet::size().

00512 {
00513   // Forward sensitivities are more efficient for Nq > Np
00514   if (qoi_indices.size(*this) > parameters.size())
00515     forward_qoi_parameter_sensitivity(qoi_indices, parameters, sensitivities);
00516   // Adjoint sensitivities are more efficient for Np > Nq,
00517   // and an adjoint may be more reusable than a forward
00518   // solution sensitivity in the Np == Nq case.
00519   else
00520     adjoint_qoi_parameter_sensitivity(qoi_indices, parameters, sensitivities);
00521 }

void libMesh::System::re_update (  )  [virtual, inherited]

Re-update the local values when the mesh has changed. This method takes the data updated by update() and makes it up-to-date on the current mesh.

Reimplemented in libMesh::TransientSystem< RBConstruction >.

Definition at line 431 of file system.C.

References libMesh::System::current_local_solution, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::n_vars(), and libMesh::System::solution.

00432 {
00433   parallel_only();
00434 
00435   // If this system is empty... don't do anything!
00436   if(!this->n_vars())
00437     return;
00438 
00439   const std::vector<dof_id_type>& send_list = this->get_dof_map().get_send_list ();
00440 
00441   // Check sizes
00442   libmesh_assert_equal_to (current_local_solution->size(), solution->size());
00443   // Not true with ghosted vectors
00444   // libmesh_assert_equal_to (current_local_solution->local_size(), solution->size());
00445   // libmesh_assert (!send_list.empty());
00446   libmesh_assert_less_equal (send_list.size(), solution->size());
00447 
00448   // Create current_local_solution from solution.  This will
00449   // put a local copy of solution into current_local_solution.
00450   solution->localize (*current_local_solution, send_list);
00451 }

void libMesh::System::read_header ( Xdr &  io,
const std::string &  version,
const bool  read_header = true,
const bool  read_additional_data = true,
const bool  read_legacy_format = false 
) [inherited]

Reads the basic data header for this System.

Definition at line 113 of file system_io.C.

References libMesh::System::_additional_data_written, libMesh::System::_written_var_indices, libMesh::System::add_variable(), libMesh::System::add_vector(), libMesh::Parallel::Communicator::broadcast(), libMesh::System::clear(), libMesh::CommWorld, libMesh::Xdr::data(), libMesh::FEType::family, libMesh::System::get_mesh(), libMesh::FEType::inf_map, libMesh::MeshBase::mesh_dimension(), libMeshEnums::MONOMIAL, libMesh::on_command_line(), libMesh::FEType::order, libMesh::out, libMesh::processor_id(), libMesh::FEType::radial_family, libMesh::FEType::radial_order, libMesh::Xdr::reading(), libMesh::System::variable_number(), libMesh::Xdr::version(), and libMeshEnums::XYZ.

Referenced by libMesh::EquationSystems::_read_impl().

00118 {
00119   // This method implements the input of a
00120   // System object, embedded in the output of
00121   // an EquationSystems<T_sys>.  This warrants some
00122   // documentation.  The output file essentially
00123   // consists of 5 sections:
00124   //
00125   // for this system
00126   //
00127   //   5.) The number of variables in the system (unsigned int)
00128   //
00129   //   for each variable in the system
00130   //
00131   //     6.) The name of the variable (string)
00132   //
00133   //     6.1.) Variable subdmains
00134   //
00135   //     7.) Combined in an FEType:
00136   //         - The approximation order(s) of the variable
00137   //           (Order Enum, cast to int/s)
00138   //         - The finite element family/ies of the variable
00139   //           (FEFamily Enum, cast to int/s)
00140   //
00141   //   end variable loop
00142   //
00143   //   8.) The number of additional vectors (unsigned int),
00144   //
00145   //     for each additional vector in the system object
00146   //
00147   //     9.) the name of the additional vector  (string)
00148   //
00149   // end system
00150   libmesh_assert (io.reading());
00151 
00152   // Possibly clear data structures and start from scratch.
00153   if (read_header_in)
00154     this->clear ();
00155 
00156   // Figure out if we need to read infinite element information.
00157   // This will be true if the version string contains " with infinite elements"
00158   const bool read_ifem_info =
00159     (version.rfind(" with infinite elements") < version.size()) ||
00160     libMesh::on_command_line ("--read_ifem_systems");
00161 
00162 
00163   {
00164     // 5.)
00165     // Read the number of variables in the system
00166     unsigned int nv=0;
00167     if (libMesh::processor_id() == 0)
00168       io.data (nv);
00169     CommWorld.broadcast(nv);
00170 
00171     _written_var_indices.clear();
00172     _written_var_indices.resize(nv, 0);
00173 
00174     for (unsigned int var=0; var<nv; var++)
00175       {
00176         // 6.)
00177         // Read the name of the var-th variable
00178         std::string var_name;
00179         if (libMesh::processor_id() == 0)
00180           io.data (var_name);
00181         CommWorld.broadcast(var_name);
00182 
00183         // 6.1.)
00184         std::set<subdomain_id_type> domains;
00185         if (io.version() >= LIBMESH_VERSION_ID(0,7,2))
00186         {
00187           std::vector<subdomain_id_type> domain_array;
00188           if (libMesh::processor_id() == 0)
00189             io.data (domain_array);
00190           for (std::vector<subdomain_id_type>::iterator it = domain_array.begin(); it != domain_array.end(); ++it)
00191                   domains.insert(*it);
00192         }
00193         CommWorld.broadcast(domains);
00194 
00195         // 7.)
00196         // Read the approximation order(s) of the var-th variable
00197         int order=0;
00198         if (libMesh::processor_id() == 0)
00199           io.data (order);
00200         CommWorld.broadcast(order);
00201 
00202 
00203         // do the same for infinite element radial_order
00204         int rad_order=0;
00205         if (read_ifem_info)
00206           {
00207             if (libMesh::processor_id() == 0)
00208               io.data(rad_order);
00209             CommWorld.broadcast(rad_order);
00210           }
00211 
00212 
00213         // Read the finite element type of the var-th variable
00214         int fam=0;
00215         if (libMesh::processor_id() == 0)
00216           io.data (fam);
00217         CommWorld.broadcast(fam);
00218         FEType type;
00219         type.order  = static_cast<Order>(order);
00220         type.family = static_cast<FEFamily>(fam);
00221 
00222         // Check for incompatibilities.  The shape function indexing was
00223         // changed for the monomial and xyz finite element families to
00224         // simplify extension to arbitrary p.  The consequence is that
00225         // old restart files will not be read correctly.  This is expected
00226         // to be an unlikely occurance, but catch it anyway.
00227         if (read_legacy_format)
00228           if ((type.family == MONOMIAL || type.family == XYZ) &&
00229               ((type.order > 2 && this->get_mesh().mesh_dimension() == 2) ||
00230                (type.order > 1 && this->get_mesh().mesh_dimension() == 3)))
00231             {
00232               libmesh_here();
00233               libMesh::out << "*****************************************************************\n"
00234                             << "* WARNING: reading a potentially incompatible restart file!!!   *\n"
00235                             << "*  contact libmesh-users@lists.sourceforge.net for more details *\n"
00236                             << "*****************************************************************"
00237                             << std::endl;
00238             }
00239 
00240         // Read additional information for infinite elements
00241         int radial_fam=0;
00242         int i_map=0;
00243         if (read_ifem_info)
00244           {
00245             if (libMesh::processor_id() == 0)
00246               io.data (radial_fam);
00247             CommWorld.broadcast(radial_fam);
00248             if (libMesh::processor_id() == 0)
00249               io.data (i_map);
00250             CommWorld.broadcast(i_map);
00251           }
00252 
00253 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
00254 
00255         type.radial_order  = static_cast<Order>(rad_order);
00256         type.radial_family = static_cast<FEFamily>(radial_fam);
00257         type.inf_map       = static_cast<InfMapType>(i_map);
00258 
00259 #endif
00260 
00261         if (read_header_in)
00262         {
00263           if (domains.empty())
00264             _written_var_indices[var] = this->add_variable (var_name, type);
00265           else
00266             _written_var_indices[var] = this->add_variable (var_name, type, &domains);
00267         }
00268     else
00269           _written_var_indices[var] = this->variable_number(var_name);
00270       }
00271   }
00272 
00273   // 8.)
00274   // Read the number of additional vectors.
00275   unsigned int nvecs=0;
00276   if (libMesh::processor_id() == 0)
00277     io.data (nvecs);
00278   CommWorld.broadcast(nvecs);
00279 
00280   // If nvecs > 0, this means that write_additional_data
00281   // was true when this file was written.  We will need to
00282   // make use of this fact later.
00283   if (nvecs > 0)
00284     this->_additional_data_written = true;
00285 
00286   for (unsigned int vec=0; vec<nvecs; vec++)
00287     {
00288       // 9.)
00289       // Read the name of the vec-th additional vector
00290       std::string vec_name;
00291       if (libMesh::processor_id() == 0)
00292         io.data (vec_name);
00293       CommWorld.broadcast(vec_name);
00294 
00295       if (read_additional_data)
00296         {
00297           // Systems now can handle adding post-initialization vectors
00298 //        libmesh_assert(this->_can_add_vectors);
00299           // Some systems may have added their own vectors already
00300 //        libmesh_assert_equal_to (this->_vectors.count(vec_name), 0);
00301 
00302           this->add_vector(vec_name);
00303         }
00304     }
00305 }

void libMesh::System::read_legacy_data ( Xdr &  io,
const bool  read_additional_data = true 
) [inherited]

Reads additional data, namely vectors, for this System.

Definition at line 309 of file system_io.C.

References libMesh::System::_additional_data_written, libMesh::System::_vectors, libMesh::System::_written_var_indices, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::Parallel::Communicator::broadcast(), libMesh::CommWorld, libMesh::Xdr::data(), end, libMesh::System::get_mesh(), libMesh::DofObject::invalid_id, libMesh::System::n_dofs(), libMesh::System::n_vars(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), libMesh::System::number(), libMesh::processor_id(), libMesh::Xdr::reading(), libMesh::System::solution, and libMesh::zero.

00311 {
00312   libmesh_deprecated();
00313 
00314   // This method implements the output of the vectors
00315   // contained in this System object, embedded in the
00316   // output of an EquationSystems<T_sys>.
00317   //
00318   //   10.) The global solution vector, re-ordered to be node-major
00319   //       (More on this later.)
00320   //
00321   //      for each additional vector in the object
00322   //
00323   //      11.) The global additional vector, re-ordered to be
00324   //           node-major (More on this later.)
00325   libmesh_assert (io.reading());
00326 
00327   // read and reordering buffers
00328   std::vector<Number> global_vector;
00329   std::vector<Number> reordered_vector;
00330 
00331   // 10.)
00332   // Read and set the solution vector
00333   {
00334     if (libMesh::processor_id() == 0)
00335       io.data (global_vector);
00336     CommWorld.broadcast(global_vector);
00337 
00338     // Remember that the stored vector is node-major.
00339     // We need to put it into whatever application-specific
00340     // ordering we may have using the dof_map.
00341     reordered_vector.resize(global_vector.size());
00342 
00343     //libMesh::out << "global_vector.size()=" << global_vector.size() << std::endl;
00344     //libMesh::out << "this->n_dofs()=" << this->n_dofs() << std::endl;
00345 
00346     libmesh_assert_equal_to (global_vector.size(), this->n_dofs());
00347 
00348     dof_id_type cnt=0;
00349 
00350     const unsigned int sys = this->number();
00351     const unsigned int nv  = this->_written_var_indices.size();
00352     libmesh_assert_less_equal (nv, this->n_vars());
00353 
00354     for (unsigned int data_var=0; data_var<nv; data_var++)
00355       {
00356         const unsigned int var = _written_var_indices[data_var];
00357 
00358         // First reorder the nodal DOF values
00359         {
00360           MeshBase::node_iterator
00361             it  = this->get_mesh().nodes_begin(),
00362             end = this->get_mesh().nodes_end();
00363 
00364           for (; it != end; ++it)
00365             for (unsigned int index=0; index<(*it)->n_comp(sys,var); index++)
00366               {
00367                 libmesh_assert_not_equal_to ((*it)->dof_number(sys, var, index),
00368                                              DofObject::invalid_id);
00369 
00370                 libmesh_assert_less (cnt, global_vector.size());
00371 
00372                 reordered_vector[(*it)->dof_number(sys, var, index)] =
00373                   global_vector[cnt++];
00374               }
00375         }
00376 
00377         // Then reorder the element DOF values
00378         {
00379           MeshBase::element_iterator
00380             it  = this->get_mesh().active_elements_begin(),
00381             end = this->get_mesh().active_elements_end();
00382 
00383           for (; it != end; ++it)
00384             for (unsigned int index=0; index<(*it)->n_comp(sys,var); index++)
00385               {
00386                 libmesh_assert_not_equal_to ((*it)->dof_number(sys, var, index),
00387                                              DofObject::invalid_id);
00388 
00389                 libmesh_assert_less (cnt, global_vector.size());
00390 
00391                 reordered_vector[(*it)->dof_number(sys, var, index)] =
00392                   global_vector[cnt++];
00393               }
00394         }
00395       }
00396 
00397     *(this->solution) = reordered_vector;
00398   }
00399 
00400   // For each additional vector, simply go through the list.
00401   // ONLY attempt to do this IF additional data was actually
00402   // written to the file for this system (controlled by the
00403   // _additional_data_written flag).
00404   if (this->_additional_data_written)
00405     {
00406       std::map<std::string, NumericVector<Number>* >::iterator
00407         pos = this->_vectors.begin();
00408 
00409       for (; pos != this->_vectors.end(); ++pos)
00410         {
00411           // 11.)
00412           // Read the values of the vec-th additional vector.
00413           // Prior do _not_ clear, but fill with zero, since the
00414           // additional vectors _have_ to have the same size
00415           // as the solution vector
00416           std::fill (global_vector.begin(), global_vector.end(), libMesh::zero);
00417 
00418           if (libMesh::processor_id() == 0)
00419             io.data (global_vector);
00420           CommWorld.broadcast(global_vector);
00421 
00422           // If read_additional_data==true, then we will keep this vector, otherwise
00423           // we are going to throw it away.
00424           if (read_additional_data)
00425             {
00426               // Remember that the stored vector is node-major.
00427               // We need to put it into whatever application-specific
00428               // ordering we may have using the dof_map.
00429               std::fill (reordered_vector.begin(),
00430                          reordered_vector.end(),
00431                          libMesh::zero);
00432 
00433               reordered_vector.resize(global_vector.size());
00434 
00435               libmesh_assert_equal_to (global_vector.size(), this->n_dofs());
00436 
00437               dof_id_type cnt=0;
00438 
00439               const unsigned int sys = this->number();
00440               const unsigned int nv  = this->_written_var_indices.size();
00441               libmesh_assert_less_equal (nv, this->n_vars());
00442 
00443               for (unsigned int data_var=0; data_var<nv; data_var++)
00444                 {
00445                   const unsigned int var = _written_var_indices[data_var];
00446                   // First reorder the nodal DOF values
00447                   {
00448                     MeshBase::node_iterator
00449                       it  = this->get_mesh().nodes_begin(),
00450                       end = this->get_mesh().nodes_end();
00451 
00452                     for (; it!=end; ++it)
00453                       for (unsigned int index=0; index<(*it)->n_comp(sys,var); index++)
00454                         {
00455                           libmesh_assert_not_equal_to ((*it)->dof_number(sys, var, index),
00456                                                        DofObject::invalid_id);
00457 
00458                           libmesh_assert_less (cnt, global_vector.size());
00459 
00460                           reordered_vector[(*it)->dof_number(sys, var, index)] =
00461                             global_vector[cnt++];
00462                         }
00463                   }
00464 
00465                   // Then reorder the element DOF values
00466                   {
00467                     MeshBase::element_iterator
00468                       it  = this->get_mesh().active_elements_begin(),
00469                       end = this->get_mesh().active_elements_end();
00470 
00471                     for (; it!=end; ++it)
00472                       for (unsigned int index=0; index<(*it)->n_comp(sys,var); index++)
00473                         {
00474                           libmesh_assert_not_equal_to ((*it)->dof_number(sys, var, index),
00475                                                        DofObject::invalid_id);
00476 
00477                           libmesh_assert_less (cnt, global_vector.size());
00478 
00479                           reordered_vector[(*it)->dof_number(sys, var, index)] =
00480                             global_vector[cnt++];
00481                         }
00482                   }
00483                 }
00484 
00485               // use the overloaded operator=(std::vector) to assign the values
00486               *(pos->second) = reordered_vector;
00487             }
00488         }
00489     } // end if (_additional_data_written)
00490 }

void libMesh::System::read_parallel_data ( Xdr &  io,
const bool  read_additional_data 
) [inherited]

Reads additional data, namely vectors, for this System. This method may safely be called on a distributed-memory mesh. This method will read an individual file for each processor in the simulation where the local solution components for that processor are stored.

This method implements the output of the vectors contained in this System object, embedded in the output of an EquationSystems<T_sys>.

9.) The global solution vector, re-ordered to be node-major (More on this later.)

for each additional vector in the object

10.) The global additional vector, re-ordered to be node-major (More on this later.)

Note that the actual IO is handled through the Xdr class (to be renamed later?) which provides a uniform interface to both the XDR (eXternal Data Representation) interface and standard ASCII output. Thus this one section of code will read XDR or ASCII files with no changes.

Definition at line 494 of file system_io.C.

References libMesh::System::_vectors, libMesh::System::_written_var_indices, libMesh::Xdr::data(), libMesh::FEType::family, libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::DofObject::invalid_id, libMesh::Xdr::is_open(), libMesh::n_processors(), libMesh::System::n_vars(), libMesh::System::number(), libMesh::processor_id(), libMesh::Xdr::reading(), libMeshEnums::SCALAR, libMesh::DofMap::SCALAR_dof_indices(), libMesh::System::solution, libMesh::Variable::type(), and libMesh::System::variable().

00496 {
00516   // PerfLog pl("IO Performance",false);
00517   // pl.push("read_parallel_data");
00518   dof_id_type total_read_size = 0;
00519   
00520   libmesh_assert (io.reading());
00521   libmesh_assert (io.is_open());
00522 
00523   // build the ordered nodes and element maps.
00524   // when writing/reading parallel files we need to iterate
00525   // over our nodes/elements in order of increasing global id().
00526   // however, this is not guaranteed to be ordering we obtain
00527   // by using the node_iterators/element_iterators directly.
00528   // so build a set, sorted by id(), that provides the ordering.
00529   // further, for memory economy build the set but then transfer
00530   // its contents to vectors, which will be sorted.
00531   std::vector<const DofObject*> ordered_nodes, ordered_elements;
00532   {
00533     std::set<const DofObject*, CompareDofObjectsByID>
00534       ordered_nodes_set (this->get_mesh().local_nodes_begin(),
00535                          this->get_mesh().local_nodes_end());
00536 
00537       ordered_nodes.insert(ordered_nodes.end(),
00538                            ordered_nodes_set.begin(),
00539                            ordered_nodes_set.end());
00540   }
00541   {
00542     std::set<const DofObject*, CompareDofObjectsByID>
00543       ordered_elements_set (this->get_mesh().local_elements_begin(),
00544                             this->get_mesh().local_elements_end());
00545 
00546       ordered_elements.insert(ordered_elements.end(),
00547                               ordered_elements_set.begin(),
00548                               ordered_elements_set.end());
00549   }
00550 
00551   std::vector<Number> io_buffer;
00552 
00553   // 9.)
00554   //
00555   // Actually read the solution components
00556   // for the ith system to disk
00557   io.data(io_buffer);
00558 
00559   total_read_size += io_buffer.size();
00560 
00561   const unsigned int sys_num = this->number();
00562   const unsigned int nv      = this->_written_var_indices.size();
00563   libmesh_assert_less_equal (nv, this->n_vars());
00564 
00565   dof_id_type cnt=0;
00566 
00567   // Loop over each non-SCALAR variable and each node, and read out the value.
00568   for (unsigned int data_var=0; data_var<nv; data_var++)
00569     {
00570       const unsigned int var = _written_var_indices[data_var];
00571       if(this->variable(var).type().family != SCALAR)
00572         {
00573           // First read the node DOF values
00574           for (std::vector<const DofObject*>::const_iterator
00575                it = ordered_nodes.begin(); it != ordered_nodes.end(); ++it)
00576             for (unsigned int comp=0; comp<(*it)->n_comp(sys_num, var); comp++)
00577               {
00578                 libmesh_assert_not_equal_to ((*it)->dof_number(sys_num, var, comp),
00579                                              DofObject::invalid_id);
00580                 libmesh_assert_less (cnt, io_buffer.size());
00581                 this->solution->set((*it)->dof_number(sys_num, var, comp), io_buffer[cnt++]);
00582               }
00583 
00584           // Then read the element DOF values
00585           for (std::vector<const DofObject*>::const_iterator
00586              it = ordered_elements.begin(); it != ordered_elements.end(); ++it)
00587             for (unsigned int comp=0; comp<(*it)->n_comp(sys_num, var); comp++)
00588             {
00589               libmesh_assert_not_equal_to ((*it)->dof_number(sys_num, var, comp),
00590                                            DofObject::invalid_id);
00591               libmesh_assert_less (cnt, io_buffer.size());
00592               this->solution->set((*it)->dof_number(sys_num, var, comp), io_buffer[cnt++]);
00593             }
00594         }
00595     }
00596 
00597   // Finally, read the SCALAR variables on the last processor
00598   for (unsigned int data_var=0; data_var<nv; data_var++)
00599     {
00600       const unsigned int var = _written_var_indices[data_var];
00601       if(this->variable(var).type().family == SCALAR)
00602         {
00603           if (libMesh::processor_id() == (libMesh::n_processors()-1))
00604             {
00605               const DofMap& dof_map = this->get_dof_map();
00606               std::vector<dof_id_type> SCALAR_dofs;
00607               dof_map.SCALAR_dof_indices(SCALAR_dofs, var);
00608 
00609               for(unsigned int i=0; i<SCALAR_dofs.size(); i++)
00610                 {
00611                   this->solution->set( SCALAR_dofs[i], io_buffer[cnt++] );
00612                 }
00613             }
00614         }
00615     }
00616 
00617   // And we're done setting solution entries
00618   this->solution->close();
00619 
00620   // Only read additional vectors if wanted
00621   if (read_additional_data)
00622     {
00623       std::map<std::string, NumericVector<Number>* >::const_iterator
00624         pos = _vectors.begin();
00625 
00626       for(; pos != this->_vectors.end(); ++pos)
00627         {
00628           cnt=0;
00629           io_buffer.clear();
00630 
00631           // 10.)
00632           //
00633           // Actually read the additional vector components
00634           // for the ith system to disk
00635           io.data(io_buffer);
00636 
00637           total_read_size += io_buffer.size();
00638 
00639           // Loop over each non-SCALAR variable and each node, and read out the value.
00640           for (unsigned int data_var=0; data_var<nv; data_var++)
00641             {
00642               const unsigned int var = _written_var_indices[data_var];
00643               if(this->variable(var).type().family != SCALAR)
00644                 {
00645                   // First read the node DOF values
00646                   for (std::vector<const DofObject*>::const_iterator
00647                        it = ordered_nodes.begin(); it != ordered_nodes.end(); ++it)
00648                     for (unsigned int comp=0; comp<(*it)->n_comp(sys_num, var); comp++)
00649                       {
00650                         libmesh_assert_not_equal_to ((*it)->dof_number(sys_num, var, comp),
00651                                                      DofObject::invalid_id);
00652                         libmesh_assert_less (cnt, io_buffer.size());
00653                         pos->second->set((*it)->dof_number(sys_num, var, comp), io_buffer[cnt++]);
00654                       }
00655 
00656                   // Then read the element DOF values
00657                   for (std::vector<const DofObject*>::const_iterator
00658                        it = ordered_elements.begin(); it != ordered_elements.end(); ++it)
00659                     for (unsigned int comp=0; comp<(*it)->n_comp(sys_num, var); comp++)
00660                       {
00661                         libmesh_assert_not_equal_to ((*it)->dof_number(sys_num, var, comp),
00662                                                      DofObject::invalid_id);
00663                         libmesh_assert_less (cnt, io_buffer.size());
00664                         pos->second->set((*it)->dof_number(sys_num, var, comp), io_buffer[cnt++]);
00665                       }
00666                 }
00667             }
00668 
00669           // Finally, read the SCALAR variables on the last processor
00670           for (unsigned int data_var=0; data_var<nv; data_var++)
00671             {
00672               const unsigned int var = _written_var_indices[data_var];
00673               if(this->variable(var).type().family == SCALAR)
00674                 {
00675                   if (libMesh::processor_id() == (libMesh::n_processors()-1))
00676                     {
00677                       const DofMap& dof_map = this->get_dof_map();
00678                       std::vector<dof_id_type> SCALAR_dofs;
00679                       dof_map.SCALAR_dof_indices(SCALAR_dofs, var);
00680 
00681                       for(unsigned int i=0; i<SCALAR_dofs.size(); i++)
00682                         {
00683                           pos->second->set( SCALAR_dofs[i], io_buffer[cnt++] );
00684                         }
00685                     }
00686                 }
00687             }
00688 
00689           // And we're done setting entries for this variable
00690           pos->second->close();
00691         }
00692     }
00693 
00694   // const Real 
00695   //   dt   = pl.get_elapsed_time(),
00696   //   rate = total_read_size*sizeof(Number)/dt; 
00697 
00698   // std::cerr << "Read " << total_read_size << " \"Number\" values\n"
00699   //        << " Elapsed time = " << dt << '\n'
00700   //        << " Rate = " << rate/1.e6 << "(MB/sec)\n\n";
00701 
00702   // pl.pop("read_parallel_data");
00703 }

void libMesh::RBParametrized::read_parameter_ranges_from_file ( const std::string &  file_name,
const bool  read_binary 
) [inherited]

Read in the parameter ranges from file. Initialize parameters to the "minimum" parameter values.

virtual void libMesh::RBConstruction::read_riesz_representors_from_files ( const std::string &  riesz_representors_dir,
const bool  write_binary_residual_representors 
) [virtual, inherited]

Read in all the Riesz representor data from files. directory_name specifies which directory to read from. io_flag specifies whether we read in all data, or only a basis (in)dependent subset.

Reimplemented in libMesh::TransientRBConstruction.

void libMesh::System::read_serialized_data ( Xdr &  io,
const bool  read_additional_data = true 
) [inherited]

Reads additional data, namely vectors, for this System. This method may safely be called on a distributed-memory mesh.

Definition at line 707 of file system_io.C.

References libMesh::System::_vectors, libMesh::Xdr::comment(), libMesh::processor_id(), libMesh::System::read_serialized_vector(), and libMesh::System::solution.

00709 {
00710   // This method implements the input of the vectors
00711   // contained in this System object, embedded in the
00712   // output of an EquationSystems<T_sys>.
00713   //
00714   //   10.) The global solution vector, re-ordered to be node-major
00715   //       (More on this later.)
00716   //
00717   //      for each additional vector in the object
00718   //
00719   //      11.) The global additional vector, re-ordered to be
00720   //          node-major (More on this later.)
00721   parallel_only();
00722   std::string comment;
00723 
00724   // PerfLog pl("IO Performance",false);
00725   // pl.push("read_serialized_data");
00726   std::size_t total_read_size = 0;
00727 
00728   // 10.)
00729   // Read the global solution vector
00730   {
00731     total_read_size +=
00732       this->read_serialized_vector(io, *this->solution);
00733 
00734     // get the comment
00735     if (libMesh::processor_id() == 0)
00736       io.comment (comment);
00737   }
00738 
00739   // 11.)
00740   // Only read additional vectors if wanted
00741   if (read_additional_data)
00742     {
00743       std::map<std::string, NumericVector<Number>* >::const_iterator
00744         pos = _vectors.begin();
00745 
00746       for(; pos != this->_vectors.end(); ++pos)
00747         {
00748           total_read_size +=
00749             this->read_serialized_vector(io, *pos->second);
00750 
00751           // get the comment
00752           if (libMesh::processor_id() == 0)
00753             io.comment (comment);
00754 
00755         }
00756     }
00757   
00758   // const Real 
00759   //   dt   = pl.get_elapsed_time(),
00760   //   rate = total_read_size*sizeof(Number)/dt; 
00761 
00762   // std::cout << "Read " << total_read_size << " \"Number\" values\n"
00763   //        << " Elapsed time = " << dt << '\n'
00764   //        << " Rate = " << rate/1.e6 << "(MB/sec)\n\n";
00765 
00766   // pl.pop("read_serialized_data");
00767 }

dof_id_type libMesh::System::read_serialized_vectors ( Xdr &  io,
const std::vector< NumericVector< Number > * > &  vectors 
) const [inherited]

Read a number of identically distributed vectors. This method allows for optimization for the multiple vector case by only communicating the metadata once.

Definition at line 2164 of file system_io.C.

References libMesh::Xdr::data(), libMesh::FEType::family, libMesh::System::get_mesh(), libMesh::MeshBase::n_elem(), libMesh::MeshTools::n_elem(), libMesh::MeshBase::n_nodes(), n_nodes, libMesh::System::n_vars(), libMesh::processor_id(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::Xdr::reading(), libMeshEnums::SCALAR, libMesh::Variable::type(), and libMesh::System::variable().

02166 {
02167   parallel_only();
02168 
02169   // Error checking
02170 // #ifndef NDEBUG
02171 //   // In parallel we better be reading a parallel vector -- if not
02172 //   // we will not set all of its components below!!
02173 //   if (libMesh::n_processors() > 1)
02174 //     {
02175 //       libmesh_assert (vec.type() == PARALLEL ||
02176 //                    vec.type() == GHOSTED);
02177 //     }
02178 // #endif
02179 
02180   libmesh_assert (io.reading());
02181 
02182   // sizes
02183   unsigned int num_vecs=0;
02184   dof_id_type vector_length=0;
02185 
02186   if (libMesh::processor_id() == 0) 
02187     {
02188       // Get the number of vectors
02189       io.data(num_vecs);
02190       // Get the buffer size
02191       io.data(vector_length);
02192 
02193       libmesh_assert_equal_to (num_vecs, vectors.size());
02194       
02195       if (num_vecs != 0)
02196         {
02197           libmesh_assert_not_equal_to (vectors[0], 0);
02198           libmesh_assert_equal_to     (vectors[0]->size(), vector_length);
02199         }
02200     }
02201 
02202   // no need to actually communicate these.
02203   // CommWorld.broadcast(num_vecs);
02204   // CommWorld.broadcast(vector_length);
02205   
02206   // Cache these - they are not free!
02207   const dof_id_type
02208     n_nodes = this->get_mesh().n_nodes(),
02209     n_elem  = this->get_mesh().n_elem();  
02210 
02211   dof_id_type read_length = 0.;
02212   
02213   //---------------------------------
02214   // Collect the values for all nodes
02215   read_length +=
02216     this->read_serialized_blocked_dof_objects (n_nodes,
02217                                                this->get_mesh().local_nodes_begin(),
02218                                                this->get_mesh().local_nodes_end(),
02219                                                io,
02220                                                vectors);
02221 
02222   //------------------------------------
02223   // Collect the values for all elements
02224   read_length +=
02225     this->read_serialized_blocked_dof_objects (n_elem,
02226                                                this->get_mesh().local_elements_begin(),
02227                                                this->get_mesh().local_elements_end(),
02228                                                io,
02229                                                vectors);
02230 
02231   //-------------------------------------------
02232   // Finally loop over all the SCALAR variables
02233   for (unsigned int vec=0; vec<vectors.size(); vec++)
02234     for (unsigned int var=0; var<this->n_vars(); var++)
02235       if(this->variable(var).type().family == SCALAR)
02236         {
02237           libmesh_assert_not_equal_to (vectors[vec], 0);
02238           
02239           read_length +=
02240             this->read_SCALAR_dofs (var, io, *vectors[vec]);
02241         }
02242 
02243   //---------------------------------------
02244   // last step - must close all the vectors
02245   for (unsigned int vec=0; vec<vectors.size(); vec++)
02246     {
02247       libmesh_assert_not_equal_to (vectors[vec], 0);      
02248       vectors[vec]->close();
02249     }
02250       
02251   return read_length;
02252 }

virtual void libMesh::RBConstruction::recompute_all_residual_terms ( const bool  compute_inner_products = true  )  [virtual, inherited]

This function computes all of the residual representors, can be useful when restarting a basis training computation. If compute_inner_products is false, we just compute the residual Riesz representors, whereas if true, we also compute all the corresponding inner product terms.

void libMesh::LinearImplicitSystem::reinit (  )  [virtual, inherited]

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

Reimplemented from libMesh::ImplicitSystem.

Reimplemented in libMesh::NewmarkSystem, and libMesh::TransientSystem< RBConstruction >.

Definition at line 85 of file linear_implicit_system.C.

References libMesh::LinearImplicitSystem::linear_solver, and libMesh::ImplicitSystem::reinit().

00086 {
00087   // re-initialize the linear solver interface
00088   linear_solver->clear();
00089 
00090   // initialize parent data
00091   Parent::reinit();
00092 }

void libMesh::LinearImplicitSystem::release_linear_solver ( LinearSolver< Number > *   )  const [virtual, inherited]

Releases a pointer to a linear solver acquired by this->get_linear_solver()

Reimplemented from libMesh::ImplicitSystem.

Definition at line 367 of file linear_implicit_system.C.

00368 {
00369 }

void libMesh::System::remove_vector ( const std::string &  vec_name  )  [inherited]

Removes the additional vector vec_name from this system

Definition at line 711 of file system.C.

References libMesh::System::_vector_projections, libMesh::System::_vector_types, libMesh::System::_vectors, and libMesh::System::have_vector().

00712 {
00713   //Return if the vector does not exist
00714   if ( !(this->have_vector(vec_name)) )
00715     return;
00716 
00717   _vectors[vec_name]->clear();
00718   delete _vectors[vec_name];
00719   _vectors[vec_name] = NULL;
00720 
00721   _vectors.erase(vec_name);
00722   _vector_projections.erase(vec_name);
00723   _vector_types.erase(vec_name);
00724 }

SparseMatrix< Number > * libMesh::ImplicitSystem::request_matrix ( const std::string &  mat_name  )  [inherited]
Returns:
a writable pointer to this system's additional matrix named mat_name, or returns NULL if no matrix by that name exists.

Definition at line 245 of file implicit_system.C.

References libMesh::ImplicitSystem::_matrices.

00246 {
00247   // Make sure the matrix exists
00248   matrices_iterator pos = _matrices.find (mat_name);
00249 
00250   if (pos == _matrices.end())
00251     return NULL;
00252 
00253   return pos->second;
00254 }

const SparseMatrix< Number > * libMesh::ImplicitSystem::request_matrix ( const std::string &  mat_name  )  const [inherited]
Returns:
a const pointer to this system's additional matrix named mat_name, or returns NULL if no matrix by that name exists.

Definition at line 232 of file implicit_system.C.

References libMesh::ImplicitSystem::_matrices.

Referenced by libMesh::ImplicitSystem::sensitivity_solve(), libMesh::NewtonSolver::solve(), and libMesh::LinearImplicitSystem::solve().

00233 {
00234   // Make sure the matrix exists
00235   const_matrices_iterator pos = _matrices.find (mat_name);
00236 
00237   if (pos == _matrices.end())
00238     return NULL;
00239 
00240   return pos->second;
00241 }

NumericVector< Number > * libMesh::System::request_vector ( const unsigned int  vec_num  )  [inherited]
Returns:
a writeable pointer to this system's additional vector number vec_num (where the vectors are counted starting with 0), or returns NULL if the system has no such vector.

Definition at line 767 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00768 {
00769   vectors_iterator v = vectors_begin();
00770   vectors_iterator v_end = vectors_end();
00771   unsigned int num = 0;
00772   while((num<vec_num) && (v!=v_end))
00773     {
00774       num++;
00775       ++v;
00776     }
00777   if (v==v_end)
00778     return NULL;
00779   return v->second;
00780 }

const NumericVector< Number > * libMesh::System::request_vector ( const unsigned int  vec_num  )  const [inherited]
Returns:
a const pointer to this system's additional vector number vec_num (where the vectors are counted starting with 0), or returns NULL if the system has no such vector.

Definition at line 750 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00751 {
00752   const_vectors_iterator v = vectors_begin();
00753   const_vectors_iterator v_end = vectors_end();
00754   unsigned int num = 0;
00755   while((num<vec_num) && (v!=v_end))
00756     {
00757       num++;
00758       ++v;
00759     }
00760   if (v==v_end)
00761     return NULL;
00762   return v->second;
00763 }

NumericVector< Number > * libMesh::System::request_vector ( const std::string &  vec_name  )  [inherited]
Returns:
a pointer to the vector if this System has a vector associated with the given name, NULL otherwise.

Definition at line 738 of file system.C.

References libMesh::System::_vectors.

00739 {
00740   vectors_iterator pos = _vectors.find(vec_name);
00741 
00742   if (pos == _vectors.end())
00743     return NULL;
00744 
00745   return pos->second;
00746 }

const NumericVector< Number > * libMesh::System::request_vector ( const std::string &  vec_name  )  const [inherited]
Returns:
a const pointer to the vector if this System has a vector associated with the given name, NULL otherwise.

Definition at line 726 of file system.C.

References libMesh::System::_vectors.

00727 {
00728   const_vectors_iterator pos = _vectors.find(vec_name);
00729 
00730   if (pos == _vectors.end())
00731     return NULL;
00732 
00733   return pos->second;
00734 }

void libMesh::RBConstructionBase< LinearImplicitSystem >::reset_alternative_solver ( AutoPtr< LinearSolver< Number > > &  ls,
const std::pair< std::string, std::string > &  orig 
) [inherited]

Resets the PC (and iterative solver, if desired) in the passed-in LinearSolver object to the values specified in the pair of strings passed as the second argument. If the "alternative_solver" string, defined below, is "unchanged", this function does nothing.

void libMesh::LinearImplicitSystem::restrict_solve_to ( const SystemSubset *  subset,
const SubsetSolveMode  subset_solve_mode = SUBSET_ZERO 
) [virtual, inherited]

After calling this method, any solve will be limited to the given subset. To disable this mode, call this method with subset being a NULL pointer.

Reimplemented from libMesh::System.

Definition at line 96 of file linear_implicit_system.C.

References libMesh::LinearImplicitSystem::_subset, libMesh::LinearImplicitSystem::_subset_solve_mode, and libMesh::SystemSubset::get_system().

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

00098 {
00099   _subset = subset;
00100   _subset_solve_mode = subset_solve_mode;
00101   if(subset!=NULL)
00102     {
00103       libmesh_assert_equal_to (&subset->get_system(), this);
00104     }
00105 }

void libMesh::System::restrict_vectors (  )  [virtual, inherited]

Restrict vectors after the mesh has coarsened

Definition at line 316 of file system.C.

References libMesh::System::_dof_map, libMesh::System::_solution_projection, libMesh::System::_vector_projections, libMesh::System::_vector_types, libMesh::System::_vectors, libMesh::System::current_local_solution, libMesh::err, libMeshEnums::GHOSTED, libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), libMesh::System::project_vector(), and libMesh::System::solution.

Referenced by libMesh::System::prolong_vectors(), and libMesh::EquationSystems::reinit().

00317 {
00318 #ifdef LIBMESH_ENABLE_AMR
00319   // Restrict the _vectors on the coarsened cells
00320   for (vectors_iterator pos = _vectors.begin(); pos != _vectors.end(); ++pos)
00321     {
00322       NumericVector<Number>* v = pos->second;
00323 
00324       if (_vector_projections[pos->first])
00325         this->project_vector (*v);
00326       else
00327       {
00328         ParallelType type = _vector_types[pos->first];
00329 
00330         if(type == GHOSTED)
00331         {
00332 #ifdef LIBMESH_ENABLE_GHOSTED
00333           pos->second->init (this->n_dofs(), this->n_local_dofs(),
00334                          _dof_map->get_send_list(), false,
00335                          GHOSTED);
00336 #else
00337           libMesh::err << "Cannot initialize ghosted vectors when they are not enabled." << std::endl;
00338           libmesh_error();
00339 #endif
00340         }
00341         else
00342           pos->second->init (this->n_dofs(), this->n_local_dofs(), false, type);
00343       }
00344     }
00345 
00346   const std::vector<dof_id_type>& send_list = _dof_map->get_send_list ();
00347 
00348   // Restrict the solution on the coarsened cells
00349   if (_solution_projection)
00350     this->project_vector (*solution);
00351 
00352 #ifdef LIBMESH_ENABLE_GHOSTED
00353   current_local_solution->init(this->n_dofs(),
00354                                this->n_local_dofs(), send_list,
00355                                false, GHOSTED);
00356 #else
00357   current_local_solution->init(this->n_dofs());
00358 #endif
00359 
00360   if (_solution_projection)
00361     solution->localize (*current_local_solution, send_list);
00362 
00363 #endif // LIBMESH_ENABLE_AMR
00364 }

std::pair< unsigned int, Real > libMesh::ImplicitSystem::sensitivity_solve ( const ParameterVector &  parameters  )  [virtual, inherited]

Assembles & solves the linear system(s) (dR/du)*u_p = -dR/dp, for those parameters contained within parameters.

Returns a pair with the total number of linear iterations performed and the (sum of the) final residual norms

Reimplemented from libMesh::System.

Definition at line 321 of file implicit_system.C.

References libMesh::System::add_sensitivity_solution(), libMesh::System::assemble_before_solve, libMesh::ImplicitSystem::assemble_residual_derivatives(), libMesh::ImplicitSystem::assembly(), libMesh::SparseMatrix< T >::close(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_dof_map(), libMesh::ImplicitSystem::get_linear_solve_parameters(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::System::get_sensitivity_rhs(), libMesh::System::get_sensitivity_solution(), libMesh::ImplicitSystem::matrix, libMesh::ImplicitSystem::release_linear_solver(), libMesh::ImplicitSystem::request_matrix(), libMesh::ParameterVector::size(), and libMesh::LinearSolver< T >::solve().

00322 {
00323   // Log how long the linear solve takes.
00324   START_LOG("sensitivity_solve()", "ImplicitSystem");
00325 
00326   // The forward system should now already be solved.
00327   // Now assemble the corresponding sensitivity system.
00328 
00329   if (this->assemble_before_solve)
00330     {
00331       // Build the Jacobian
00332       this->assembly(false, true);
00333       this->matrix->close();
00334 
00335       // Reset and build the RHS from the residual derivatives
00336       this->assemble_residual_derivatives(parameters);
00337     }
00338 
00339   // The sensitivity problem is linear
00340   LinearSolver<Number> *linear_solver = this->get_linear_solver();
00341 
00342   // Our iteration counts and residuals will be sums of the individual
00343   // results
00344   std::pair<unsigned int, Real> solver_params =
00345     this->get_linear_solve_parameters();
00346   std::pair<unsigned int, Real> totalrval = std::make_pair(0,0.0);
00347 
00348   // Solve the linear system.
00349   SparseMatrix<Number> *pc = this->request_matrix("Preconditioner");
00350   for (unsigned int p=0; p != parameters.size(); ++p)
00351     {
00352       std::pair<unsigned int, Real> rval =
00353         linear_solver->solve (*matrix, pc,
00354                               this->add_sensitivity_solution(p),
00355                               this->get_sensitivity_rhs(p),
00356                               solver_params.second,
00357                               solver_params.first);
00358 
00359       totalrval.first  += rval.first;
00360       totalrval.second += rval.second;
00361     }
00362 
00363   // The linear solver may not have fit our constraints exactly
00364 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00365   for (unsigned int p=0; p != parameters.size(); ++p)
00366     this->get_dof_map().enforce_constraints_exactly
00367       (*this, &this->get_sensitivity_solution(p),
00368        /* homogeneous = */ true);
00369 #endif
00370 
00371   this->release_linear_solver(linear_solver);
00372 
00373   // Stop logging the nonlinear solve
00374   STOP_LOG("sensitivity_solve()", "ImplicitSystem");
00375 
00376   return totalrval;
00377 }

void libMesh::System::set_adjoint_already_solved ( bool  setting  )  [inline, inherited]

Setter for the adjoint_already_solved boolean

Definition at line 365 of file system.h.

References libMesh::System::adjoint_already_solved.

00366   { adjoint_already_solved = setting;}

std::pair<std::string,std::string> libMesh::RBConstructionBase< LinearImplicitSystem >::set_alternative_solver ( AutoPtr< LinearSolver< Number > > &  ls  )  [inherited]

Changes the current PC (and iterative solver, if desired) in the passed-in LinearSolver object to an alternative solver specified by the alternative_solver string stored in this class. You might use this to e.g. switch to a sparse direct solver for the multiple RHS solves executed during the update_residual_terms function. The return strings are names of the original PC and KSP objects, you can reset these using the reset_alternative_solver() function below.

void libMesh::System::set_basic_system_only (  )  [inline, inherited]

Sets the system to be "basic only": i.e. advanced system components such as ImplicitSystem matrices may not be initialized. This is useful for efficiency in certain utility programs that never use System::solve(). This method must be called after the System or derived class is created but before it is initialized; e.g. from within EquationSystems::read()

Definition at line 1890 of file system.h.

References libMesh::System::_basic_system_only.

Referenced by libMesh::EquationSystems::_read_impl().

01891 {
01892   _basic_system_only = true;
01893 }

void libMesh::RBConstruction::set_constraint_assembly ( ElemAssembly &  constraint_assembly_in  )  [inherited]

Set the rb_assembly_expansion object.

virtual void libMesh::RBConstruction::set_context_solution_vec ( NumericVector< Number > &  vec  )  [protected, virtual, inherited]

Set current_local_solution = vec so that we can access vec from FEMContext during assembly. Override in subclasses if different behavior is required.

void libMesh::RBConstructionBase< LinearImplicitSystem >::set_deterministic_training_parameter_name ( const std::string  name  )  [inherited]

Set the name of the parameter that we will generate deterministic training parameters for. Defaults to "NONE".

void libMesh::RBConstructionBase< LinearImplicitSystem >::set_deterministic_training_parameter_repeats ( unsigned int  repeats  )  [inherited]

Set the number of times each sample of the deterministic training parameter is repeated.

void libMesh::RBConstruction::set_inner_product_assembly ( ElemAssembly &  inner_product_assembly_in  )  [inherited]

Set the rb_assembly_expansion object.

virtual void libMesh::RBConstruction::set_Nmax ( unsigned int  Nmax  )  [virtual, inherited]
void libMesh::RBParametrized::set_parameters ( const RBParameters &  params  )  [inherited]

Set the current parameters to params

void libMesh::RBConstructionBase< LinearImplicitSystem >::set_params_from_training_set ( unsigned int  index  )  [protected, inherited]

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

virtual void libMesh::RBConstructionBase< LinearImplicitSystem >::set_params_from_training_set_and_broadcast ( unsigned int  index  )  [protected, virtual, inherited]

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

void libMesh::RBConstruction::set_quiet_mode ( bool  quiet_mode_in  )  [inline, inherited]

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

Definition at line 190 of file rb_construction.h.

References libMesh::RBConstruction::quiet_mode.

00191     { this->quiet_mode = quiet_mode_in; }

void libMesh::RBConstruction::set_rb_assembly_expansion ( RBAssemblyExpansion &  rb_assembly_expansion_in  )  [inherited]

Set the rb_assembly_expansion object.

void libMesh::RBConstruction::set_rb_evaluation ( RBEvaluation &  rb_eval_in  )  [inherited]

Set the RBEvaluation object.

void libMesh::RBConstructionBase< LinearImplicitSystem >::set_training_random_seed ( unsigned int  seed  )  [inherited]

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

void libMesh::RBConstruction::set_training_tolerance ( Real  new_training_tolerance  )  [inline, inherited]

Get/set the tolerance for the basis training.

Definition at line 174 of file rb_construction.h.

References libMesh::RBConstruction::training_tolerance.

00175     {this->training_tolerance = new_training_tolerance; }

void libMesh::System::set_vector_preservation ( const std::string &  vec_name,
bool  preserve 
) [inherited]

Allows one to set the boolean controlling whether the vector identified by vec_name should be "preserved": projected to new meshes, saved, etc.

Definition at line 887 of file system.C.

References libMesh::System::_vector_projections.

00889 {
00890   _vector_projections[vec_name] = preserve;
00891 }

void libMesh::LinearImplicitSystem::solve (  )  [virtual, inherited]

Assembles & solves the linear system A*x=b.

Reimplemented from libMesh::ExplicitSystem.

Reimplemented in libMesh::FrequencySystem.

Definition at line 109 of file linear_implicit_system.C.

References libMesh::LinearImplicitSystem::_final_linear_residual, libMesh::LinearImplicitSystem::_n_linear_iterations, libMesh::LinearImplicitSystem::_shell_matrix, libMesh::LinearImplicitSystem::_subset, libMesh::LinearImplicitSystem::_subset_solve_mode, libMesh::LinearImplicitSystem::assemble(), libMesh::System::assemble_before_solve, libMesh::SystemSubset::dof_ids(), libMesh::Parameters::get(), libMesh::System::get_equation_systems(), libMesh::LinearImplicitSystem::linear_solver, libMesh::ImplicitSystem::matrix, libMesh::EquationSystems::parameters, libMesh::Real, libMesh::ImplicitSystem::request_matrix(), libMesh::ExplicitSystem::rhs, libMesh::System::solution, and libMesh::System::update().

00110 {
00111   if (this->assemble_before_solve)
00112     // Assemble the linear system
00113     this->assemble ();
00114 
00115   // Log how long the linear solve takes.
00116   // This gets done by the LinearSolver classes now [RHS]
00117   // START_LOG("solve()", "System");
00118 
00119   // Get a reference to the EquationSystems
00120   const EquationSystems& es =
00121     this->get_equation_systems();
00122 
00123   // Get the user-specifiied linear solver tolerance
00124   const Real tol            =
00125     es.parameters.get<Real>("linear solver tolerance");
00126 
00127   // Get the user-specified maximum # of linear solver iterations
00128   const unsigned int maxits =
00129     es.parameters.get<unsigned int>("linear solver maximum iterations");
00130 
00131   if(_subset!=NULL)
00132     {
00133       linear_solver->restrict_solve_to(&_subset->dof_ids(),_subset_solve_mode);
00134     }
00135 
00136   // Solve the linear system.  Several cases:
00137   std::pair<unsigned int, Real> rval = std::make_pair(0,0.0);
00138   if(_shell_matrix)
00139     // 1.) Shell matrix with or without user-supplied preconditioner.
00140     rval = linear_solver->solve(*_shell_matrix, this->request_matrix("Preconditioner"), *solution, *rhs, tol, maxits);
00141   else
00142     // 2.) No shell matrix, with or without user-supplied preconditioner
00143     rval = linear_solver->solve (*matrix, this->request_matrix("Preconditioner"), *solution, *rhs, tol, maxits);
00144 
00145   if(_subset!=NULL)
00146     {
00147       linear_solver->restrict_solve_to(NULL);
00148     }
00149 
00150   // Store the number of linear iterations required to
00151   // solve and the final residual.
00152   _n_linear_iterations   = rval.first;
00153   _final_linear_residual = rval.second;
00154 
00155   // Stop logging the linear solve
00156   // This gets done by the LinearSolver classes now [RHS]
00157   // STOP_LOG("solve()", "System");
00158 
00159   // Update the system after the solve
00160   this->update();
00161 }

sys_type& libMesh::RBConstruction::system (  )  [inline, inherited]
Returns:
a clever pointer to the system.

Reimplemented from libMesh::RBConstructionBase< LinearImplicitSystem >.

Reimplemented in libMesh::TransientSystem< RBConstruction >.

Definition at line 113 of file rb_construction.h.

00113 { return *this; }

virtual std::string libMesh::RBConstruction::system_type (  )  const [virtual, inherited]
Returns:
a string indicating the type of the system.

Reimplemented from libMesh::LinearImplicitSystem.

Reimplemented in libMesh::TransientSystem< RBConstruction >.

virtual Real libMesh::RBEIMConstruction::train_reduced_basis ( const std::string &  directory_name = "offline_data",
const bool  resize_rb_eval_data = true 
) [virtual]

Override train_reduced_basis to first initialize _parametrized_functions_in_training_set.

Reimplemented from libMesh::RBConstruction.

virtual void libMesh::RBConstruction::truth_assembly (  )  [protected, virtual, inherited]

Assemble the truth matrix and right-hand side for current_parameters.

Reimplemented in libMesh::TransientRBConstruction.

virtual Real libMesh::RBEIMConstruction::truth_solve ( int  plot_solution  )  [virtual]

Load the truth representation of the parametrized function at the current parameters into the solution vector. The truth representation is the projection of parametrized_function into the finite element space. If plot_solution > 0 the solution will be plotted to an output file.

Reimplemented from libMesh::RBConstruction.

void libMesh::System::update (  )  [virtual, inherited]

Update the local values to reflect the solution on neighboring processors.

Definition at line 410 of file system.C.

References libMesh::System::_dof_map, libMesh::System::current_local_solution, and libMesh::System::solution.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::FEMSystem::assemble_qoi(), libMesh::FEMSystem::assemble_qoi_derivative(), libMesh::NonlinearImplicitSystem::assembly(), libMesh::FEMSystem::assembly(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::GMVIO::copy_nodal_solution(), DMFunction_libMesh(), DMJacobian_libMesh(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::FEMSystem::mesh_position_get(), libMesh::FEMSystem::postprocess(), libMesh::NonlinearImplicitSystem::solve(), libMesh::NewtonSolver::solve(), libMesh::LinearImplicitSystem::solve(), and libMesh::ExplicitSystem::solve().

00411 {
00412   libmesh_assert(solution->closed());
00413 
00414   const std::vector<dof_id_type>& send_list = _dof_map->get_send_list ();
00415 
00416   // Check sizes
00417   libmesh_assert_equal_to (current_local_solution->size(), solution->size());
00418 // More processors than elements => empty send_list
00419 //  libmesh_assert (!send_list.empty());
00420   libmesh_assert_less_equal (send_list.size(), solution->size());
00421 
00422   // Create current_local_solution from solution.  This will
00423   // put a local copy of solution into current_local_solution.
00424   // Only the necessary values (specified by the send_list)
00425   // are copied to minimize communication
00426   solution->localize (*current_local_solution, send_list);
00427 }

void libMesh::System::update_global_solution ( std::vector< Number > &  global_soln,
const unsigned int  dest_proc 
) const [inherited]

Fill the input vector global_soln so that it contains the global solution on processor dest_proc. Requires communication with all other processors.

Definition at line 665 of file system.C.

References libMesh::System::solution.

00667 {
00668   global_soln.resize        (solution->size());
00669 
00670   solution->localize_to_one (global_soln, dest_proc);
00671 }

void libMesh::System::update_global_solution ( std::vector< Number > &  global_soln  )  const [inherited]

Fill the input vector global_soln so that it contains the global solution on all processors. Requires communication with all other processors.

Definition at line 656 of file system.C.

References libMesh::System::solution.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_solution_vector(), libMesh::ExactErrorEstimator::estimate_error(), and libMesh::EquationSystems::get_solution().

00657 {
00658   global_soln.resize (solution->size());
00659 
00660   solution->localize (global_soln);
00661 }

void libMesh::RBConstruction::update_greedy_param_list (  )  [protected, inherited]

Update the list of Greedily chosen parameters with current_parameters.

virtual void libMesh::RBEIMConstruction::update_RB_system_matrices (  )  [protected, virtual]

Compute the reduced basis matrices for the current basis. Overload to update the inner product matrix that is used to compute the best fit to parametrized_function.

Reimplemented from libMesh::RBConstruction.

virtual void libMesh::RBConstruction::update_residual_terms ( bool  compute_inner_products = true  )  [protected, virtual, inherited]

Compute the terms that are combined `online' to determine the dual norm of the residual. By default, inner product terms are also computed, but you can turn this feature off e.g. if you are already reading in that data from files.

Reimplemented in libMesh::TransientRBConstruction.

virtual void libMesh::RBEIMConstruction::update_system (  )  [protected, virtual]

Update the system after enriching the RB space; this calls a series of functions to update the system properly.

Reimplemented from libMesh::RBConstruction.

void libMesh::System::user_assembly (  )  [virtual, inherited]

Calls user's attached assembly function, or is overloaded by the user in derived classes.

Definition at line 1909 of file system.C.

References libMesh::System::_assemble_system_function, libMesh::System::_assemble_system_object, libMesh::System::_equation_systems, libMesh::System::Assembly::assemble(), and libMesh::System::name().

Referenced by libMesh::System::assemble().

01910 {
01911   // Call the user-provided assembly function,
01912   // if it was provided
01913   if (_assemble_system_function != NULL)
01914     this->_assemble_system_function (_equation_systems, this->name());
01915 
01916   // ...or the user-provided assembly object.
01917   else if (_assemble_system_object != NULL)
01918     this->_assemble_system_object->assemble();
01919 }

void libMesh::System::user_constrain (  )  [virtual, inherited]

Calls user's attached constraint function, or is overloaded by the user in derived classes.

Definition at line 1923 of file system.C.

References libMesh::System::_constrain_system_function, libMesh::System::_constrain_system_object, libMesh::System::_equation_systems, libMesh::System::Constraint::constrain(), and libMesh::System::name().

Referenced by libMesh::EquationSystems::allgather(), libMesh::System::init_data(), and libMesh::EquationSystems::reinit().

01924 {
01925   // Call the user-provided constraint function,
01926   // if it was provided
01927   if (_constrain_system_function!= NULL)
01928     this->_constrain_system_function(_equation_systems, this->name());
01929 
01930   // ...or the user-provided constraint object.
01931   else if (_constrain_system_object != NULL)
01932     this->_constrain_system_object->constrain();
01933 }

void libMesh::System::user_initialization (  )  [virtual, inherited]

Calls user's attached initialization function, or is overloaded by the user in derived classes.

Definition at line 1895 of file system.C.

References libMesh::System::_equation_systems, libMesh::System::_init_system_function, libMesh::System::_init_system_object, libMesh::System::Initialization::initialize(), and libMesh::System::name().

Referenced by libMesh::System::init(), and libMesh::NewmarkSystem::initial_conditions().

01896 {
01897   // Call the user-provided intialization function,
01898   // if it was provided
01899   if (_init_system_function != NULL)
01900     this->_init_system_function (_equation_systems, this->name());
01901 
01902   // ...or the user-provided initialization object.
01903   else if (_init_system_object != NULL)
01904     this->_init_system_object->initialize();
01905 }

void libMesh::System::user_QOI ( const QoISet &  qoi_indices  )  [virtual, inherited]

Calls user's attached quantity of interest function, or is overloaded by the user in derived classes.

Definition at line 1937 of file system.C.

References libMesh::System::_equation_systems, libMesh::System::_qoi_evaluate_function, libMesh::System::_qoi_evaluate_object, libMesh::System::name(), and libMesh::System::QOI::qoi().

Referenced by libMesh::System::assemble_qoi().

01938 {
01939   // Call the user-provided quantity of interest function,
01940   // if it was provided
01941   if (_qoi_evaluate_function != NULL)
01942     this->_qoi_evaluate_function(_equation_systems, this->name(), qoi_indices);
01943 
01944   // ...or the user-provided QOI function object.
01945   else if (_qoi_evaluate_object != NULL)
01946     this->_qoi_evaluate_object->qoi(qoi_indices);
01947 }

void libMesh::System::user_QOI_derivative ( const QoISet &  qoi_indices  )  [virtual, inherited]

Calls user's attached quantity of interest derivative function, or is overloaded by the user in derived classes.

Definition at line 1951 of file system.C.

References libMesh::System::_equation_systems, libMesh::System::_qoi_evaluate_derivative_function, libMesh::System::_qoi_evaluate_derivative_object, libMesh::System::name(), and libMesh::System::QOIDerivative::qoi_derivative().

Referenced by libMesh::System::assemble_qoi_derivative().

01952 {
01953   // Call the user-provided quantity of interest derivative,
01954   // if it was provided
01955   if (_qoi_evaluate_derivative_function != NULL)
01956     this->_qoi_evaluate_derivative_function(_equation_systems, this->name(), qoi_indices);
01957 
01958   // ...or the user-provided QOI derivative function object.
01959   else if (_qoi_evaluate_derivative_object != NULL)
01960     this->_qoi_evaluate_derivative_object->qoi_derivative(qoi_indices);
01961 }

const Variable & libMesh::System::variable ( unsigned int  var  )  const [inline, inherited]

Return a constant reference to Variable var.

Definition at line 1926 of file system.h.

References libMesh::System::_variables.

Referenced by libMesh::EquationSystems::build_solution_vector(), libMesh::WrappedFunction< Output >::component(), libMesh::EquationSystems::get_solution(), libMesh::WrappedFunction< Output >::operator()(), libMesh::System::project_vector(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_vector(), and libMesh::System::write_serialized_vectors().

01927 {
01928   libmesh_assert_less (i, _variables.size());
01929 
01930   return _variables[i];
01931 }

const VariableGroup & libMesh::System::variable_group ( unsigned int  vg  )  const [inline, inherited]

Return a constant reference to VariableGroup vg.

Definition at line 1936 of file system.h.

References libMesh::System::_variable_groups.

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

01937 {
01938   libmesh_assert_less (vg, _variable_groups.size());
01939 
01940   return _variable_groups[vg];
01941 }

const std::string & libMesh::System::variable_name ( const unsigned int  i  )  const [inline, inherited]
Returns:
the name of variable i.

Definition at line 1946 of file system.h.

References libMesh::System::_variables.

Referenced by libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::KellyErrorEstimator::boundary_side_integration(), libMesh::DiscontinuityMeasure::boundary_side_integration(), libMesh::WrappedFunction< Output >::component(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ExactSolution::ExactSolution(), libMesh::EquationSystems::get_solution(), libMesh::WrappedFunction< Output >::operator()(), and libMesh::System::write_header().

01947 {
01948   libmesh_assert_less (i, _variables.size());
01949 
01950   return _variables[i].name();
01951 }

unsigned short int libMesh::System::variable_number ( const std::string &  var  )  const [inherited]
Returns:
the variable number assoicated with the user-specified variable named var.

Definition at line 1239 of file system.C.

References libMesh::System::_variable_numbers, libMesh::System::_variables, and libMesh::err.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::GMVIO::copy_nodal_solution(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::System::read_header(), libMesh::System::variable_scalar_number(), libMesh::System::variable_type(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

01240 {
01241   // Make sure the variable exists
01242   std::map<std::string, unsigned short int>::const_iterator
01243     pos = _variable_numbers.find(var);
01244 
01245   if (pos == _variable_numbers.end())
01246     {
01247       libMesh::err << "ERROR: variable "
01248                     << var
01249                     << " does not exist in this system!"
01250                     << std::endl;
01251       libmesh_error();
01252     }
01253   libmesh_assert_equal_to (_variables[pos->second].name(), var);
01254 
01255   return pos->second;
01256 }

unsigned int libMesh::System::variable_scalar_number ( unsigned int  var_num,
unsigned int  component 
) const [inline, inherited]
Returns:
an index, starting from 0 for the first component of the first variable, and incrementing for each component of each (potentially vector-valued) variable in the system in order. For systems with only scalar-valued variables, this will be the same as var_num

Irony: currently our only non-scalar-valued variable type is SCALAR.

Definition at line 1967 of file system.h.

References libMesh::System::_variables.

01969 {
01970   return _variables[var_num].first_scalar_number() + component;
01971 }

unsigned int libMesh::System::variable_scalar_number ( const std::string &  var,
unsigned int  component 
) const [inline, inherited]
Returns:
an index, starting from 0 for the first component of the first variable, and incrementing for each component of each (potentially vector-valued) variable in the system in order. For systems with only scalar-valued variables, this will be the same as variable_number(var)

Irony: currently our only non-scalar-valued variable type is SCALAR.

Definition at line 1957 of file system.h.

References libMesh::System::variable_number().

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::WrappedFunction< Output >::component(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::WrappedFunction< Output >::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), libMesh::ProjectSolution::operator()(), and libMesh::System::project_vector().

01959 {
01960   return variable_scalar_number(this->variable_number(var), component);
01961 }

const FEType & libMesh::System::variable_type ( const std::string &  var  )  const [inline, inherited]
Returns:
the finite element type for variable var.

Definition at line 1986 of file system.h.

References libMesh::System::_variables, and libMesh::System::variable_number().

01987 {
01988   return _variables[this->variable_number(var)].type();
01989 }

const FEType & libMesh::System::variable_type ( const unsigned int  i  )  const [inline, inherited]
Returns:
the finite element type variable number i.

Definition at line 1976 of file system.h.

References libMesh::System::_variables.

Referenced by libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_solution_vector(), libMesh::GMVIO::copy_nodal_solution(), libMesh::FEMContext::FEMContext(), libMesh::EquationSystems::get_solution(), libMesh::System::write_header(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

01977 {
01978   libmesh_assert_less (i, _variables.size());
01979 
01980   return _variables[i].type();
01981 }

const std::string & libMesh::System::vector_name ( const NumericVector< Number > &  vec_reference  )  const [inherited]
Returns:
the name of a system vector, given a reference to that vector

Definition at line 868 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00869 {
00870   const_vectors_iterator v = vectors_begin();
00871   const_vectors_iterator v_end = vectors_end();
00872 
00873   for(; v != v_end; ++v)
00874     {
00875       // Check if the current vector is the one whose name we want
00876       if(&vec_reference == v->second)
00877         break; // exit loop if it is
00878     }
00879 
00880   // Before returning, make sure we didnt loop till the end and not find any match
00881   libmesh_assert (v != v_end);
00882 
00883   // Return the string associated with the current vector
00884   return v->first;
00885 }

const std::string & libMesh::System::vector_name ( const unsigned int  vec_num  )  const [inherited]
Returns:
the name of this system's additional vector number vec_num (where the vectors are counted starting with 0).

Definition at line 854 of file system.C.

References libMesh::System::vectors_begin(), and libMesh::System::vectors_end().

00855 {
00856   const_vectors_iterator v = vectors_begin();
00857   const_vectors_iterator v_end = vectors_end();
00858   unsigned int num = 0;
00859   while((num<vec_num) && (v!=v_end))
00860     {
00861       num++;
00862       ++v;
00863     }
00864   libmesh_assert (v != v_end);
00865   return v->first;
00866 }

bool libMesh::System::vector_preservation ( const std::string &  vec_name  )  const [inherited]
Returns:
the boolean describing whether the vector identified by vec_name should be "preserved": projected to new meshes, saved, etc.

Definition at line 895 of file system.C.

References libMesh::System::_vector_projections.

Referenced by libMesh::MemorySolutionHistory::store().

00896 {
00897   if (_vector_projections.find(vec_name) == _vector_projections.end())
00898     return false;
00899 
00900   return _vector_projections.find(vec_name)->second;
00901 }

System::const_vectors_iterator libMesh::System::vectors_begin (  )  const [inline, inherited]

Beginning of vectors container

Definition at line 2044 of file system.h.

References libMesh::System::_vectors.

02045 {
02046   return _vectors.begin();
02047 }

System::vectors_iterator libMesh::System::vectors_begin (  )  [inline, inherited]

Beginning of vectors container

Definition at line 2038 of file system.h.

References libMesh::System::_vectors.

Referenced by libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_vector(), libMesh::System::request_vector(), libMesh::MemorySolutionHistory::store(), libMesh::VTKIO::system_vectors_to_vtk(), and libMesh::System::vector_name().

02039 {
02040   return _vectors.begin();
02041 }

System::const_vectors_iterator libMesh::System::vectors_end (  )  const [inline, inherited]

End of vectors container

Definition at line 2056 of file system.h.

References libMesh::System::_vectors.

02057 {
02058   return _vectors.end();
02059 }

System::vectors_iterator libMesh::System::vectors_end (  )  [inline, inherited]

End of vectors container

Definition at line 2050 of file system.h.

References libMesh::System::_vectors.

Referenced by libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_vector(), libMesh::System::request_vector(), libMesh::MemorySolutionHistory::store(), libMesh::VTKIO::system_vectors_to_vtk(), and libMesh::System::vector_name().

02051 {
02052   return _vectors.end();
02053 }

std::pair< unsigned int, Real > libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve ( const ParameterVector &  parameters,
const ParameterVector &  weights,
const QoISet &  qoi_indices = QoISet() 
) [virtual, inherited]

Assembles & solves the linear system(s) (dR/du)^T*z_w = sum(w_p*(d^2q/dudp - d^2R/dudp*z)), for those parameters p contained within parameters, weighted by the values w_p found within weights.

Assumes that adjoint_solve has already calculated z for each qoi in qoi_indices.

Returns a pair with the total number of linear iterations performed and the (sum of the) final residual norms

Reimplemented from libMesh::System.

Definition at line 437 of file implicit_system.C.

References libMesh::System::add_weighted_sensitivity_adjoint_solution(), libMesh::ExplicitSystem::assemble_qoi_derivative(), libMesh::ImplicitSystem::assembly(), libMesh::NumericVector< T >::close(), libMesh::SparseMatrix< T >::close(), libMesh::ParameterVector::deep_copy(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_adjoint_rhs(), libMesh::System::get_adjoint_solution(), libMesh::System::get_dof_map(), libMesh::ImplicitSystem::get_linear_solve_parameters(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::SparseMatrix< T >::get_transpose(), libMesh::System::get_weighted_sensitivity_adjoint_solution(), libMesh::QoISet::has_index(), libMesh::ImplicitSystem::matrix, libMesh::System::qoi, libMesh::Real, libMesh::AutoPtr< Tp >::release(), libMesh::ImplicitSystem::release_linear_solver(), libMesh::ExplicitSystem::rhs, libMesh::LinearSolver< T >::solve(), libMesh::TOLERANCE, libMesh::ParameterVector::value_copy(), libMesh::SparseMatrix< T >::vector_mult_add(), and libMesh::NumericVector< T >::zero_clone().

00440 {
00441   // Log how long the linear solve takes.
00442   START_LOG("weighted_sensitivity_adjoint_solve()", "ImplicitSystem");
00443 
00444   // We currently get partial derivatives via central differencing
00445   const Real delta_p = TOLERANCE;
00446 
00447   // The forward system should now already be solved.
00448   // The adjoint system should now already be solved.
00449   // Now we're assembling a weighted sum of adjoint-adjoint systems:
00450   //
00451   // dR/du (u, sum_l(w_l*z^l)) = sum_l(w_l*(Q''_ul - R''_ul (u, z)))
00452 
00453   // We'll assemble the rhs first, because the R'' term will require
00454   // perturbing the jacobian
00455 
00456   // We'll use temporary rhs vectors, because we haven't (yet) found
00457   // any good reasons why users might want to save these:
00458 
00459   std::vector<NumericVector<Number> *> temprhs(this->qoi.size());
00460   for (unsigned int i=0; i != this->qoi.size(); ++i)
00461     if (qoi_indices.has_index(i))
00462       temprhs[i] = this->rhs->zero_clone().release();
00463 
00464   // We approximate the _l partial derivatives via a central
00465   // differencing perturbation in the w_l direction:
00466   //
00467   // sum_l(w_l*v_l) ~= (v(p + dp*w_l*e_l) - v(p - dp*w_l*e_l))/(2*dp)
00468 
00469   // PETSc doesn't implement SGEMX, so neither does NumericVector,
00470   // so we want to avoid calculating f -= R'*z.  We'll thus evaluate
00471   // the above equation by first adding -v(p+dp...), then multiplying
00472   // the intermediate result vectors by -1, then adding -v(p-dp...),
00473   // then finally dividing by 2*dp.
00474 
00475   ParameterVector oldparameters, parameterperturbation;
00476   parameters.deep_copy(oldparameters);
00477   weights.deep_copy(parameterperturbation);
00478   parameterperturbation *= delta_p;
00479   parameters += parameterperturbation;
00480 
00481   this->assembly(false, true);
00482   this->matrix->close();
00483 
00484   // Take the discrete adjoint, so that we can calculate R_u(u,z) with
00485   // a matrix-vector product of R_u and z.
00486   matrix->get_transpose(*matrix);
00487 
00488   this->assemble_qoi_derivative(qoi_indices);
00489   for (unsigned int i=0; i != this->qoi.size(); ++i)
00490     if (qoi_indices.has_index(i))
00491       {
00492         this->get_adjoint_rhs(i).close();
00493         *(temprhs[i]) -= this->get_adjoint_rhs(i);
00494         this->matrix->vector_mult_add(*(temprhs[i]), this->get_adjoint_solution(i));
00495         *(temprhs[i]) *= -1.0;
00496       }
00497 
00498   oldparameters.value_copy(parameters);
00499   parameterperturbation *= -1.0;
00500   parameters += parameterperturbation;
00501 
00502   this->assembly(false, true);
00503   this->matrix->close();
00504   matrix->get_transpose(*matrix);
00505 
00506   this->assemble_qoi_derivative(qoi_indices);
00507   for (unsigned int i=0; i != this->qoi.size(); ++i)
00508     if (qoi_indices.has_index(i))
00509       {
00510         this->get_adjoint_rhs(i).close();
00511         *(temprhs[i]) -= this->get_adjoint_rhs(i);
00512         this->matrix->vector_mult_add(*(temprhs[i]), this->get_adjoint_solution(i));
00513         *(temprhs[i]) /= (2.0*delta_p);
00514       }
00515 
00516   // Finally, assemble the jacobian at the non-perturbed parameter
00517   // values.  Ignore assemble_before_solve; if we had a good
00518   // non-perturbed matrix before we've already overwritten it.
00519   oldparameters.value_copy(parameters);
00520 
00521   // if (this->assemble_before_solve)
00522     {
00523       // Build the Jacobian
00524       this->assembly(false, true);
00525       this->matrix->close();
00526 
00527       // Take the discrete adjoint
00528       matrix->get_transpose(*matrix);
00529     }
00530 
00531   // The weighted adjoint-adjoint problem is linear
00532   LinearSolver<Number> *linear_solver = this->get_linear_solver();
00533 
00534   // Our iteration counts and residuals will be sums of the individual
00535   // results
00536   std::pair<unsigned int, Real> solver_params =
00537     this->get_linear_solve_parameters();
00538   std::pair<unsigned int, Real> totalrval = std::make_pair(0,0.0);
00539 
00540   for (unsigned int i=0; i != this->qoi.size(); ++i)
00541     if (qoi_indices.has_index(i))
00542       {
00543         const std::pair<unsigned int, Real> rval =
00544           linear_solver->solve (*matrix, this->add_weighted_sensitivity_adjoint_solution(i),
00545                                  *(temprhs[i]),
00546                                  solver_params.second,
00547                                  solver_params.first);
00548 
00549         totalrval.first  += rval.first;
00550         totalrval.second += rval.second;
00551       }
00552 
00553   this->release_linear_solver(linear_solver);
00554 
00555   for (unsigned int i=0; i != this->qoi.size(); ++i)
00556     if (qoi_indices.has_index(i))
00557       delete temprhs[i];
00558 
00559   // The linear solver may not have fit our constraints exactly
00560 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00561   for (unsigned int i=0; i != this->qoi.size(); ++i)
00562     if (qoi_indices.has_index(i))
00563       this->get_dof_map().enforce_constraints_exactly
00564         (*this, &this->get_weighted_sensitivity_adjoint_solution(i),
00565          /* homogeneous = */ true);
00566 #endif
00567 
00568   // Stop logging the nonlinear solve
00569   STOP_LOG("weighted_sensitivity_adjoint_solve()", "ImplicitSystem");
00570 
00571   return totalrval;
00572 }

std::pair< unsigned int, Real > libMesh::ImplicitSystem::weighted_sensitivity_solve ( const ParameterVector &  parameters,
const ParameterVector &  weights 
) [virtual, inherited]

Assembles & solves the linear system(s) (dR/du)*u_w = sum(w_p*-dR/dp), for those parameters p contained within parameters weighted by the values w_p found within weights.

Returns a pair with the total number of linear iterations performed and the (sum of the) final residual norms

Reimplemented from libMesh::System.

Definition at line 577 of file implicit_system.C.

References libMesh::System::add_weighted_sensitivity_solution(), libMesh::ImplicitSystem::assembly(), libMesh::NumericVector< T >::clone(), libMesh::SparseMatrix< T >::close(), libMesh::NumericVector< T >::close(), libMesh::ParameterVector::deep_copy(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_dof_map(), libMesh::ImplicitSystem::get_linear_solve_parameters(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::System::get_weighted_sensitivity_solution(), libMesh::ImplicitSystem::matrix, libMesh::Real, libMesh::ImplicitSystem::release_linear_solver(), libMesh::ExplicitSystem::rhs, libMesh::LinearSolver< T >::solve(), libMesh::TOLERANCE, and libMesh::ParameterVector::value_copy().

00579 {
00580   // Log how long the linear solve takes.
00581   START_LOG("weighted_sensitivity_solve()", "ImplicitSystem");
00582 
00583   // We currently get partial derivatives via central differencing
00584   const Real delta_p = TOLERANCE;
00585 
00586   // The forward system should now already be solved.
00587 
00588   // Now we're assembling a weighted sum of sensitivity systems:
00589   //
00590   // dR/du (u, v)(sum(w_l*u'_l)) = -sum_l(w_l*R'_l (u, v)) forall v
00591 
00592   // We'll assemble the rhs first, because the R' term will require
00593   // perturbing the system, and some applications may not be able to
00594   // assemble a perturbed residual without simultaneously constructing
00595   // a perturbed jacobian.
00596 
00597   // We approximate the _l partial derivatives via a central
00598   // differencing perturbation in the w_l direction:
00599   //
00600   // sum_l(w_l*v_l) ~= (v(p + dp*w_l*e_l) - v(p - dp*w_l*e_l))/(2*dp)
00601 
00602   ParameterVector oldparameters, parameterperturbation;
00603   parameters.deep_copy(oldparameters);
00604   weights.deep_copy(parameterperturbation);
00605   parameterperturbation *= delta_p;
00606   parameters += parameterperturbation;
00607 
00608   this->assembly(true, false);
00609   this->rhs->close();
00610 
00611   AutoPtr<NumericVector<Number> > temprhs = this->rhs->clone();
00612 
00613   oldparameters.value_copy(parameters);
00614   parameterperturbation *= -1.0;
00615   parameters += parameterperturbation;
00616 
00617   this->assembly(true, false);
00618   this->rhs->close();
00619 
00620   *temprhs -= *(this->rhs);
00621   *temprhs /= (2.0*delta_p);
00622 
00623   // Finally, assemble the jacobian at the non-perturbed parameter
00624   // values
00625   oldparameters.value_copy(parameters);
00626 
00627   // Build the Jacobian
00628   this->assembly(false, true);
00629   this->matrix->close();
00630 
00631   // The weighted sensitivity problem is linear
00632   LinearSolver<Number> *linear_solver = this->get_linear_solver();
00633 
00634   std::pair<unsigned int, Real> solver_params =
00635     this->get_linear_solve_parameters();
00636 
00637   const std::pair<unsigned int, Real> rval =
00638     linear_solver->solve (*matrix, this->add_weighted_sensitivity_solution(),
00639                           *temprhs,
00640                           solver_params.second,
00641                           solver_params.first);
00642 
00643   this->release_linear_solver(linear_solver);
00644 
00645   // The linear solver may not have fit our constraints exactly
00646 #ifdef LIBMESH_ENABLE_CONSTRAINTS
00647   this->get_dof_map().enforce_constraints_exactly
00648     (*this, &this->get_weighted_sensitivity_solution(),
00649      /* homogeneous = */ true);
00650 #endif
00651 
00652   // Stop logging the nonlinear solve
00653   STOP_LOG("weighted_sensitivity_solve()", "ImplicitSystem");
00654 
00655   return rval;
00656 }

void libMesh::System::write_header ( Xdr &  io,
const std::string &  version,
const bool  write_additional_data 
) const [inherited]

Writes the basic data header for this System.

This method implements the output of a System object, embedded in the output of an EquationSystems<T_sys>. This warrants some documentation. The output of this part consists of 5 sections:

for this system

5.) The number of variables in the system (unsigned int)

for each variable in the system

6.) The name of the variable (string)

6.1.) subdomain where the variable lives

7.) Combined in an FEType:

  • The approximation order(s) of the variable (Order Enum, cast to int/s)
  • The finite element family/ies of the variable (FEFamily Enum, cast to int/s)

end variable loop

8.) The number of additional vectors (unsigned int),

for each additional vector in the system object

9.) the name of the additional vector (string)

end system

Definition at line 1240 of file system_io.C.

References libMesh::System::_vectors, libMesh::Variable::active_subdomains(), libMesh::Xdr::data(), libMesh::FEType::family, libMesh::System::get_mesh(), libMesh::FEType::inf_map, libMesh::System::n_vars(), libMesh::System::n_vectors(), libMesh::System::name(), libMesh::FEType::order, libMesh::MeshBase::processor_id(), libMesh::FEType::radial_family, libMesh::FEType::radial_order, libMesh::System::variable(), libMesh::System::variable_name(), libMesh::System::variable_type(), and libMesh::Xdr::writing().

01243 {
01277   libmesh_assert (io.writing());
01278 
01279 
01280   // Only write the header information
01281   // if we are processor 0.
01282   if (this->get_mesh().processor_id() != 0)
01283     return;
01284 
01285   std::string comment;
01286   char buf[80];
01287 
01288   // 5.)
01289   // Write the number of variables in the system
01290 
01291   {
01292     // set up the comment
01293     comment = "# No. of Variables in System \"";
01294     comment += this->name();
01295     comment += "\"";
01296 
01297     unsigned int nv = this->n_vars();
01298     io.data (nv, comment.c_str());
01299   }
01300 
01301 
01302   for (unsigned int var=0; var<this->n_vars(); var++)
01303     {
01304       // 6.)
01305       // Write the name of the var-th variable
01306       {
01307         // set up the comment
01308         comment  = "#   Name, Variable No. ";
01309         std::sprintf(buf, "%d", var);
01310         comment += buf;
01311         comment += ", System \"";
01312         comment += this->name();
01313         comment += "\"";
01314 
01315         std::string var_name = this->variable_name(var);
01316         io.data (var_name, comment.c_str());
01317       }
01318 
01319       // 6.1.) Variable subdomains
01320       {
01321         // set up the comment
01322                 comment  = "#     Subdomains, Variable \"";
01323                 std::sprintf(buf, "%s", this->variable_name(var).c_str());
01324                 comment += buf;
</