DUBeat  1.0.1
High-order discontinuous Galerkin methods and applications to cardiac electrophysiology
ModelDG< basis > Class Template Referenceabstract

Class representing the resolution of problems using discontinuous Galerkin methods. More...

#include <model_DG.hpp>

Inheritance diagram for ModelDG< basis >:
Collaboration diagram for ModelDG< basis >:

Public Member Functions

 ModelDG (std::string model_name)
 Constructor. More...
 
 ModelDG (ModelDG< basis > &ModelDG)=default
 Default copy constructor. More...
 
 ModelDG (const ModelDG< basis > &ModelDG)=default
 Default const copy constructor. More...
 
 ModelDG (ModelDG< basis > &&ModelDG)=default
 Default move constructor. More...
 
virtual void declare_parameters (lifex::ParamHandler &params) const override
 Declare main parameters. More...
 
virtual void parse_parameters (lifex::ParamHandler &params) override
 Parse parameters from .prm file. More...
 
unsigned int get_dofs_per_cell () const
 Return the number of degrees of freedom per element. More...
 
virtual void run () override
 Run the simulation. More...
 
virtual ~ModelDG ()=default
 Destructor. More...
 

Protected Member Functions

virtual void setup_system ()
 Setup of the problem before the resolution. More...
 
void make_sparsity_pattern (const DoFHandlerDG< basis > &dof, dealii::DynamicSparsityPattern &sparsity, const dealii::AffineConstraints< double > &constraints=dealii::AffineConstraints< double >(), const bool keep_constrained_dofs=true, const dealii::types::subdomain_id subdomain_id=dealii::numbers::invalid_subdomain_id)
 Creation of the sparsity pattern to assign to the system matrix before assembling. More...
 
void initialize_solution (lifex::LinAlg::MPI::Vector &solution_owned, lifex::LinAlg::MPI::Vector &solution)
 To inizialize the solutions using the deal.II reinit. More...
 
virtual void assemble_system ()=0
 Assembly of the linear system, pure virtual. More...
 
void create_mesh ()
 Load the mesh from the default path. More...
 
void create_mesh (std::string mesh_path)
 Load the mesh from a user-defined path. More...
 
void solve_system ()
 System solving. More...
 
void compute_errors (const lifex::LinAlg::MPI::Vector &solution_owned, const lifex::LinAlg::MPI::Vector &solution_ex_owned, const std::shared_ptr< dealii::Function< lifex::dim >> &u_ex, const std::shared_ptr< dealii::Function< lifex::dim >> &grad_u_ex, const char *solution_name=(char *)"u") const
 To compute the \(L^\infty\) error, the \(L^2\) error, the \(H^1\) error and the \(DG\) error at the end of system solving, it exploits the DGComputeErrors<basis> class. More...
 
virtual void output_results (std::string output_name="solution") const
 Output of results. More...
 
virtual void conversion_to_fem (lifex::LinAlg::MPI::Vector &sol_owned)
 To convert a discretized solution from modal to nodal basis (does nothing if problem is already in nodal basis), in-place version. More...
 
virtual lifex::LinAlg::MPI::Vector conversion_to_fem (const lifex::LinAlg::MPI::Vector &sol_owned) const
 To convert a discretized solution from modal to nodal basis (does nothing if problem is already in nodal basis), const version. More...
 
virtual void conversion_to_fem (lifex::LinAlg::MPI::Vector &sol_owned, const std::string fem_file_path, const unsigned int degree_fem=1, const double scaling_factor=1) const
 
virtual void conversion_to_dub (lifex::LinAlg::MPI::Vector &sol_owned)
 To convert a discretized solution in Dubiner basis (only for problems using Dubiner basis). More...
 
virtual lifex::LinAlg::MPI::Vector conversion_to_dub (const lifex::LinAlg::MPI::Vector &sol_owned) const
 To convert a discretized solution in Dubiner basis (only for problems using Dubiner basis). More...
 
virtual void conversion_to_dub (lifex::LinAlg::MPI::Vector &sol_owned, const std::string fem_file_path, const unsigned int degree_fem=1, const double scaling_factor=1) const
 
