HArD::Core2D
Hybrid Arbitrary Degree::Core 2D - Library to implement 2D schemes with edge and cell polynomials as unknowns
Classes | Typedefs | Functions | Variables
HMM

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 More...
 
using HArDCore2D::HMM_StefanPME_Transient::source_function_type = std::function< double(const double &, const VectorRd &, const Cell *)>
 type for source More...
 
using HArDCore2D::HMM_StefanPME_Transient::grad_function_type = std::function< VectorRd(const double &, const VectorRd &, const Cell *)>
 type for gradient More...
 
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) More...
 
template<typename T >
using HArDCore2D::SpatialFunctionType = std::function< T(const VectorRd &)>
 Generic type for function that only depends on spatial coordinate. More...
 
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. More...
 
template<typename T >
using HArDCore2D::TemporalSpatialFunctionType = std::function< T(const double &, const VectorRd &)>
 Generic type for function that depends on space and time. More...
 
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. More...
 
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. More...
 
template<typename T >
using HArDCore2D::DoubleVector = std::vector< std::vector< T > >
 
using HArDCore2D::StochStefanPME::solution_function_type = TemporalSpatialFunctionType< double >
 type for solution More...
 
using HArDCore2D::StochStefanPME::source_function_type = PiecewiseSpatialFunctionType< double >
 type for source (at a given fixed time) More...
 
using HArDCore2D::StochStefanPME::grad_function_type = PiecewiseTemporalSpatialFunctionType< VectorRd >
 type for gradient More...
 
using HArDCore2D::StochStefanPME::lapl_function_type = PiecewiseTemporalSpatialFunctionType< double >
 type for laplacian More...
 
using HArDCore2D::StochStefanPME::tensor_function_type = PiecewiseSpatialFunctionType< Eigen::Matrix2d >
 type for diffusion tensor (does not depend on time) More...
 
using HArDCore2D::StochStefanPME::WienerType = SpatialFunctionType< double >
 type for the Wiener process at a fixed time (only variation in space) More...
 
using HArDCore2D::StochStefanPME::GradWienerType = SpatialFunctionType< VectorRd >
 type for gradient of Wiener process More...
 
using HArDCore2D::StochStefanPME::ReactionType = SpatialFunctionType< double >
 Type for reaction term mu^2. More...
 

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. More...
 
UVector HArDCore2D::HMM_StefanPME_Transient::iterate (const double tps, const double dt, const UVector &Xn)
 Execute one time iteration. More...
 
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) ) More...
 
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. More...
 
double HArDCore2D::HMM_StefanPME_Transient::Lp_MassLumped (const UVector &Xh, double p) const
 Mass-lumped Lp norm of a function given by a vector. More...
 
double HArDCore2D::HMM_StefanPME_Transient::EnergyNorm (const UVector &Xh) const
 Discrete energy norm (associated to the diffusion operator) More...
 
double HArDCore2D::HMM_StefanPME_Transient::get_assembly_time () const
 cpu time to assemble the scheme More...
 
double HArDCore2D::HMM_StefanPME_Transient::get_solving_time () const
 cpu time to solve the scheme More...
 
double HArDCore2D::HMM_StefanPME_Transient::get_itime (size_t idx) const
 various intermediate assembly times More...
 
double HArDCore2D::HMM_StefanPME_Transient::get_solving_error () const
 residual after solving the scheme More...
 
size_t HArDCore2D::HMM_StefanPME_Transient::get_nb_newton () const
 number of Newton iterations More...
 
Eigen::MatrixXd HArDCore2D::HMM_StefanPME_Transient::get_MassT (size_t iT) const
 Mass matrix in cell iT. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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) ) More...
 
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. More...
 
double HArDCore2D::StochStefanPME::Lp_MassLumped (const UVector &Xh, double p) const
 Mass-lumped Lp norm of a function given by a vector. More...
 
double HArDCore2D::StochStefanPME::EnergyNorm (const UVector &Xh) const
 Discrete energy norm (associated to the diffusion operator) More...
 
double HArDCore2D::StochStefanPME::Integral (const UVector &Xh) const
 Integral of a function given by a vector. More...
 
double HArDCore2D::StochStefanPME::get_assembly_time () const
 cpu time to assemble the scheme More...
 
double HArDCore2D::StochStefanPME::get_solution_time () const
 cpu time to solve the linear systems More...
 
double HArDCore2D::StochStefanPME::get_itime (size_t idx) const
 various intermediate assembly times More...
 
double HArDCore2D::StochStefanPME::get_solving_error () const
 residual after solving the scheme More...
 
Eigen::MatrixXd HArDCore2D::StochStefanPME::get_MassT (size_t iT) const
 Mass matrix in cell iT. More...
 

Variables

static SpatialFunctionType< double > HArDCore2D::zero_scalar_function = [](const VectorRd & x)->double { return 0.;}
 Zero spatial function. More...
 
static const double HArDCore2D::PI = boost::math::constants::pi<double>()
 Type for vector of vector. More...
 
static std::function< double(const size_t &, const size_t &)> HArDCore2D::mui
 
