|
HArD::Core2D
Hybrid Arbitrary Degree::Core 2D - Library to implement 2D schemes with edge and cell polynomials as unknowns
|
Hybrid Mimetic Mixed method. More...
Classes | |
| class | HArDCore2D::HMM_StefanPME_Transient |
| The vector Xh manipulated in the resolution has mixed components, corresponding either to the unknown u or to \(\zeta(u)\), depending on the choice of weight of mass-lumping for the cell/edge unknowns. If no weight is put on the edges (resp. the cells), then the edge (resp. cell) unknowns represent \(\zeta(u)\). Otherwise, they represent u. More... | |
| class | HArDCore2D::StochStefanPME |
| The vector Xh manipulated in the resolution has mixed components, corresponding either to the unknown \(u\) or to \(\zeta(e^W u)\), depending on the choice of weight of mass-lumping for the cell/edge unknowns. If no weight is put on the edges (resp. the cells), then the edge (resp. cell) unknowns represent \(\zeta(e^W u)\). Otherwise, they represent u. More... | |
Typedefs | |
| using | HArDCore2D::HMM_StefanPME_Transient::solution_function_type = std::function< double(const double &, const VectorRd &)> |
| type for solution | |
| using | HArDCore2D::HMM_StefanPME_Transient::source_function_type = std::function< double(const double &, const VectorRd &, const Cell *)> |
| type for source | |
| using | HArDCore2D::HMM_StefanPME_Transient::grad_function_type = std::function< VectorRd(const double &, const VectorRd &, const Cell *)> |
| type for gradient | |
| using | HArDCore2D::HMM_StefanPME_Transient::tensor_function_type = std::function< Eigen::Matrix2d(const double, const double, const Cell *)> |
| type for diffusion tensor (does not depend on time) | |
| template<typename T > | |
| using | HArDCore2D::SpatialFunctionType = std::function< T(const VectorRd &)> |
| Generic type for function that only depends on spatial coordinate. | |
| template<typename T > | |
| using | HArDCore2D::PiecewiseSpatialFunctionType = std::function< T(const VectorRd &, const Cell *)> |
| Generic type for function that depends on space, with formula changing from one cell to the next. | |
| template<typename T > | |
| using | HArDCore2D::TemporalSpatialFunctionType = std::function< T(const double &, const VectorRd &)> |
| Generic type for function that depends on space and time. | |
| template<typename T > | |
| using | HArDCore2D::PiecewiseTemporalSpatialFunctionType = std::function< T(const double &, const VectorRd &, const Cell *)> |
| Generic type for function that depends on space and time, with formula changing from one cell to the next. | |
| template<typename T > | |
| using | HArDCore2D::EigFunctionType = std::function< T(const size_t &, const size_t &, const VectorRd &)> |
| Type for "eigenfunctions" e_i, depending on two indices k,l. | |
| template<typename T > | |
| using | HArDCore2D::DoubleVector = std::vector< std::vector< T > > |
| using | HArDCore2D::StochStefanPME::solution_function_type = TemporalSpatialFunctionType< double > |
| type for solution | |
| using | HArDCore2D::StochStefanPME::source_function_type = PiecewiseSpatialFunctionType< double > |
| type for source (at a given fixed time) | |
| using | HArDCore2D::StochStefanPME::grad_function_type = PiecewiseTemporalSpatialFunctionType< VectorRd > |
| type for gradient | |
| using | HArDCore2D::StochStefanPME::lapl_function_type = PiecewiseTemporalSpatialFunctionType< double > |
| type for laplacian | |
| using | HArDCore2D::StochStefanPME::tensor_function_type = PiecewiseSpatialFunctionType< Eigen::Matrix2d > |
| type for diffusion tensor (does not depend on time) | |
| using | HArDCore2D::StochStefanPME::WienerType = SpatialFunctionType< double > |
| type for the Wiener process at a fixed time (only variation in space) | |
| using | HArDCore2D::StochStefanPME::GradWienerType = SpatialFunctionType< VectorRd > |
| type for gradient of Wiener process | |
| using | HArDCore2D::StochStefanPME::ReactionType = SpatialFunctionType< double > |
| Type for reaction term mu^2. | |
Functions | |
| HArDCore2D::HMM_StefanPME_Transient::HMM_StefanPME_Transient (HybridCore &hmm, tensor_function_type kappa, source_function_type source, BoundaryConditions BC, solution_function_type exact_solution, grad_function_type grad_exact_solution, TestCaseNonLinearity::nonlinearity_function_type zeta, double weight, std::string solver_type, std::ostream &output=std::cout) | |
| Constructor of the class. | |
| UVector | HArDCore2D::HMM_StefanPME_Transient::iterate (const double tps, const double dt, const UVector &Xn) |
| Execute one time iteration. | |
| Eigen::VectorXd | HArDCore2D::HMM_StefanPME_Transient::apply_nonlinearity (const Eigen::VectorXd &Y, const std::string type) const |
| Compute non-linearity on vector (depends if weight=0, weight=1 or weight\in (0,1) ) | |
| UVector | HArDCore2D::HMM_StefanPME_Transient::apply_nonlinearity (const UVector &Y, const std::string type) const |
| double | HArDCore2D::HMM_StefanPME_Transient::L2_MassLumped (const UVector &Xh) const |
| Mass-lumped L2 norm of a function given by a vector. | |
| double | HArDCore2D::HMM_StefanPME_Transient::Lp_MassLumped (const UVector &Xh, double p) const |
| Mass-lumped Lp norm of a function given by a vector. | |
| double | HArDCore2D::HMM_StefanPME_Transient::EnergyNorm (const UVector &Xh) const |
| Discrete energy norm (associated to the diffusion operator) | |
| double | HArDCore2D::HMM_StefanPME_Transient::get_assembly_time () const |
| cpu time to assemble the scheme | |
| double | HArDCore2D::HMM_StefanPME_Transient::get_solving_time () const |
| cpu time to solve the scheme | |
| double | HArDCore2D::HMM_StefanPME_Transient::get_itime (size_t idx) const |
| various intermediate assembly times | |
| double | HArDCore2D::HMM_StefanPME_Transient::get_solving_error () const |
| residual after solving the scheme | |
| size_t | HArDCore2D::HMM_StefanPME_Transient::get_nb_newton () const |
| number of Newton iterations | |
| Eigen::MatrixXd | HArDCore2D::HMM_StefanPME_Transient::get_MassT (size_t iT) const |
| Mass matrix in cell iT. | |
| static double | HArDCore2D::StDev (const Eigen::VectorXd &v) |
| HArDCore2D::StochStefanPME::StochStefanPME (HybridCore &hmm, tensor_function_type kappa, TestCaseNonLinearity::nonlinearity_function_type zeta, BoundaryConditions BC, solution_function_type exact_solution, grad_function_type grad_exact_solution, solution_function_type time_der_exact_solution, lapl_function_type minus_Lapl_exact_solution, double weight, std::string solver_type, std::ostream &output=std::cout) | |
| Constructor of the class. | |
| StochStefanPME::source_function_type | HArDCore2D::StochStefanPME::compute_source (const bool &use_exact_source, const double &t, const WienerType &W, const GradWienerType &grad_W, const WienerType &Delta_exp_W, const ReactionType &mu_squared) |
| Compute the source to zero or the exact one based on the selected solution and Wiener process. | |
| std::pair< UVector, size_t > | HArDCore2D::StochStefanPME::iterate (const double tps, const double dt, const UVector &Xn, const UVector &Ih_W, const WienerType &W, const ReactionType &mu_squared, const source_function_type &source) |
| Execute one time iteration: return the calculated vector, and the number of Newton iterations. | |
| void | HArDCore2D::StochStefanPME::SetWienerProcess (const int &caseW, const double &t, WienerType &W, UVector &I_W, GradWienerType &grad_W, WienerType &Delta_exp_W, ReactionType &mu_squared, bool verbose=false) |
| Set the Wiener process, deterministic case. | |
| void | HArDCore2D::StochStefanPME::SimulateWienerProcess (const double &dt, DoubleVector< double > &Wtime, WienerType &W, const DoubleVector< Eigen::VectorXd > &Ih_eigs, UVector &I_W, GradWienerType &grad_W, WienerType &Delta_exp_W, ReactionType &mu_squared) |
| Set the Wiener process, random case. | |
| Eigen::VectorXd | HArDCore2D::StochStefanPME::apply_nonlinearity_eW (const std::string type, const Eigen::VectorXd &Y, const Eigen::VectorXd &I_W) const |
| Compute non-linearity and e^W on vector Y (depends if weight=0, weight=1 or weight\in (0,1) ) | |
| UVector | HArDCore2D::StochStefanPME::apply_nonlinearity_eW (const std::string type, const UVector &Y, const Eigen::VectorXd &I_W) const |
| double | HArDCore2D::StochStefanPME::L2_MassLumped (const UVector &Xh) const |
| Mass-lumped L2 norm of a function given by a vector. | |
| double | HArDCore2D::StochStefanPME::Lp_MassLumped (const UVector &Xh, double p) const |
| Mass-lumped Lp norm of a function given by a vector. | |
| double | HArDCore2D::StochStefanPME::EnergyNorm (const UVector &Xh) const |
| Discrete energy norm (associated to the diffusion operator) | |
| double | HArDCore2D::StochStefanPME::Integral (const UVector &Xh) const |
| Integral of a function given by a vector. | |
| double | HArDCore2D::StochStefanPME::get_assembly_time () const |
| cpu time to assemble the scheme | |
| double | HArDCore2D::StochStefanPME::get_solution_time () const |
| cpu time to solve the linear systems | |
| double | HArDCore2D::StochStefanPME::get_itime (size_t idx) const |
| various intermediate assembly times | |
| double | HArDCore2D::StochStefanPME::get_solving_error () const |
| residual after solving the scheme | |
| Eigen::MatrixXd | HArDCore2D::StochStefanPME::get_MassT (size_t iT) const |
| Mass matrix in cell iT. | |
Variables | |
| static SpatialFunctionType< double > | HArDCore2D::zero_scalar_function = [](const VectorRd & x)->double { return 0.;} |
| Zero spatial function. | |
| static const double | HArDCore2D::PI = boost::math::constants::pi<double>() |
| Type for vector of vector. | |
| static std::function< double(const size_t &, const size_t &)> | HArDCore2D::mui |
| static EigFunctionType< double > | HArDCore2D::mui_ei |
| static EigFunctionType< VectorRd > | HArDCore2D::grad_mui_ei |
| static EigFunctionType< double > | HArDCore2D::Delta_mui_ei |
Hybrid Mimetic Mixed method.
| using HArDCore2D::DoubleVector = typedef std::vector<std::vector<T> > |
| using HArDCore2D::EigFunctionType = typedef std::function<T(const size_t&, const size_t&, const VectorRd&)> |
Type for "eigenfunctions" e_i, depending on two indices k,l.
| using HArDCore2D::HMM_StefanPME_Transient::grad_function_type = std::function<VectorRd(const double&, const VectorRd&, const Cell*)> |
type for gradient
| using HArDCore2D::StochStefanPME::grad_function_type = PiecewiseTemporalSpatialFunctionType<VectorRd> |
type for gradient
type for gradient of Wiener process
type for laplacian
| using HArDCore2D::PiecewiseSpatialFunctionType = typedef std::function<T(const VectorRd&, const Cell*)> |
Generic type for function that depends on space, with formula changing from one cell to the next.
| using HArDCore2D::PiecewiseTemporalSpatialFunctionType = typedef std::function<T(const double&, const VectorRd&, const Cell*)> |
Generic type for function that depends on space and time, with formula changing from one cell to the next.
| using HArDCore2D::StochStefanPME::ReactionType = SpatialFunctionType<double> |
Type for reaction term mu^2.
| using HArDCore2D::HMM_StefanPME_Transient::solution_function_type = std::function<double(const double&, const VectorRd&)> |
type for solution
type for solution
| using HArDCore2D::HMM_StefanPME_Transient::source_function_type = std::function<double(const double&, const VectorRd&, const Cell*)> |
type for source
type for source (at a given fixed time)
| using HArDCore2D::SpatialFunctionType = typedef std::function<T(const VectorRd&)> |
Generic type for function that only depends on spatial coordinate.
| using HArDCore2D::TemporalSpatialFunctionType = typedef std::function<T(const double&, const VectorRd&)> |
Generic type for function that depends on space and time.
| using HArDCore2D::HMM_StefanPME_Transient::tensor_function_type = std::function<Eigen::Matrix2d(const double, const double, const Cell*)> |
type for diffusion tensor (does not depend on time)
| using HArDCore2D::StochStefanPME::tensor_function_type = PiecewiseSpatialFunctionType<Eigen::Matrix2d> |
type for diffusion tensor (does not depend on time)
| using HArDCore2D::StochStefanPME::WienerType = SpatialFunctionType<double> |
type for the Wiener process at a fixed time (only variation in space)
| Eigen::VectorXd HArDCore2D::HMM_StefanPME_Transient::apply_nonlinearity | ( | const Eigen::VectorXd & | Y, |
| const std::string | type | ||
| ) | const |
Compute non-linearity on vector (depends if weight=0, weight=1 or weight\in (0,1) )
| UVector HArDCore2D::HMM_StefanPME_Transient::apply_nonlinearity | ( | const UVector & | Y, |
| const std::string | type | ||
| ) | const |
| Eigen::VectorXd HArDCore2D::StochStefanPME::apply_nonlinearity_eW | ( | const std::string | type, |
| const Eigen::VectorXd & | Y, | ||
| const Eigen::VectorXd & | I_W | ||
| ) | const |
Compute non-linearity and e^W on vector Y (depends if weight=0, weight=1 or weight\in (0,1) )
| UVector HArDCore2D::StochStefanPME::apply_nonlinearity_eW | ( | const std::string | type, |
| const UVector & | Y, | ||
| const Eigen::VectorXd & | I_W | ||
| ) | const |
| StochStefanPME::source_function_type HArDCore2D::StochStefanPME::compute_source | ( | const bool & | use_exact_source, |
| const double & | t, | ||
| const WienerType & | W, | ||
| const GradWienerType & | grad_W, | ||
| const WienerType & | Delta_exp_W, | ||
| const ReactionType & | mu_squared | ||
| ) |
Compute the source to zero or the exact one based on the selected solution and Wiener process.
| use_exact_source | use exact source or zero |
| t | Time at which to compute the source term |
| W | Wiener process |
| grad_W | Gradient of Wiener process |
| Delta_exp_W | Laplacian of exponential of Wiener process |
| mu_squared | mu^2 for the reaction term |
| double HArDCore2D::HMM_StefanPME_Transient::EnergyNorm | ( | const UVector & | Xh | ) | const |
Discrete energy norm (associated to the diffusion operator)
| double HArDCore2D::StochStefanPME::EnergyNorm | ( | const UVector & | Xh | ) | const |
Discrete energy norm (associated to the diffusion operator)
|
inline |
cpu time to assemble the scheme
|
inline |
cpu time to assemble the scheme
|
inline |
various intermediate assembly times
|
inline |
various intermediate assembly times
|
inline |
Mass matrix in cell iT.
|
inline |
Mass matrix in cell iT.
|
inline |
number of Newton iterations
|
inline |
cpu time to solve the linear systems
|
inline |
residual after solving the scheme
|
inline |
residual after solving the scheme
|
inline |
cpu time to solve the scheme
| HArDCore2D::HMM_StefanPME_Transient::HMM_StefanPME_Transient | ( | HybridCore & | hmm, |
| tensor_function_type | kappa, | ||
| source_function_type | source, | ||
| BoundaryConditions | BC, | ||
| solution_function_type | exact_solution, | ||
| grad_function_type | grad_exact_solution, | ||
| TestCaseNonLinearity::nonlinearity_function_type | zeta, | ||
| double | weight, | ||
| std::string | solver_type, | ||
| std::ostream & | output = std::cout |
||
| ) |
Constructor of the class.
| hmm | instance of the hybridcore class (with basis functions, etc.) |
| kappa | diffusion tensor |
| source | source term |
| BC | type of boundary conditions |
| exact_solution | exact solution |
| grad_exact_solution | gradient of the exact solution |
| zeta | function describing the nonlinearity |
| weight | proportion of weight (between 0 and 1) of mass-lumping on the edges |
| solver_type | type of solver to use for the global system (bicgstab at the moment) |
| output | optional argument for output of messages |
| double HArDCore2D::StochStefanPME::Integral | ( | const UVector & | Xh | ) | const |
Integral of a function given by a vector.
| UVector HArDCore2D::HMM_StefanPME_Transient::iterate | ( | const double | tps, |
| const double | dt, | ||
| const UVector & | Xn | ||
| ) |
Execute one time iteration.
| tps | current time (for computing source) |
| dt | time step |
| Xn | Unkown at previous time step |
| std::pair< UVector, size_t > HArDCore2D::StochStefanPME::iterate | ( | const double | tps, |
| const double | dt, | ||
| const UVector & | Xn, | ||
| const UVector & | Ih_W, | ||
| const WienerType & | W, | ||
| const ReactionType & | mu_squared, | ||
| const source_function_type & | source | ||
| ) |
Execute one time iteration: return the calculated vector, and the number of Newton iterations.
| tps | current time (for computing source) |
| dt | time step |
| Xn | Unkown at previous time step |
| Ih_W | Interpolate of Wiener process at time tps |
| W | Wiener process at time tps, to fix the BCs |
| mu_squared | Reaction term |
| source | source term at time tps |
| double HArDCore2D::HMM_StefanPME_Transient::L2_MassLumped | ( | const UVector & | Xh | ) | const |
Mass-lumped L2 norm of a function given by a vector.
| double HArDCore2D::StochStefanPME::L2_MassLumped | ( | const UVector & | Xh | ) | const |
Mass-lumped L2 norm of a function given by a vector.
| double HArDCore2D::HMM_StefanPME_Transient::Lp_MassLumped | ( | const UVector & | Xh, |
| double | p | ||
| ) | const |
Mass-lumped Lp norm of a function given by a vector.
| double HArDCore2D::StochStefanPME::Lp_MassLumped | ( | const UVector & | Xh, |
| double | p | ||
| ) | const |
Mass-lumped Lp norm of a function given by a vector.
| void HArDCore2D::StochStefanPME::SetWienerProcess | ( | const int & | caseW, |
| const double & | t, | ||
| WienerType & | W, | ||
| UVector & | I_W, | ||
| GradWienerType & | grad_W, | ||
| WienerType & | Delta_exp_W, | ||
| ReactionType & | mu_squared, | ||
| bool | verbose = false |
||
| ) |
Set the Wiener process, deterministic case.
| caseW | Type of Wiener process (can be deterministic...) |
| t | Time at which to compute the Wiener process |
| W | To store the Wiener process |
| I_W | Interpolate of the calculated Wiener process |
| grad_W | To store its gradient |
| Delta_exp_W | To store the Laplacian of its exponential |
| mu_squared | To store mu^2 |
| verbose | Determine if text should be output |
| void HArDCore2D::StochStefanPME::SimulateWienerProcess | ( | const double & | dt, |
| DoubleVector< double > & | Wtime, | ||
| WienerType & | W, | ||
| const DoubleVector< Eigen::VectorXd > & | Ih_eigs, | ||
| UVector & | I_W, | ||
| GradWienerType & | grad_W, | ||
| WienerType & | Delta_exp_W, | ||
| ReactionType & | mu_squared | ||
| ) |
Set the Wiener process, random case.
Simulate the Wiener process, by random increment of the temporal component from the current value.
| dt | Time increment to update the process |
| Wtime | Current values of the time component W_i of W, will be updated |
| W | Computed complete Wiener process |
| Ih_eigs | Interpolates of the eigenfunctions e_i to compute the interpolate of the Wiener process |
| I_W | To store the interpolate of the Wiener process |
| grad_W | To store the gradient of the Wiener process |
| Delta_exp_W | To store the Laplacian of the exponential of the Wiener process |
| mu_squared | Reaction term mu^2 |
|
inlinestatic |
| HArDCore2D::StochStefanPME::StochStefanPME | ( | HybridCore & | hmm, |
| tensor_function_type | kappa, | ||
| TestCaseNonLinearity::nonlinearity_function_type | zeta, | ||
| BoundaryConditions | BC, | ||
| solution_function_type | exact_solution, | ||
| grad_function_type | grad_exact_solution, | ||
| solution_function_type | time_der_exact_solution, | ||
| lapl_function_type | minus_Lapl_exact_solution, | ||
| double | weight, | ||
| std::string | solver_type, | ||
| std::ostream & | output = std::cout |
||
| ) |
Constructor of the class.
| hmm | instance of the hybridcore class (with basis functions, etc.) |
| kappa | diffusion tensor |
| zeta | function describing the nonlinearity |
| BC | type of boundary conditions |
| exact_solution | exact solution |
| grad_exact_solution | gradient of the exact solution |
| time_der_exact_solution | time derivative of exact solution |
| minus_Lapl_exact_solution | minus Laplacian of the exact solution |
| weight | proportion of weight (between 0 and 1) of mass-lumping on the edges |
| solver_type | type of solver to use for the global system |
| output | optional argument for output of messages |
|
static |
|
static |
|
static |
|
static |
|
static |
Type for vector of vector.
|
static |
Zero spatial function.