void discretize_analytical_solution (const std::shared_ptr< dealii::Function< lifex::dim >> &u_analytical, lifex::LinAlg::MPI::Vector &sol_owned)
 Conversion of an analytical solution from FEM to basis coefficients. More...
 
void conversion_to_fem (lifex::LinAlg::MPI::Vector &sol_owned)
 Conversion of a discretized solution from Dubiner coefficients to FEM coefficients. More...
 
lifex::LinAlg::MPI::Vector conversion_to_fem (const lifex::LinAlg::MPI::Vector &sol_owned) const
 Conversion to FEM coefficients, const version. More...
 
void conversion_to_fem (lifex::LinAlg::MPI::Vector &sol_owned, const std::string fem_file_path, const unsigned int degree_fem, const double scaling_factor) const
 
void conversion_to_dub (lifex::LinAlg::MPI::Vector &sol_owned)
 Conversion of a discretized solution from FEM coefficients to Dubiner coefficients. More...
 
lifex::LinAlg::MPI::Vector conversion_to_dub (const lifex::LinAlg::MPI::Vector &sol_owned) const
 Conversion to DUB coefficients, const version. More...
 
void conversion_to_dub (lifex::LinAlg::MPI::Vector &sol_owned, const std::string fem_file_path, const unsigned int degree_fem, const double scaling_factor) const
 
void discretize_analytical_solution (const std::shared_ptr< dealii::Function< lifex::dim >> &u_analytical, lifex::LinAlg::MPI::Vector &sol_owned)
 Conversion of an analytical solution from FEM to basis coefficients. More...
 
void discretize_analytical_solution (const std::shared_ptr< dealii::Function< lifex::dim >> &u_analytical, lifex::LinAlg::MPI::Vector &sol_owned)
 Conversion of an analytical solution from FEM to basis coefficients. More...
 

Protected Attributes

const std::string model_name
 Name of the class/problem. More...
 
unsigned int prm_fe_degree
 Polynomials degree. More...
 
unsigned int prm_n_refinements
 Mesh refinement level (>=1). More...
 
double prm_penalty_coeff
 DG Penalty coefficient. More...
 
double prm_stability_coeff
 DG stabilty coefficient. More...
 
double scaling_factor = 1
 Scaling factor. More...
 
std::shared_ptr< lifex::utils::MeshHandler > triangulation
 Triangulation (internal use for useful already implemented methods). More...
 
unsigned int dofs_per_cell
 Number of degrees of freedom per cell. More...
 
DoFHandlerDG< basis > dof_handler
 DoFHandler (internal use for useful already implemented methods). More...
 
std::shared_ptr< DUBFEMHandler< basis > > dub_fem_values
 Member used for conversions between analytical, nodal and modal representations of the solutions. More...
 
std::unique_ptr< AssembleDG< basis > > assemble
 Matrix assembler. More...
 
lifex::utils::LinearSolverHandler< lifex::LinAlg::MPI::Vector > linear_solver
 Linear solver handler. More...
 
lifex::utils::PreconditionerHandler preconditioner
 Preconditioner handler. More...
 
lifex::LinAlg::MPI::SparseMatrix matrix
 Distributed matrix of the linear system. More...
 
lifex::LinAlg::MPI::Vector rhs
 Distributed right hand side vector of the linear system. More...
 
lifex::LinAlg::MPI::Vector solution_owned
 Distributed solution vector, without ghost entries. More...
 
lifex::LinAlg::MPI::Vector solution
 Distributed solution vector, with ghost entries. More...
 
lifex::LinAlg::MPI::Vector solution_ex_owned
 Distributed exact solution vector, without ghost entries. More...
 
lifex::LinAlg::MPI::Vector solution_ex
 Distributed exact solution vector, without ghost entries. More...
 
std::shared_ptr< lifex::utils::FunctionDirichlet > u_ex
 Pointer to exact solution function. More...
 
std::shared_ptr< dealii::Function< lifex::dim > > grad_u_ex
 Pointer to exact gradient solution Function. More...
 
