Feel++
0.92.0
|
Certified Reduced Basis Model class. More...
#include <crbmodel.hpp>
Public Types | |
Typedefs | |
typedef ModelType | model_type |
model type | |
typedef boost::shared_ptr < ModelType > | model_ptrtype |
typedef model_type::value_type | value_type |
value_type | |
typedef ModelType::mesh_type | mesh_type |
mesh type | |
typedef ModelType::mesh_ptrtype | mesh_ptrtype |
mesh shared_ptr | |
typedef ModelType::space_type | space_type |
space_type | |
typedef model_type::functionspace_type | functionspace_type |
function space type | |
typedef model_type::functionspace_ptrtype | functionspace_ptrtype |
typedef model_type::element_type | element_type |
element of the functionspace type | |
typedef model_type::element_ptrtype | element_ptrtype |
typedef model_type::backend_type | backend_type |
typedef boost::shared_ptr < backend_type > | backend_ptrtype |
typedef model_type::sparse_matrix_ptrtype | sparse_matrix_ptrtype |
typedef model_type::vector_ptrtype | vector_ptrtype |
typedef model_type::vector_type | vector_type |
typedef model_type::eigen_matrix_type | eigen_matrix_type |
typedef model_type::parameterspace_type | parameterspace_type |
typedef model_type::parameterspace_ptrtype | parameterspace_ptrtype |
typedef model_type::parameter_type | parameter_type |
typedef model_type::parameter_ptrtype | parameter_ptrtype |
typedef model_type::sampling_type | sampling_type |
typedef model_type::sampling_ptrtype | sampling_ptrtype |
typedef Eigen::VectorXd | theta_vector_type |
typedef boost::tuple < sparse_matrix_ptrtype, sparse_matrix_ptrtype, std::vector< vector_ptrtype > > | offline_merge_type |
typedef boost::tuple < std::vector < sparse_matrix_ptrtype > , std::vector < sparse_matrix_ptrtype > , std::vector< std::vector < vector_ptrtype > > > | affine_decomposition_type |
typedef boost::tuple < theta_vector_type, theta_vector_type, std::vector < theta_vector_type > > | thetaq_type |
Public Member Functions | |
Constructors, destructor | |
CRBModel () | |
CRBModel (po::variables_map const &vm, CRBModelMode mode=CRBModelMode::PFEM) | |
CRBModel (model_ptrtype &model) | |
CRBModel (CRBModel const &o) | |
virtual | ~CRBModel () |
destructor | |
FEELPP_DONT_INLINE void | init () |
initialize the model (mesh, function space, operators, matrices, ...) | |
Operator overloads | |
CRBModel & | operator= (CRBModel const &o) |
copy operator | |
Accessors | |
po::variables_map | vm () const |
virtual sparse_matrix_ptrtype | newMatrix () const |
sparse_matrix_ptrtype const & | innerProduct () const |
Returns the matrix associated with the ![]() | |
sparse_matrix_ptrtype const & | h1 () const |
Returns the matrix associated with the ![]() | |
functionspace_ptrtype | functionSpace () const |
Returns the function space. | |
size_type | Qa () const |
return the number of ![]() | |
size_type | Qm () const |
return the number of ![]() | |
size_type | Qm (mpl::bool_< true >) const |
size_type | Qm (mpl::bool_< false >) const |
size_type | Nl () const |
return the number of outputs | |
size_type | Ql (int l) const |
return the number of ![]() | |
parameterspace_ptrtype | parameterSpace () const |
return the parameter space | |
Mutators | |
void | setMeshSize (double s) |
Methods | |
thetaq_type | computeThetaq (parameter_type const &mu, double time=0) |
compute the thetaq given mu | |
thetaq_type | computeThetaq (parameter_type const &mu, mpl::bool_< true >, double time=0) |
thetaq_type | computeThetaq (parameter_type const &mu, mpl::bool_< false >, double time=0) |
offline_merge_type | update (parameter_type const &mu, double time=0) |
update the model wrt mu | |
affine_decomposition_type | computeAffineDecomposition () |
Compute the affine decomposition of the various forms. | |
affine_decomposition_type | computeAffineDecomposition (mpl::bool_< true >) |
affine_decomposition_type | computeAffineDecomposition (mpl::bool_< false >) |
value_type | h1 (element_type const &xi_i, element_type const &xi_j) |
the inner product ![]() | |
value_type | h1 (element_type const &xi_i) |
the inner product ![]() | |
sparse_matrix_ptrtype | Aq (uint16_type q, bool transpose=false) const |
Returns the matrix Aq [q] of the affine decomposition of the bilinear form. | |
sparse_matrix_ptrtype | Mq (uint16_type q, bool transpose=false) const |
Returns the matrix Mq [q] of the affine decomposition of the bilinear form (time dependent) | |
value_type | Aq (uint16_type q, element_type const &xi_i, element_type const &xi_j, bool transpose=false) |
the inner product ![]() | |
value_type | Mq (uint16_type q, element_type const &xi_i, element_type const &xi_j, bool transpose=false) |
the inner product ![]() | |
theta_vector_type const & | thetaAq () const |
Returns the vector coefficients. | |
theta_vector_type const & | thetaMq () const |
Returns the vector coefficients. | |
theta_vector_type const & | thetaMq (mpl::bool_< true >) const |
theta_vector_type const & | thetaMq (mpl::bool_< false >) const |
value_type | thetaAq (int q) const |
Returns the value of the ![]() ![]() | |
value_type | thetaMq (int q) const |
Returns the value of the ![]() ![]() | |
value_type | thetaMq (int q, mpl::bool_< true >) const |
value_type | thetaMq (int q, mpl::bool_< false >) const |
std::vector< theta_vector_type > const & | thetaL () const |
Returns the vector coefficients. | |
value_type | thetaL (int l, int q) const |
Returns the value of the ![]() ![]() | |
vector_ptrtype | Fq (uint16_type l, uint16_type q) const |
the vector Fq [q] of the affine decomposition of the right hand side | |
value_type | Fq (uint16_type l, uint16_type q, element_type const &xi) |
the inner product ![]() | |
value_type | Fq (uint16_type l, uint16_type q, element_ptrtype const &xi) |
the inner product ![]() | |
double | scalarProduct (vector_type const &X, vector_type const &Y) |
double | scalarProduct (vector_ptrtype const &X, vector_ptrtype const &Y) |
double | scalarProductForPod (vector_type const &X, vector_type const &Y) |
double | scalarProductForPod (vector_type const &X, vector_type const &Y, mpl::bool_< true >) |
double | scalarProductForPod (vector_type const &X, vector_type const &Y, mpl::bool_< false >) |
double | scalarProductForPod (vector_ptrtype const &X, vector_ptrtype const &Y) |
double | scalarProductForPod (vector_ptrtype const &X, vector_ptrtype const &Y, mpl::bool_< true >) |
double | scalarProductForPod (vector_ptrtype const &X, vector_ptrtype const &Y, mpl::bool_< false >) |
void | solve (parameter_type const &mu) |
void | solve (parameter_type const &mu, element_ptrtype &u) |
void | solve (parameter_type const &mu, element_ptrtype &u, vector_ptrtype const &L, bool transpose=false) |
void | l2solve (vector_ptrtype &u, vector_ptrtype const &f) |
boost::tuple< int, value_type > | solve (sparse_matrix_ptrtype const &A, vector_ptrtype &x, vector_ptrtype const &f, bool tranpose=false) |
solve ![]() | |
void | run (const double *X, unsigned long N, double *Y, unsigned long P) |
export a vector of elements | |
std::map< double, boost::tuple < double, double, double, int, double, double > > | run () |
std::map< double, boost::tuple < double, double, double, int, double, double > > | run (mpl::bool_< true >) |
std::map< double, boost::tuple < double, double, double, int, double, double > > | run (mpl::bool_< false >) |
std::map< double, boost::tuple < double, double, double, int, double, double > > | runq () |
value_type | output (int output_index, parameter_type const &mu) |
int | computeNumberOfSnapshots () |
int | computeNumberOfSnapshots (mpl::bool_< true >) |
int | computeNumberOfSnapshots (mpl::bool_< false >) |
double | timeStep () |
double | timeStep (mpl::bool_< true >) |
double | timeStep (mpl::bool_< false >) |
double | timeInitial () |
double | timeInitial (mpl::bool_< true >) |
double | timeInitial (mpl::bool_< false >) |
double | timeFinal () |
double | timeFinal (mpl::bool_< true >) |
double | timeFinal (mpl::bool_< false >) |
int | timeOrder () |
int | timeOrder (mpl::bool_< true >) |
int | timeOrder (mpl::bool_< false >) |
bool | isSteady () |
bool | isSteady (mpl::bool_< true >) |
bool | isSteady (mpl::bool_< false >) |
void | initializationField (element_ptrtype &initial_field, parameter_type const &mu) |
void | initializationField (element_ptrtype &initial_field, parameter_type const &mu, mpl::bool_< true >) |
void | initializationField (element_ptrtype &initial_field, parameter_type const &mu, mpl::bool_< false >) |
Static Public Attributes | |
Constants | |
static const uint16_type | ParameterSpaceDimension = ModelType::ParameterSpaceDimension |
static const bool | is_time_dependent = ModelType::is_time_dependent |
Protected Attributes | |
std::vector < sparse_matrix_ptrtype > | M_Aq |
affine decomposition terms for the left hand side | |
std::vector < sparse_matrix_ptrtype > | M_Mq |
affine decomposition terms ( time dependent ) | |
std::vector< std::vector < vector_ptrtype > > | M_Fq |
affine decomposition terms for the right hand side |
Certified Reduced Basis Model class.
This class implements the requirements over a model to be usable by the certified reduced basis method.
ModelType | the type of the finite element method model |
The FEM model type should derive from this class and fill the vector of matrices M_Aq and vector of vectors M_Fq
Feel::CRBModel< ModelType >::CRBModel | ( | model_ptrtype & | model | ) | [inline] |
model | the model to be used |
References Feel::CRBModel< ModelType >::init().
Feel::CRBModel< ModelType >::CRBModel | ( | CRBModel< ModelType > const & | o | ) | [inline] |
copy constructor
References Feel::CRBModel< ModelType >::init().
sparse_matrix_ptrtype Feel::CRBModel< ModelType >::Aq | ( | uint16_type | q, |
bool | transpose = false |
||
) | const [inline] |
Returns the matrix Aq
[q] of the affine decomposition of the bilinear form.
q | the index of the component in the affine decomposition |
transpose | transpose A_q |
Aq
[q] of the affine decomposition of the bilinear form References Feel::CRBModel< ModelType >::M_Aq.
value_type Feel::CRBModel< ModelType >::Aq | ( | uint16_type | q, |
element_type const & | xi_i, | ||
element_type const & | xi_j, | ||
bool | transpose = false |
||
) | [inline] |
the inner product
q | the index of the component in the affine decomposition |
xi_i | an element of the function space |
xi_j | an element of the function space |
transpose | transpose A_q |
References Feel::CRBModel< ModelType >::M_Aq.
affine_decomposition_type Feel::CRBModel< ModelType >::computeAffineDecomposition | ( | ) | [inline] |
Compute the affine decomposition of the various forms.
This function defined in the M_model
assembles the parameter independant part of the affine decomposition of the bilinear and linear forms.
vector_ptrtype Feel::CRBModel< ModelType >::Fq | ( | uint16_type | l, |
uint16_type | q | ||
) | const [inline] |
the vector Fq
[q] of the affine decomposition of the right hand side
References Feel::CRBModel< ModelType >::M_Fq.
value_type Feel::CRBModel< ModelType >::Fq | ( | uint16_type | l, |
uint16_type | q, | ||
element_type const & | xi | ||
) | [inline] |
the inner product
Denote the algebraic representation of the linear form associated with the right hand side.
q | the index of the component in the affine decomposition |
xi | an element of the function space |
References Feel::inner_product(), and Feel::CRBModel< ModelType >::M_Fq.
value_type Feel::CRBModel< ModelType >::Fq | ( | uint16_type | l, |
uint16_type | q, | ||
element_ptrtype const & | xi | ||
) | [inline] |
the inner product
Denote the algebraic representation of the linear form associated with the right hand side.
q | the index of the component in the affine decomposition |
xi | an element of the function space |
References Feel::inner_product(), and Feel::CRBModel< ModelType >::M_Fq.
value_type Feel::CRBModel< ModelType >::h1 | ( | element_type const & | xi_i, |
element_type const & | xi_j | ||
) | [inline] |
the inner product
xi_i | an element of the function space |
xi_j | an element of the function space |
transpose | transpose A_q |
value_type Feel::CRBModel< ModelType >::h1 | ( | element_type const & | xi_i | ) | [inline] |
the inner product
xi_i | an element of the function space |
xi_j | an element of the function space |
transpose | transpose A_q |
void Feel::CRBModel< ModelType >::l2solve | ( | vector_ptrtype & | u, |
vector_ptrtype const & | f | ||
) | [inline] |
solve where
is the matrix associated to the
norm
sparse_matrix_ptrtype Feel::CRBModel< ModelType >::Mq | ( | uint16_type | q, |
bool | transpose = false |
||
) | const [inline] |
Returns the matrix Mq
[q] of the affine decomposition of the bilinear form (time dependent)
q | the index of the component in the affine decomposition |
transpose | transpose M_q |
Mq
[q] of the affine decomposition of the bilinear form (ime dependent) References Feel::CRBModel< ModelType >::M_Mq.
value_type Feel::CRBModel< ModelType >::Mq | ( | uint16_type | q, |
element_type const & | xi_i, | ||
element_type const & | xi_j, | ||
bool | transpose = false |
||
) | [inline] |
the inner product
q | the index of the component in the affine decomposition |
xi_i | an element of the function space |
xi_j | an element of the function space |
transpose | transpose M_q |
References Feel::CRBModel< ModelType >::M_Mq.
virtual sparse_matrix_ptrtype Feel::CRBModel< ModelType >::newMatrix | ( | ) | const [inline, virtual] |
create a new matrix
value_type Feel::CRBModel< ModelType >::output | ( | int | output_index, |
parameter_type const & | mu | ||
) | [inline] |
Given the output index output_index
and the parameter mu
, return the value of the corresponding FEM output
void Feel::CRBModel< TruthModelType >::run | ( | const double * | X, |
unsigned long | N, | ||
double * | Y, | ||
unsigned long | P | ||
) |
export a vector of elements
v | a vector of shared_ptr<> elements of the functions space run the model |
References Feel::CRBModel< ModelType >::run().
Referenced by Feel::CRBModel< ModelType >::run().
double Feel::CRBModel< ModelType >::scalarProduct | ( | vector_type const & | X, |
vector_type const & | Y | ||
) | [inline] |
returns the scalar product of the vector x and vector y
double Feel::CRBModel< ModelType >::scalarProduct | ( | vector_ptrtype const & | X, |
vector_ptrtype const & | Y | ||
) | [inline] |
returns the scalar product of the vector x and vector y
double Feel::CRBModel< ModelType >::scalarProductForPod | ( | vector_type const & | X, |
vector_type const & | Y | ||
) | [inline] |
returns the scalar product used to assemble POD matrix of the vector x and vector y
Referenced by Feel::CRBModel< ModelType >::scalarProductForPod().
double Feel::CRBModel< ModelType >::scalarProductForPod | ( | vector_ptrtype const & | X, |
vector_ptrtype const & | Y | ||
) | [inline] |
returns the scalar product used to assemble POD matrix of the vector x and vector y
References Feel::CRBModel< ModelType >::scalarProductForPod().
void Feel::CRBModel< ModelType >::setMeshSize | ( | double | s | ) | [inline] |
set the mesh characteristic length to s
void Feel::CRBModel< ModelType >::solve | ( | parameter_type const & | mu | ) | [inline] |
solve the model for a given parameter mu
void Feel::CRBModel< ModelType >::solve | ( | parameter_type const & | mu, |
element_ptrtype & | u | ||
) | [inline] |
solve the model for a given parameter mu
void Feel::CRBModel< ModelType >::solve | ( | parameter_type const & | mu, |
element_ptrtype & | u, | ||
vector_ptrtype const & | L, | ||
bool | transpose = false |
||
) | [inline] |
solve the model for a given parameter mu
and L
as right hand side
transpose | if true solve the transposed(dual) problem |
boost::tuple<int, value_type> Feel::CRBModel< ModelType >::solve | ( | sparse_matrix_ptrtype const & | A, |
vector_ptrtype & | x, | ||
vector_ptrtype const & | f, | ||
bool | tranpose = false |
||
) | [inline] |
solve
tranpose
is true then solve for A | matrix |
x | solution vector |
f | right hand side vector |
transpose | if is true solve for ![]() ![]() |
po::variables_map Feel::CRBModel< ModelType >::vm | ( | ) | const [inline] |
variables_map