DUBeat
1.0.1
High-order discontinuous Galerkin methods and applications to cardiac electrophysiology
|
Class representing the resolution of problems using discontinuous Galerkin methods. More...
#include <model_DG.hpp>
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 ¶ms) const override |
Declare main parameters. More... | |
virtual void | parse_parameters (lifex::ParamHandler ¶ms) 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... | |
Class representing the resolution of problems using discontinuous Galerkin methods.
Definition at line 62 of file model_DG.hpp.
Constructor.
Definition at line 66 of file model_DG.hpp.
Default copy constructor.
Default const copy constructor.
Default move constructor.
|
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 >.
|
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.
|
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.
|
protected |
Conversion to DUB coefficients, const version.
Definition at line 699 of file model_DG.hpp.
|
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.
|
protected |
Conversion of a discretized solution from FEM coefficients to Dubiner coefficients.
Definition at line 681 of file model_DG.hpp.
|
protected |
Definition at line 719 of file model_DG.hpp.
|
protectedvirtual |
Definition at line 709 of file model_DG.hpp.
|
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.
|
protected |
Conversion to FEM coefficients, const version.
Definition at line 634 of file model_DG.hpp.
|
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.
|
protected |
Conversion of a discretized solution from Dubiner coefficients to FEM coefficients.
Definition at line 616 of file model_DG.hpp.
|
protected |
Definition at line 654 of file model_DG.hpp.
|
protectedvirtual |
Definition at line 644 of file model_DG.hpp.
|
protected |
Load the mesh from the default path.
Definition at line 504 of file model_DG.hpp.
|
protected |
Load the mesh from a user-defined path.
Definition at line 522 of file model_DG.hpp.
|
overridevirtual |
Declare main parameters.
Reimplemented in ModelDG_t< basis >, and DUBeat::models::MonodomainFHNDG< basis >.
Definition at line 262 of file model_DG.hpp.
|
protected |
Conversion of an analytical solution from FEM to basis coefficients.
|
protected |
Conversion of an analytical solution from FEM to basis coefficients.
Specialization for FEM basis.
Definition at line 737 of file model_DG.hpp.
|
protected |
Conversion of an analytical solution from FEM to basis coefficients.
Specialization for Dubiner basis.
Definition at line 748 of file model_DG.hpp.
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.
|
protected |
To inizialize the solutions using the deal.II reinit.
Definition at line 492 of file model_DG.hpp.
|
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.
|
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.
|
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.
|
overridevirtual |
Run the simulation.
Reimplemented in ModelDG_t< basis >, and DUBeat::models::MonodomainFHNDG< basis >.
Definition at line 351 of file model_DG.hpp.
|
protectedvirtual |
Setup of the problem before the resolution.
Definition at line 382 of file model_DG.hpp.
|
protected |
System solving.
Definition at line 538 of file model_DG.hpp.
|
protected |
Matrix assembler.
Definition at line 233 of file model_DG.hpp.
|
protected |
DoFHandler (internal use for useful already implemented methods).
Definition at line 228 of file model_DG.hpp.
|
protected |
Number of degrees of freedom per cell.
Definition at line 226 of file model_DG.hpp.
|
protected |
Member used for conversions between analytical, nodal and modal representations of the solutions.
Definition at line 231 of file model_DG.hpp.
|
protected |
Known forcing term.
Definition at line 255 of file model_DG.hpp.
|
protected |
Neumann boundary conditions.
Definition at line 257 of file model_DG.hpp.
|
protected |
Pointer to exact gradient solution Function.
Definition at line 253 of file model_DG.hpp.
|
protected |
Linear solver handler.
Definition at line 235 of file model_DG.hpp.
|
protected |
Distributed matrix of the linear system.
Definition at line 239 of file model_DG.hpp.
|
protected |
Name of the class/problem.
Definition at line 212 of file model_DG.hpp.
|
protected |
Preconditioner handler.
Definition at line 237 of file model_DG.hpp.
|
protected |
Polynomials degree.
Definition at line 214 of file model_DG.hpp.
|
protected |
Mesh refinement level (>=1).
Definition at line 216 of file model_DG.hpp.
|
protected |
DG Penalty coefficient.
Definition at line 218 of file model_DG.hpp.
|
protected |
DG stabilty coefficient.
Definition at line 220 of file model_DG.hpp.
|
protected |
Distributed right hand side vector of the linear system.
Definition at line 241 of file model_DG.hpp.
|
protected |
Scaling factor.
Definition at line 222 of file model_DG.hpp.
|
protected |
Distributed solution vector, with ghost entries.
Definition at line 245 of file model_DG.hpp.
|
protected |
Distributed exact solution vector, without ghost entries.
Definition at line 249 of file model_DG.hpp.
|
protected |
Distributed exact solution vector, without ghost entries.
Definition at line 247 of file model_DG.hpp.
|
protected |
Distributed solution vector, without ghost entries.
Definition at line 243 of file model_DG.hpp.
|
protected |
Triangulation (internal use for useful already implemented methods).
Definition at line 224 of file model_DG.hpp.
|
protected |
Pointer to exact solution function.
Definition at line 251 of file model_DG.hpp.