std::shared_ptr< dealii::Function< lifex::dim > > f_ex
 Known forcing term. More...
 
std::shared_ptr< dealii::Function< lifex::dim > > g_n
 Neumann boundary conditions. More...
 

Detailed Description

template<class basis>
class ModelDG< basis >

Class representing the resolution of problems using discontinuous Galerkin methods.

Definition at line 62 of file model_DG.hpp.

Constructor & Destructor Documentation

◆ ModelDG() [1/4]

template<class basis >
ModelDG< basis >::ModelDG ( std::string  model_name)
inline

Constructor.

Definition at line 66 of file model_DG.hpp.

◆ ModelDG() [2/4]

template<class basis >
ModelDG< basis >::ModelDG ( ModelDG< basis > &  ModelDG)
default

Default copy constructor.

◆ ModelDG() [3/4]

template<class basis >
ModelDG< basis >::ModelDG ( const ModelDG< basis > &  ModelDG)
default

Default const copy constructor.

◆ ModelDG() [4/4]

template<class basis >
ModelDG< basis >::ModelDG ( ModelDG< basis > &&  ModelDG)
default

Default move constructor.

◆ ~ModelDG()

template<class basis >
virtual ModelDG< basis >::~ModelDG ( )
virtualdefault

Destructor.

Member Function Documentation

◆ assemble_system()

template<class basis >
virtual void ModelDG< basis >::assemble_system ( )
protectedpure virtual

Assembly of the linear system, pure virtual.

To specialize in relation to the model to solve.

Implemented in DUBeat::models::MonodomainFHNDG< basis >, DUBeat::models::LaplaceDG< basis >, and DUBeat::models::HeatDG< basis >.

◆ compute_errors()

template<class basis >
void ModelDG< basis >::compute_errors ( const lifex::LinAlg::MPI::Vector &  solution_owned,
const lifex::LinAlg::MPI::Vector &  solution_ex_owned,
const std::shared_ptr< dealii::Function< lifex::dim >> &  u_ex,
const std::shared_ptr< dealii::Function< lifex::dim >> &  grad_u_ex,
const char *  solution_name = (char *)"u" 
) const
protected

To compute the \(L^\infty\) error, the \(L^2\) error, the \(H^1\) error and the \(DG\) error at the end of system solving, it exploits the DGComputeErrors<basis> class.

Definition at line 547 of file model_DG.hpp.

◆ conversion_to_dub() [1/6]

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::conversion_to_dub ( const lifex::LinAlg::MPI::Vector &  sol_owned) const
protectedvirtual

To convert a discretized solution in Dubiner basis (only for problems using Dubiner basis).

Conversion to DUB coefficients, const version.

Const version.

Definition at line 690 of file model_DG.hpp.

◆ conversion_to_dub() [2/6]

lifex::LinAlg::MPI::Vector ModelDG< DUBValues< lifex::dim > >::conversion_to_dub ( const lifex::LinAlg::MPI::Vector &  sol_owned) const
protected

Conversion to DUB coefficients, const version.

Definition at line 699 of file model_DG.hpp.

◆ conversion_to_dub() [3/6]

template<class basis >
void ModelDG< basis >::conversion_to_dub ( lifex::LinAlg::MPI::Vector &  sol_owned)
protectedvirtual

To convert a discretized solution in Dubiner basis (only for problems using Dubiner basis).

Conversion of a discretized solution from FEM coefficients to Dubiner coefficients.

In-place version.

Useless if we are not using Dubiner basis functions.

Definition at line 672 of file model_DG.hpp.

◆ conversion_to_dub() [4/6]

void ModelDG< DUBValues< lifex::dim > >::conversion_to_dub ( lifex::LinAlg::MPI::Vector &  sol_owned)
protected

Conversion of a discretized solution from FEM coefficients to Dubiner coefficients.

Definition at line 681 of file model_DG.hpp.

◆ conversion_to_dub() [5/6]

void ModelDG< DUBValues< lifex::dim > >::conversion_to_dub ( lifex::LinAlg::MPI::Vector &  sol_owned,
const std::string  fem_file_path,
const unsigned int  degree_fem,
const double  scaling_factor 
) const
protected