static EigFunctionType< double > HArDCore2D::mui_ei
 
static EigFunctionType< VectorRdHArDCore2D::grad_mui_ei
 
static EigFunctionType< double > HArDCore2D::Delta_mui_ei
 

Detailed Description

Hybrid Mimetic Mixed method.

Typedef Documentation

◆ DoubleVector

template<typename T >
using HArDCore2D::DoubleVector = typedef std::vector<std::vector<T> >

◆ EigFunctionType

template<typename 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.

◆ grad_function_type [1/2]

type for gradient

◆ grad_function_type [2/2]

type for gradient

◆ GradWienerType

type for gradient of Wiener process

◆ lapl_function_type

type for laplacian

◆ PiecewiseSpatialFunctionType

template<typename T >
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.

◆ PiecewiseTemporalSpatialFunctionType

template<typename T >
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.

◆ ReactionType

Type for reaction term mu^2.

◆ solution_function_type [1/2]

type for solution

◆ solution_function_type [2/2]

type for solution

◆ source_function_type [1/2]

using HArDCore2D::HMM_StefanPME_Transient::source_function_type = std::function<double(const double&, const VectorRd&, const Cell*)>

type for source

◆ source_function_type [2/2]

type for source (at a given fixed time)

◆ SpatialFunctionType

template<typename T >
using HArDCore2D::SpatialFunctionType = typedef std::function<T(const VectorRd&)>

Generic type for function that only depends on spatial coordinate.

◆ TemporalSpatialFunctionType

template<typename T >
using HArDCore2D::TemporalSpatialFunctionType = typedef std::function<T(const double&, const VectorRd&)>

Generic type for function that depends on space and time.

◆ tensor_function_type [1/2]

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)

◆ tensor_function_type [2/2]

type for diffusion tensor (does not depend on time)

◆ WienerType

type for the Wiener process at a fixed time (only variation in space)

Function Documentation

◆ apply_nonlinearity() [1/2]

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

◆ apply_nonlinearity() [2/2]

UVector HArDCore2D::HMM_StefanPME_Transient::apply_nonlinearity ( const UVector Y,
const std::string  type 
) const

◆ apply_nonlinearity_eW() [1/2]

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

◆ apply_nonlinearity_eW() [2/2]

UVector HArDCore2D::StochStefanPME::apply_nonlinearity_eW ( const std::string  type,
const UVector Y,
const Eigen::VectorXd &  I_W 
) const

◆ compute_source()

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.

Parameters
use_exact_sourceuse exact source or zero
tTime at which to compute the source term
WWiener process
grad_WGradient of Wiener process
Delta_exp_WLaplacian of exponential of Wiener process
mu_squaredmu^2 for the reaction term

◆ EnergyNorm() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::EnergyNorm ( const UVector Xh) const

Discrete energy norm (associated to the diffusion operator)

◆ EnergyNorm() [2/2]

double HArDCore2D::StochStefanPME::EnergyNorm ( const UVector Xh) const

Discrete energy norm (associated to the diffusion operator)

◆ get_assembly_time() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::get_assembly_time ( ) const
inline

cpu time to assemble the scheme

◆ get_assembly_time() [2/2]

double HArDCore2D::StochStefanPME::get_assembly_time ( ) const
inline

cpu time to assemble the scheme

◆ get_itime() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::get_itime ( size_t  idx) const
inline

various intermediate assembly times

◆ get_itime() [2/2]

double HArDCore2D::StochStefanPME::get_itime ( size_t  idx) const
inline

various intermediate assembly times

◆ get_MassT() [1/2]

Eigen::MatrixXd HArDCore2D::HMM_StefanPME_Transient::get_MassT ( size_t  iT) const
inline

Mass matrix in cell iT.

◆ get_MassT() [2/2]

Eigen::MatrixXd HArDCore2D::StochStefanPME::get_MassT ( size_t  iT) const
inline

Mass matrix in cell iT.

◆ get_nb_newton()

size_t HArDCore2D::HMM_StefanPME_Transient::get_nb_newton ( ) const
inline

number of Newton iterations

◆ get_solution_time()

double HArDCore2D::StochStefanPME::get_solution_time ( ) const
inline

cpu time to solve the linear systems

◆ get_solving_error() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::get_solving_error ( ) const
inline

residual after solving the scheme

◆ get_solving_error() [2/2]

double HArDCore2D::StochStefanPME::get_solving_error ( ) const
inline

residual after solving the scheme

◆ get_solving_time()

double HArDCore2D::HMM_StefanPME_Transient::get_solving_time ( ) const
inline

cpu time to solve the scheme

◆ HMM_StefanPME_Transient()

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.

Parameters
hmminstance of the hybridcore class (with basis functions, etc.)
kappadiffusion tensor
sourcesource term
BCtype of boundary conditions
exact_solutionexact solution
grad_exact_solutiongradient of the exact solution
zetafunction describing the nonlinearity
weightproportion of weight (between 0 and 1) of mass-lumping on the edges
solver_typetype of solver to use for the global system (bicgstab at the moment)
outputoptional argument for output of messages

