#include <convection_crb.hpp>
The class derives from the Application class the template arguments are :
Order_s | velocity polynomial order |
Order_t | temperature polynomial order |
Order_p | pressure polynomial order |
Inherits ModelCrbBase< ParameterDefinition, FunctionSpaceDefinition, EimDefinition< ParameterDefinition, FunctionSpaceDefinition > >, and enable_shared_from_this< ConvectionCrb >.
Public Types | |
typedef boost::tuple < std::vector< std::vector < sparse_matrix_ptrtype > >, std::vector< std::vector < std::vector< vector_ptrtype > > > > | affine_decomposition_type |
typedef boost::shared_ptr < backend_type > | backend_ptrtype |
typedef Backend< double > | backend_type |
typedef Lagrange< Order_p, Scalar, Continuous, PointSetFekete > | basis_p_type |
typedef Lagrange< Order_t, Scalar, Continuous, PointSetFekete > | basis_t_type |
typedef bases< basis_u_type, basis_p_type, basis_t_type > | basis_type |
typedef Lagrange< Order_s, Vectorial, Continuous, PointSetFekete > | basis_u_type |
typedef std::vector < std::vector< double > > | beta_vector_type |
typedef boost::shared_ptr < eigen_matrix_type > | eigen_matrix_ptrtype |
typedef Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > | eigen_matrix_type |
typedef element_type::sub_element < 0 >::type | element_0_type |
typedef element_type::sub_element < 1 >::type | element_1_type |
typedef element_type::sub_element < 2 >::type | element_2_type |
typedef boost::shared_ptr < element_type > | element_ptrtype |
typedef space_type::element_type | element_type |
typedef eigen_matrix_type | ematrix_type |
typedef Simplex< CONVECTION_DIM > | entity_type |
typedef Exporter< mesh_type > | export_type |
typedef space_ptrtype | functionspace_ptrtype |
typedef space_type | functionspace_type |
typedef boost::shared_ptr < funlin_type > | funlin_ptrtype |
typedef FsFunctionalLinear < space_type > | funlin_type |
typedef super_type::funs_type | funs_type |
typedef super_type::funsd_type | funsd_type |
typedef Eigen::MatrixXd | matrixN_type |
typedef boost::shared_ptr < mesh_type > | mesh_ptrtype |
typedef Mesh< entity_type > | mesh_type |
typedef boost::shared_ptr < oplin_type > | oplin_ptrtype |
typedef OperatorLinear < space_type, space_type > | oplin_type |
typedef parameterspace_type::element_ptrtype | parameter_ptrtype |
typedef parameterspace_type::element_type | parameter_type |
typedef boost::shared_ptr < parameterspace_type > | parameterspace_ptrtype |
typedef ParameterSpace < ParameterSpaceDimension > | parameterspace_type |
typedef boost::shared_ptr < rbfunctionspace_type > | rbfunctionspace_ptrtype |
typedef ReducedBasisSpace < ConvectionCrb, mesh_type, basis_type, value_type > | rbfunctionspace_type |
typedef parameterspace_type::sampling_ptrtype | sampling_ptrtype |
typedef parameterspace_type::sampling_type | sampling_type |
typedef ConvectionCrb | self_type |
typedef boost::shared_ptr < space_type > | space_ptrtype |
typedef FunctionSpace < mesh_type, basis_type > | space_type |
typedef backend_type::sparse_matrix_ptrtype | sparse_matrix_ptrtype |
typedef backend_type::sparse_matrix_type | sparse_matrix_type |
typedef ModelCrbBase < ParameterDefinition, FunctionSpaceDefinition, EimDefinition < ParameterDefinition, FunctionSpaceDefinition > > | super_type |
typedef boost::shared_ptr < T_space_type > | T_space_ptrtype |
typedef FunctionSpace < mesh_type, basis_t_type > | T_space_type |
typedef boost::shared_ptr < U_space_type > | U_space_ptrtype |
typedef FunctionSpace < mesh_type, basis_u_type > | U_space_type |
typedef double | value_type |
numerical type is double | |
typedef backend_type::vector_ptrtype | vector_ptrtype |
typedef backend_type::vector_type | vector_type |
Public Member Functions | |
affine_decomposition_type | computeAffineDecomposition () |
boost::tuple< beta_vector_type, std::vector< beta_vector_type > > | computeBetaQm (parameter_type const &mu, double time=0) |
compute the beta coefficient for both bilinear and linear form More... | |
boost::tuple< beta_vector_type, std::vector< beta_vector_type > > | computeBetaQm (element_type const &T, parameter_type const &mu, double time=0) |
sparse_matrix_ptrtype | computeTrilinearForm (const element_type &X) |
ConvectionCrb (po::variables_map const &vm) | |
Feel::gmsh_ptrtype | createMesh () |
void | exportResults (element_type &u) |
void | exportResults (element_ptrtype &U, int t) |
void | exportResults (element_type &U, double t) |
space_ptrtype | functionSpace () |
void | initModel () |
sparse_matrix_ptrtype | innerProduct () |
sparse_matrix_ptrtype | jacobian (const element_type &X) |
void | l2solve (vector_ptrtype &u, vector_ptrtype const &f) |
int | mMaxA (int q) |
int | mMaxF (int output_index, int q) |
std::string | modelName () |
sparse_matrix_ptrtype | newMatrix () const |
vector_ptrtype | newVector () const |
int | Nl () const |
po::options_description const & | optionsDescription () const |
value_type | output (int output_index, parameter_type const &mu, element_type &unknown, bool need_to_solve=false) |
parameterspace_ptrtype | parameterSpace () const |
return the parameter space | |
int | Qa () const |
int | QaTri () const |
int | Ql (int l) const |
rbfunctionspace_ptrtype | rBFunctionSpace () |
Returns the reduced basis function space. | |
parameter_type | refParameter () |
vector_ptrtype | residual (const element_type &X) |
void | run (const double *X, unsigned long N, double *Y, unsigned long P) |
double | scalarProduct (vector_ptrtype const &X, vector_ptrtype const &Y) |
double | scalarProduct (vector_type const &x, vector_type const &y) |
void | setMeshSize (double s) |
void | solve (sparse_matrix_ptrtype &D, element_type &u, vector_ptrtype &F) |
void | solve (parameter_type const &mu, element_ptrtype &T) |
solve the model for parameter mu More... | |
element_type | solve (parameter_type const &mu) |
void | update (parameter_type const &mu) |
void | updateJacobian (const vector_ptrtype &X, sparse_matrix_ptrtype &J) |
void | updateJacobianWithoutAffineDecomposition (const vector_ptrtype &X, sparse_matrix_ptrtype &J) |
void | updateResidual (const vector_ptrtype &X, vector_ptrtype &R) |
po::variables_map const & | vm () const |
boost::tuple< beta_vector_type, std::vector< beta_vector_type > > ConvectionCrb::computeBetaQm | ( | parameter_type const & | mu, |
double | time = 0 |
||
) |
void ConvectionCrb::exportResults | ( | element_type & | u | ) |
export results to ensight format (enabled by –export cmd line options)
void ConvectionCrb::l2solve | ( | vector_ptrtype & | u, |
vector_ptrtype const & | f | ||
) |
solve \( M u = f \)
int ConvectionCrb::Nl | ( | ) | const |
there is at least one output which is the right hand side of the primal problem
double ConvectionCrb::output | ( | int | output_index, |
parameter_type const & | mu, | ||
element_type & | unknown, | ||
bool | need_to_solve = false |
||
) |
Given the output index output_index
and the parameter mu
, return the value of the corresponding FEM output
int ConvectionCrb::Ql | ( | int | l | ) | const |
l | the index of output |
q
th output term void ConvectionCrb::run | ( | const double * | X, |
unsigned long | N, | ||
double * | Y, | ||
unsigned long | P | ||
) |
specific interface for OpenTURNS
X | input vector of size N |
N | size of input vector X |
Y | input vector of size P |
P | size of input vector Y |
double ConvectionCrb::scalarProduct | ( | vector_ptrtype const & | X, |
vector_ptrtype const & | Y | ||
) |
returns the scalar product of the boost::shared_ptr vector x and boost::shared_ptr vector y
double ConvectionCrb::scalarProduct | ( | vector_type const & | x, |
vector_type const & | y | ||
) |
returns the scalar product of the vector x and vector y
void ConvectionCrb::solve | ( | parameter_type const & | mu, |
element_ptrtype & | T | ||
) |
solve the model for parameter mu
mu | the model parameter |
T | the temperature field |
References computeBetaQm(), exportResults(), and vm().
ConvectionCrb::element_type ConvectionCrb::solve | ( | parameter_type const & | mu | ) |
solve for a given parameter mu
|
inline |
get the variable map