Definition at line 719 of file model_DG.hpp.

◆ conversion_to_dub() [6/6]

template<class basis >
void ModelDG< basis >::conversion_to_dub ( lifex::LinAlg::MPI::Vector &  sol_owned,
const std::string  fem_file_path,
const unsigned int  degree_fem = 1,
const double  scaling_factor = 1 
) const
protectedvirtual

Definition at line 709 of file model_DG.hpp.

◆ conversion_to_fem() [1/6]

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::conversion_to_fem ( const lifex::LinAlg::MPI::Vector &  sol_owned) const
protectedvirtual

To convert a discretized solution from modal to nodal basis (does nothing if problem is already in nodal basis), const version.

Conversion to FEM coefficients, const version.

Definition at line 625 of file model_DG.hpp.

◆ conversion_to_fem() [2/6]

lifex::LinAlg::MPI::Vector ModelDG< DUBValues< lifex::dim > >::conversion_to_fem ( const lifex::LinAlg::MPI::Vector &  sol_owned) const
protected

Conversion to FEM coefficients, const version.

Definition at line 634 of file model_DG.hpp.

◆ conversion_to_fem() [3/6]

template<class basis >
void ModelDG< basis >::conversion_to_fem ( lifex::LinAlg::MPI::Vector &  sol_owned)
protectedvirtual

To convert a discretized solution from modal to nodal basis (does nothing if problem is already in nodal basis), in-place version.

Conversion of a discretized solution from Dubiner coefficients to FEM coefficients.

It exploits dubiner_to_fem from DUB_FEM_handler.hpp so the same considerations as in output_results hold.

Useless if we are not using Dubiner basis functions.

Definition at line 607 of file model_DG.hpp.

◆ conversion_to_fem() [4/6]

void ModelDG< DUBValues< lifex::dim > >::conversion_to_fem ( lifex::LinAlg::MPI::Vector &  sol_owned)
protected

Conversion of a discretized solution from Dubiner coefficients to FEM coefficients.

Definition at line 616 of file model_DG.hpp.

◆ conversion_to_fem() [5/6]

void ModelDG< DUBValues< lifex::dim > >::conversion_to_fem ( lifex::LinAlg::MPI::Vector &  sol_owned,
const std::string  fem_file_path,
const unsigned int  degree_fem,
const double  scaling_factor 
) const
protected

Definition at line 654 of file model_DG.hpp.

◆ conversion_to_fem() [6/6]

template<class basis >
void ModelDG< basis >::conversion_to_fem ( lifex::LinAlg::MPI::Vector &  sol_owned,
const std::string  fem_file_path,
const unsigned int  degree_fem = 1,
const double  scaling_factor = 1 
) const
protectedvirtual

Definition at line 644 of file model_DG.hpp.

◆ create_mesh() [1/2]

template<class basis >
void ModelDG< basis >::create_mesh
protected

Load the mesh from the default path.

Definition at line 504 of file model_DG.hpp.

◆ create_mesh() [2/2]

template<class basis >
void ModelDG< basis >::create_mesh ( std::string  mesh_path)
protected

Load the mesh from a user-defined path.

Definition at line 522 of file model_DG.hpp.

◆ declare_parameters()

template<class basis >
void ModelDG< basis >::declare_parameters ( lifex::ParamHandler &  params) const
overridevirtual

Declare main parameters.

Reimplemented in ModelDG_t< basis >, and DUBeat::models::MonodomainFHNDG< basis >.

Definition at line 262 of file model_DG.hpp.

◆ discretize_analytical_solution() [1/3]

template<class basis >
void ModelDG< basis >::discretize_analytical_solution ( const std::shared_ptr< dealii::Function< lifex::dim >> &  u_analytical,
lifex::LinAlg::MPI::Vector &  sol_owned 
)
protected

Conversion of an analytical solution from FEM to basis coefficients.

◆ discretize_analytical_solution() [2/3]