◆ Integral()

double HArDCore2D::StochStefanPME::Integral ( const UVector Xh) const

Integral of a function given by a vector.

◆ iterate() [1/2]

UVector HArDCore2D::HMM_StefanPME_Transient::iterate ( const double  tps,
const double  dt,
const UVector Xn 
)

Execute one time iteration.

Parameters
tpscurrent time (for computing source)
dttime step
XnUnkown at previous time step

◆ iterate() [2/2]

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.

Parameters
tpscurrent time (for computing source)
dttime step
XnUnkown at previous time step
Ih_WInterpolate of Wiener process at time tps
WWiener process at time tps, to fix the BCs
mu_squaredReaction term
sourcesource term at time tps

◆ L2_MassLumped() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::L2_MassLumped ( const UVector Xh) const

Mass-lumped L2 norm of a function given by a vector.

◆ L2_MassLumped() [2/2]

double HArDCore2D::StochStefanPME::L2_MassLumped ( const UVector Xh) const

Mass-lumped L2 norm of a function given by a vector.

◆ Lp_MassLumped() [1/2]

double HArDCore2D::HMM_StefanPME_Transient::Lp_MassLumped ( const UVector Xh,
double  p 
) const

Mass-lumped Lp norm of a function given by a vector.

◆ Lp_MassLumped() [2/2]

double HArDCore2D::StochStefanPME::Lp_MassLumped ( const UVector Xh,
double  p 
) const

Mass-lumped Lp norm of a function given by a vector.

◆ SetWienerProcess()

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.

Parameters
caseWType of Wiener process (can be deterministic...)
tTime at which to compute the Wiener process
WTo store the Wiener process
I_WInterpolate of the calculated Wiener process
grad_WTo store its gradient
Delta_exp_WTo store the Laplacian of its exponential
mu_squaredTo store mu^2
verboseDetermine if text should be output

◆ SimulateWienerProcess()

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.

Parameters
dtTime increment to update the process
WtimeCurrent values of the time component W_i of W, will be updated
WComputed complete Wiener process
Ih_eigsInterpolates of the eigenfunctions e_i to compute the interpolate of the Wiener process
I_WTo store the interpolate of the Wiener process
grad_WTo store the gradient of the Wiener process
Delta_exp_WTo store the Laplacian of the exponential of the Wiener process
mu_squaredReaction term mu^2

◆ StDev()

static double HArDCore2D::StDev ( const Eigen::VectorXd &  v)
inlinestatic

◆ StochStefanPME()

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.

Parameters
hmminstance of the hybridcore class (with basis functions, etc.)
kappadiffusion tensor
zetafunction describing the nonlinearity
BCtype of boundary conditions
exact_solutionexact solution
grad_exact_solutiongradient of the exact solution
time_der_exact_solutiontime derivative of exact solution
minus_Lapl_exact_solutionminus Laplacian of the exact solution
weightproportion of weight (between 0 and 1) of mass-lumping on the edges
solver_typetype of solver to use for the global system
outputoptional argument for output of messages

Variable Documentation

◆ Delta_mui_ei

EigFunctionType<double> HArDCore2D::Delta_mui_ei
static
Initial value:
= [](const size_t &k, const size_t &l, const VectorRd & p)->double
{
return - mui(k,l) * (std::pow(k, 2) + std::pow(l, 2)) * std::pow(PI,2) * sin(k * PI * p.x()) * sin(l * PI * p.y());
}
Eigen::Vector2d VectorRd
Definition: basis.hpp:55
static std::function< double(const size_t &, const size_t &)> mui
Definition: HMM_StochTrans_StefanPME.hpp:79
static const double PI
Type for vector of vector.
Definition: HMM_StochTrans_StefanPME.hpp:77

◆ grad_mui_ei

EigFunctionType<VectorRd> HArDCore2D::grad_mui_ei
static
Initial value:
= [](const size_t &k, const size_t &l, const VectorRd & p)->VectorRd
{
return mui(k,l) * VectorRd(
k * PI * cos(k * PI * p.x()) * sin(l * PI * p.y()),
l * PI * sin(k * PI * p.x()) * cos(l * PI * p.y())
);
}
Eigen::Matrix< double, DIMENSION, 1 > VectorRd
Definition: Polytope2D.hpp:19

◆ mui

std::function<double(const size_t &, const size_t &)> HArDCore2D::mui
static
Initial value:
= [](const size_t &k, const size_t &l)->double
{
return 1./(std::pow(k, 2) + std::pow(l, 2));
}

◆ mui_ei

EigFunctionType<double> HArDCore2D::mui_ei
static
Initial value:
= [](const size_t &k, const size_t &l, const VectorRd & p)->double
{
return mui(k,l) * sin(k * PI * p.x()) * sin(l * PI * p.y());
}

◆ PI

const double HArDCore2D::PI = boost::math::constants::pi<double>()
static

Type for vector of vector.

◆ zero_scalar_function

SpatialFunctionType<double> HArDCore2D::zero_scalar_function = [](const VectorRd & x)->double { return 0.;}
static

Zero spatial function.