void ModelDG< dealii::FE_SimplexDGP< lifex::dim > >::discretize_analytical_solution ( const std::shared_ptr< dealii::Function< lifex::dim >> &  u_analytical,
lifex::LinAlg::MPI::Vector &  sol_owned 
)
protected

Conversion of an analytical solution from FEM to basis coefficients.

Specialization for FEM basis.

Definition at line 737 of file model_DG.hpp.

◆ discretize_analytical_solution() [3/3]

void ModelDG< DUBValues< lifex::dim > >::discretize_analytical_solution ( const std::shared_ptr< dealii::Function< lifex::dim >> &  u_analytical,
lifex::LinAlg::MPI::Vector &  sol_owned 
)
protected

Conversion of an analytical solution from FEM to basis coefficients.

Specialization for Dubiner basis.

Definition at line 748 of file model_DG.hpp.

◆ get_dofs_per_cell()

template<class basis >
unsigned int ModelDG< basis >::get_dofs_per_cell

Return the number of degrees of freedom per element.

Definition at line 330 of file model_DG.hpp.

◆ initialize_solution()

template<class basis >
void ModelDG< basis >::initialize_solution ( lifex::LinAlg::MPI::Vector &  solution_owned,
lifex::LinAlg::MPI::Vector &  solution 
)
protected

To inizialize the solutions using the deal.II reinit.

Definition at line 492 of file model_DG.hpp.

◆ make_sparsity_pattern()

template<class basis >
void ModelDG< basis >::make_sparsity_pattern ( const DoFHandlerDG< basis > &  dof,
dealii::DynamicSparsityPattern &  sparsity,
const dealii::AffineConstraints< double > &  constraints = dealii::AffineConstraints<double>(),
const bool  keep_constrained_dofs = true,
const dealii::types::subdomain_id  subdomain_id = dealii::numbers::invalid_subdomain_id 
)
protected

Creation of the sparsity pattern to assign to the system matrix before assembling.

It has been readapted from the deal.II DoFTools::make_sparsity_pattern() method.

Definition at line 440 of file model_DG.hpp.

◆ output_results()

template<class basis >
void ModelDG< basis >::output_results ( std::string  output_name = "solution") const
protectedvirtual

Output of results.

Note that, since it exploits dubiner_to_fem from DUB_FEM_handler.hpp, if the polynomial order chosen to solve the model is \(>2\), the FE space considered for the output will be of order at most 2. This leads to an output vector which correspondent visualization might be less refined. See DUB_FEM_handler.hpp for more details.

Definition at line 581 of file model_DG.hpp.

◆ parse_parameters()

template<class basis >
void ModelDG< basis >::parse_parameters ( lifex::ParamHandler &  params)
overridevirtual

Parse parameters from .prm file.

Reimplemented in ModelDG_t< basis >, and DUBeat::models::MonodomainFHNDG< basis >.

Definition at line 301 of file model_DG.hpp.

◆ run()

template<class basis >
void ModelDG< basis >::run
overridevirtual

Run the simulation.

Reimplemented in ModelDG_t< basis >, and DUBeat::models::MonodomainFHNDG< basis >.

Definition at line 351 of file model_DG.hpp.

◆ setup_system()

template<class basis >
void ModelDG< basis >::setup_system
protectedvirtual

Setup of the problem before the resolution.

Definition at line 382 of file model_DG.hpp.

◆ solve_system()

template<class basis >
void ModelDG< basis >::solve_system
protected

System solving.

Definition at line 538 of file model_DG.hpp.

Member Data Documentation

◆ assemble

template<class basis >
std::unique_ptr<AssembleDG<basis> > ModelDG< basis >::assemble
protected

Matrix assembler.

Definition at line 233 of file model_DG.hpp.

◆ dof_handler

template<class basis >
DoFHandlerDG<basis> ModelDG< basis >::dof_handler
protected

DoFHandler (internal use for useful already implemented methods).

Definition at line 228 of file model_DG.hpp.

◆ dofs_per_cell

template<class basis >
unsigned int ModelDG< basis >::dofs_per_cell
protected

Number of degrees of freedom per cell.

Definition at line 226 of file model_DG.hpp.

◆ dub_fem_values

template<class basis >
std::shared_ptr<DUBFEMHandler<basis> > ModelDG< basis >::dub_fem_values
protected

Member used for conversions between analytical, nodal and modal representations of the solutions.

Definition at line 231 of file model_DG.hpp.

◆ f_ex

template<class basis >
std::shared_ptr<dealii::Function<lifex::dim> > ModelDG< basis >::f_ex
protected

Known forcing term.

Definition at line 255 of file model_DG.hpp.

◆ g_n

template<class basis >
std::shared_ptr<dealii::Function<lifex::dim> > ModelDG< basis >::g_n
protected

Neumann boundary conditions.

Definition at line 257 of file model_DG.hpp.

◆ grad_u_ex

template<class basis >
std::shared_ptr<dealii::Function<lifex::dim> > ModelDG< basis >::grad_u_ex
protected

Pointer to exact gradient solution Function.

Definition at line 253 of file model_DG.hpp.

◆ linear_solver

template<class basis >
lifex::utils::LinearSolverHandler<lifex::LinAlg::MPI::Vector> ModelDG< basis >::linear_solver
protected

Linear solver handler.

Definition at line 235 of file model_DG.hpp.

◆ matrix

template<class basis >
lifex::LinAlg::MPI::SparseMatrix ModelDG< basis >::matrix
protected

Distributed matrix of the linear system.

Definition at line 239 of file model_DG.hpp.

◆ model_name

template<class basis >
const std::string ModelDG< basis >::model_name
protected

Name of the class/problem.

Definition at line 212 of file model_DG.hpp.

◆ preconditioner

template<class basis >
lifex::utils::PreconditionerHandler ModelDG< basis >::preconditioner
protected

Preconditioner handler.

Definition at line 237 of file model_DG.hpp.

◆ prm_fe_degree

template<class basis >
unsigned int ModelDG< basis >::prm_fe_degree
protected

Polynomials degree.

Definition at line 214 of file model_DG.hpp.

◆ prm_n_refinements

template<class basis >
unsigned int ModelDG< basis >::prm_n_refinements
protected

Mesh refinement level (>=1).

Definition at line 216 of file model_DG.hpp.

◆ prm_penalty_coeff

template<class basis >
double ModelDG< basis >::prm_penalty_coeff
protected

DG Penalty coefficient.

Definition at line 218 of file model_DG.hpp.

◆ prm_stability_coeff

template<class basis >
double ModelDG< basis >::prm_stability_coeff
protected

DG stabilty coefficient.

Definition at line 220 of file model_DG.hpp.

◆ rhs

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::rhs
protected

Distributed right hand side vector of the linear system.

Definition at line 241 of file model_DG.hpp.

◆ scaling_factor

template<class basis >
double ModelDG< basis >::scaling_factor = 1
protected

Scaling factor.

Definition at line 222 of file model_DG.hpp.

◆ solution

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::solution
protected

Distributed solution vector, with ghost entries.

Definition at line 245 of file model_DG.hpp.

◆ solution_ex

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::solution_ex
protected

Distributed exact solution vector, without ghost entries.

Definition at line 249 of file model_DG.hpp.

◆ solution_ex_owned

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::solution_ex_owned
protected

Distributed exact solution vector, without ghost entries.

Definition at line 247 of file model_DG.hpp.

◆ solution_owned

template<class basis >
lifex::LinAlg::MPI::Vector ModelDG< basis >::solution_owned
protected

Distributed solution vector, without ghost entries.

Definition at line 243 of file model_DG.hpp.

◆ triangulation

template<class basis >
std::shared_ptr<lifex::utils::MeshHandler> ModelDG< basis >::triangulation
protected

Triangulation (internal use for useful already implemented methods).

Definition at line 224 of file model_DG.hpp.

◆ u_ex

template<class basis >
std::shared_ptr<lifex::utils::FunctionDirichlet> ModelDG< basis >::u_ex
protected

Pointer to exact solution function.

Definition at line 251 of file model_DG.hpp.


The documentation for this class was generated from the following file: