CEPS  24.01
Cardiac ElectroPhysiology Simulator
AbstractDiscretization Class Referenceabstract

Detailed Description

Abstract Class for all numerical method (FE, FD, FV etc)

Definition at line 60 of file AbstractDiscretization.hpp.

#include <AbstractDiscretization.hpp>

Inheritance diagram for AbstractDiscretization:
[legend]

Public Member Functions

 AbstractDiscretization ()=delete
 Deleted default constructor. More...
 
 AbstractDiscretization (const AbstractDiscretization &that)=delete
 Deleted Copy constructor (objects are too big, and not really useful) More...
 
AbstractDiscretizationoperator= (const AbstractDiscretization &)=delete
 Deleted assignment operator. More...
 
 ~AbstractDiscretization () override
 Destructor. More...
 
AbstractPdeProblemgetProblem () const
 Return link to pde problem. More...
 
DistributedFactorygetDofsFactory () const
 Get the stored DistributedFactory. More...
 
DistributedInfos< DegreeOfFreedom * > * getDistributedDofs () const
 Get the stored Degree Of Freedom repartition. More...
 
const DegreeOfFreedomTreegetDegreeOfFreedomTree () const
 The sorted degrees of freedom. More...
 
const CepsMap< CepsUnknownIndex, CepsVector< DegreeOfFreedom * > > & getDegreesOfFreedomForUnknown ()
 Dofs sorted by unknown. More...
 
const CepsVector< DegreeOfFreedom * > & getDegreesOfFreedomForUnknown (CepsUnknownIndex uId)
 Dofs for given unknown. More...
 
const CepsVector< DegreeOfFreedom * > & getDegreesOfFreedomForUnknown (Unknown *u)
 Dofs sorted by unknown. More...
 
GeometrygetGeometry () const
 Get the stored Geometry. More...
 
virtual void registerSpatialUnkown (Unknown *u)
 Changes the location of unknown u to be adapted to the discrectization for example CepsLocationFlag::Point for finite elements, Cell for finite volumes. More...
 
DVecPtr newDofVector () const
 Get a new vector from the factory. More...
 
DHVecPtr newDofHaloVector () const
 Get a new vector from the factory, with halo data. More...
 
DMatPtr newDofMatrix () const
 Get a new matrix from the factory. More...
 
DegreeOfFreedomgetDof (CepsGlobalIndex globalDofId)
 Get the dof with given globalID, nullptr if not found. More...
 
const DegreeOfFreedomgetDof (CepsGlobalIndex globalDofId) const
 Get the dof with given globalID, nullptr if not found. More...
 
CepsBool isDofSpatial (CepsGlobalIndex globalDofId) const
 Check if the Unknow is a spatial one (typically CepsLocationFlag != ZeroD) More...
 
virtual CepsIndex getDofSpatialId (const DegreeOfFreedom *dof) const =0
 Return the spatial id (CellId or NodeId) for this row number. More...
 
CepsVector< DegreeOfFreedom * > getOwnedDofs () const
 Get the degrees of freedom owned by process. More...
 
CepsVector< DegreeOfFreedom * > getHaloDofs () const
 Get the degrees of freedom owned by other processes. More...
 
virtual void extractValuesForUnknown (const CepsUnknownIndex &uId, DVecPtr vec, CepsVector< CepsIndex > &indices, CepsVector< CepsReal > &values, CepsBool sendToMaster=false, CepsBool returnGeomIndices=false)=0
 Slicing method that extracts the data corresponding to a specific unknown. More...
 
void evaluateFunctionOnDofs (ceps::Function< CepsReal(CepsStandardArgs)> *func, DHVecPtr v, CepsReal t=0.) const
 Fills vector v with values of function func at owned and halo dofs and time t. More...
 
void findClosestDof (const CepsReal3D &x, Unknown *u, CepsDofGlobalIndex &dofId, CepsProcId &owner)
 Set dofId to that if the closest to position x, also sets owner. More...
 
CepsReal getZeroDValue (Unknown *u, DHVecPtr sol)
 Get the value of a 0D unknown into a solution vector. More...
 
virtual CepsReal dotProduct (DHVecPtr u, DHVecPtr v, CepsBool boundary=false, const CepsSet< CepsAttribute > &attrs={}, const CepsVector< Unknown * > &unknowns={})=0
 Returns the dot product of two vectors of degrees of freedom. More...
 
CepsReal l1Norm (DHVecPtr v, CepsBool boundary=false, const CepsSet< CepsAttribute > &attrs={}, const CepsVector< Unknown * > &unknowns={})
 L1-norm of a given vector. More...
 
CepsReal l2Norm (DHVecPtr v, CepsBool boundary=false, const CepsSet< CepsAttribute > &attrs={}, const CepsVector< Unknown * > &unknowns={})
 L2-norm of a given vector. Computation of dotproduct depends on the discretization. More...
 
virtual CepsReal h1Norm (DHVecPtr v, CepsBool boundary=false, const CepsSet< CepsAttribute > &attrs={}, const CepsVector< Unknown * > &unknowns={})=0
 H1-norm of a given vector. Computation depends on the discretization. More...
 
- Public Member Functions inherited from CepsObject
 CepsObject ()=default
 default constructor More...
 
 CepsObject (const CepsObject &)=default
 Copy constructor. More...
 
virtual ~CepsObject ()=default
 Destructor. More...
 
CepsObjectoperator= (const CepsObject &)=default
 Assignment operator. More...
 
CepsObjecttoBaseObject ()
 Returns a pointer to CepsObject class. More...
 
const CepsObjecttoBaseObject () const
 Returns a pointer to CepsObject class, const version. More...
 
ProfilergetProfiler () const
 Access to profiler. More...
 

Protected Member Functions

 AbstractDiscretization (AbstractPdeProblem *problem)
 Constructor to build an object from PDE problem. More...
 
void buildDofs ()
 Builds the degrees of freedom data structures from the PDE data. More...
 
virtual void buildDofsTree ()=0
 Builds the degrees of freedom tree from the PDE data. More...
 
void buildDofsDistributedInfo ()
 Builds the degrees of freedom distribution from the dofs tree. More...
 
void createDofsForZeroDUnknowns ()
 Adds one dof for each zeroD unknown, must be called after all other dofs are distributed. More...
 
virtual void setSpatialUnknownsDofsInteractions ()=0
 Build dofs neighboring list for spatial unknowns. Different in FE (same element dofs) and VF, FD dofs on neighbor geom elements. More...
 
virtual void setZeroDUnknownsDofsInteractions ()
 Updates neighboring lists with interaction with 0D unknowns. More...
 
void setupDofsFactory ()
 Prepare the factory that will generate vectors of dofs. More...
 
void extractValuesForUnknownInternal (const CepsUnknownIndex &uId, DVecPtr vec, CepsVector< CepsIndex > &indices, CepsVector< CepsReal > &values, CepsBool sendToMaster, CepsBool returnGeomIndices, std::function< CepsBool(DegreeOfFreedom *)> &selector)
 Slicing method that extracts the data corresponding to a specific unknown. More...
 

Protected Attributes

AbstractPdeProblemm_problem
 [not owned] Pointer on the pde problem More...
 
Geometrym_geometry
 [not owned] Pointer on the geometry More...
 
DistributedInfos< DegreeOfFreedom * > * m_distributedDofs
 [owned] Pointer on the dofs structure More...
 
DegreeOfFreedomTree m_dofsTree
 [owned] Pointer on the distributed degree of freedom tree More...
 
DistributedFactorym_dofsFactory
 [owned] The factory that builds vectors and matrices More...
 
CepsBool m_dofsBuilt
 Tells if the class has setup all the dofs. More...
 
CepsMap< CepsUnknownIndex, CepsVector< DegreeOfFreedom * > > m_dofsForUnknown
 Reverse mapping of dofs, with unknown index as key. More...
 
CepsLocationFlag m_spatialUnknownsLocation
 Where unknowns are located. Must be set in child classes. More...
 
CepsVector< CepsMap< CepsDofGlobalIndex, CepsSet< CepsDofGlobalIndex > > > m_extraAdjacencyToSend
 A map to correct the missing adjacency that can occur in rare occasions in 3D. More...
 
CepsMap< CepsDofGlobalIndex, CepsSet< CepsDofGlobalIndex > > m_extraAdjacencyToRecv
 A map to correct the missing adjacency that can occur in rare occasions in 3D. More...
 

Additional Inherited Members

- Static Protected Attributes inherited from CepsObject
static Profiler m_profiler
 The same profiler for each big object. More...
 

Constructor & Destructor Documentation

◆ AbstractDiscretization() [1/3]

AbstractDiscretization::AbstractDiscretization ( )
delete

Deleted default constructor.

◆ AbstractDiscretization() [2/3]

AbstractDiscretization::AbstractDiscretization ( const AbstractDiscretization that)
delete

Deleted Copy constructor (objects are too big, and not really useful)

◆ ~AbstractDiscretization()

AbstractDiscretization::~AbstractDiscretization ( )
override

Destructor.

Definition at line 52 of file AbstractDiscretization.cpp.

◆ AbstractDiscretization() [3/3]

AbstractDiscretization::AbstractDiscretization ( AbstractPdeProblem problem)
protected

Constructor to build an object from PDE problem.

Definition at line 34 of file AbstractDiscretization.cpp.

Member Function Documentation

◆ buildDofs()

void AbstractDiscretization::buildDofs ( )
protected

Builds the degrees of freedom data structures from the PDE data.

Definition at line 264 of file AbstractDiscretization.cpp.

◆ buildDofsDistributedInfo()

void AbstractDiscretization::buildDofsDistributedInfo ( )
protected

Builds the degrees of freedom distribution from the dofs tree.

Definition at line 288 of file AbstractDiscretization.cpp.

◆ buildDofsTree()

virtual void AbstractDiscretization::buildDofsTree ( )
protectedpure virtual

Builds the degrees of freedom tree from the PDE data.

Implemented in FiniteElements.

◆ createDofsForZeroDUnknowns()

void AbstractDiscretization::createDofsForZeroDUnknowns ( )
protected

Adds one dof for each zeroD unknown, must be called after all other dofs are distributed.

Definition at line 338 of file AbstractDiscretization.cpp.

◆ dotProduct()

virtual CepsReal AbstractDiscretization::dotProduct ( DHVecPtr  u,
DHVecPtr  v,
CepsBool  boundary = false,
const CepsSet< CepsAttribute > &  attrs = {},
const CepsVector< Unknown * > &  unknowns = {} 
)
pure virtual

Returns the dot product of two vectors of degrees of freedom.

Implemented in FiniteElements.

◆ evaluateFunctionOnDofs()

void AbstractDiscretization::evaluateFunctionOnDofs ( ceps::Function< CepsReal(CepsStandardArgs)> *  func,
DHVecPtr  v,
CepsReal  t = 0. 
) const

Fills vector v with values of function func at owned and halo dofs and time t.

Definition at line 169 of file AbstractDiscretization.cpp.

◆ extractValuesForUnknown()

virtual void AbstractDiscretization::extractValuesForUnknown ( const CepsUnknownIndex uId,
DVecPtr  vec,
CepsVector< CepsIndex > &  indices,
CepsVector< CepsReal > &  values,
CepsBool  sendToMaster = false,
CepsBool  returnGeomIndices = false 
)
pure virtual

Slicing method that extracts the data corresponding to a specific unknown.

Parameters
[in]uIdunknown identifier
[in]vecdistributed vector from which to extract
[out]indicesglobal row indices of the extracted dofs
[out]valuesextracted values
[in]sendToMasterdata will be shared with master CPU (e.g. for output)
[in]returnGeomIndicesfill indices with geomIndices instead of dofIndices (outputs) /!\ Only dofs defined on geom points will be parsed !

Implemented in FiniteElements.

◆ extractValuesForUnknownInternal()

void AbstractDiscretization::extractValuesForUnknownInternal ( const CepsUnknownIndex uId,
DVecPtr  vec,
CepsVector< CepsIndex > &  indices,
CepsVector< CepsReal > &  values,
CepsBool  sendToMaster,
CepsBool  returnGeomIndices,
std::function< CepsBool(DegreeOfFreedom *)> &  selector 
)
protected

Slicing method that extracts the data corresponding to a specific unknown.

Parameters
[in]uIdunknown identifier
[in]vecdistributed vector from which to extract
[out]indicesglobal row indices of the extracted dofs
[out]valuesextracted values
[in]sendToMasterdata will be shared with master CPU (e.g. for output)
[in]returnGeomIndicesfill indices with geomIndices instead of dofIndices (outputs)
[in]selectorfunction that tells how to keep dofs /!\ Only dofs defined on geom points will be parsed !

Definition at line 544 of file AbstractDiscretization.cpp.

◆ findClosestDof()

void AbstractDiscretization::findClosestDof ( const CepsReal3D x,
Unknown u,
CepsDofGlobalIndex dofId,
CepsProcId owner 
)

Set dofId to that if the closest to position x, also sets owner.

Definition at line 189 of file AbstractDiscretization.cpp.

◆ getDegreeOfFreedomTree()

const DegreeOfFreedomTree & AbstractDiscretization::getDegreeOfFreedomTree ( ) const

The sorted degrees of freedom.

Definition at line 83 of file AbstractDiscretization.cpp.

◆ getDegreesOfFreedomForUnknown() [1/3]

const CepsMap< CepsUnknownIndex, CepsVector< DegreeOfFreedom * > > & AbstractDiscretization::getDegreesOfFreedomForUnknown ( )

Dofs sorted by unknown.

Definition at line 89 of file AbstractDiscretization.cpp.

◆ getDegreesOfFreedomForUnknown() [2/3]

const CepsVector< DegreeOfFreedom * > & AbstractDiscretization::getDegreesOfFreedomForUnknown ( CepsUnknownIndex  uId)

Dofs for given unknown.

Definition at line 95 of file AbstractDiscretization.cpp.

◆ getDegreesOfFreedomForUnknown() [3/3]

const CepsVector< DegreeOfFreedom * > & AbstractDiscretization::getDegreesOfFreedomForUnknown ( Unknown u)

Dofs sorted by unknown.

Definition at line 101 of file AbstractDiscretization.cpp.

◆ getDistributedDofs()

DistributedInfos< DegreeOfFreedom * > * AbstractDiscretization::getDistributedDofs ( ) const

Get the stored Degree Of Freedom repartition.

Definition at line 77 of file AbstractDiscretization.cpp.

◆ getDof() [1/2]

DegreeOfFreedom * AbstractDiscretization::getDof ( CepsGlobalIndex  globalDofId)

Get the dof with given globalID, nullptr if not found.

Definition at line 139 of file AbstractDiscretization.cpp.

◆ getDof() [2/2]

const DegreeOfFreedom * AbstractDiscretization::getDof ( CepsGlobalIndex  globalDofId) const

Get the dof with given globalID, nullptr if not found.

Definition at line 148 of file AbstractDiscretization.cpp.

◆ getDofsFactory()

DistributedFactory * AbstractDiscretization::getDofsFactory ( ) const

Get the stored DistributedFactory.

Definition at line 71 of file AbstractDiscretization.cpp.

◆ getDofSpatialId()

virtual CepsIndex AbstractDiscretization::getDofSpatialId ( const DegreeOfFreedom dof) const
pure virtual

Return the spatial id (CellId or NodeId) for this row number.

Implemented in FiniteElements.

◆ getGeometry()

Geometry * AbstractDiscretization::getGeometry ( ) const

Get the stored Geometry.

Definition at line 107 of file AbstractDiscretization.cpp.

◆ getHaloDofs()

CepsVector< DegreeOfFreedom * > AbstractDiscretization::getHaloDofs ( ) const

Get the degrees of freedom owned by other processes.

Definition at line 163 of file AbstractDiscretization.cpp.

◆ getOwnedDofs()

CepsVector< DegreeOfFreedom * > AbstractDiscretization::getOwnedDofs ( ) const

Get the degrees of freedom owned by process.

Definition at line 157 of file AbstractDiscretization.cpp.

◆ getProblem()

AbstractPdeProblem * AbstractDiscretization::getProblem ( ) const

Return link to pde problem.

Definition at line 65 of file AbstractDiscretization.cpp.

◆ getZeroDValue()

CepsReal AbstractDiscretization::getZeroDValue ( Unknown u,
DHVecPtr  sol 
)

Get the value of a 0D unknown into a solution vector.

Definition at line 221 of file AbstractDiscretization.cpp.

◆ h1Norm()

virtual CepsReal AbstractDiscretization::h1Norm ( DHVecPtr  v,
CepsBool  boundary = false,
const CepsSet< CepsAttribute > &  attrs = {},
const CepsVector< Unknown * > &  unknowns = {} 
)
pure virtual

H1-norm of a given vector. Computation depends on the discretization.

Implemented in FiniteElements.

◆ isDofSpatial()

CepsBool AbstractDiscretization::isDofSpatial ( CepsGlobalIndex  globalDofId) const

Check if the Unknow is a spatial one (typically CepsLocationFlag != ZeroD)

◆ l1Norm()

CepsReal AbstractDiscretization::l1Norm ( DHVecPtr  v,
CepsBool  boundary = false,
const CepsSet< CepsAttribute > &  attrs = {},
const CepsVector< Unknown * > &  unknowns = {} 
)

L1-norm of a given vector.

Definition at line 236 of file AbstractDiscretization.cpp.

◆ l2Norm()

CepsReal AbstractDiscretization::l2Norm ( DHVecPtr  v,
CepsBool  boundary = false,
const CepsSet< CepsAttribute > &  attrs = {},
const CepsVector< Unknown * > &  unknowns = {} 
)

L2-norm of a given vector. Computation of dotproduct depends on the discretization.

Definition at line 251 of file AbstractDiscretization.cpp.

◆ newDofHaloVector()

DHVecPtr AbstractDiscretization::newDofHaloVector ( ) const

Get a new vector from the factory, with halo data.

Definition at line 127 of file AbstractDiscretization.cpp.

◆ newDofMatrix()

DMatPtr AbstractDiscretization::newDofMatrix ( ) const

Get a new matrix from the factory.

Definition at line 133 of file AbstractDiscretization.cpp.

◆ newDofVector()

DVecPtr AbstractDiscretization::newDofVector ( ) const

Get a new vector from the factory.

Definition at line 121 of file AbstractDiscretization.cpp.

◆ operator=()

AbstractDiscretization& AbstractDiscretization::operator= ( const AbstractDiscretization )
delete

Deleted assignment operator.

◆ registerSpatialUnkown()

void AbstractDiscretization::registerSpatialUnkown ( Unknown u)
virtual

Changes the location of unknown u to be adapted to the discrectization for example CepsLocationFlag::Point for finite elements, Cell for finite volumes.

Definition at line 113 of file AbstractDiscretization.cpp.

◆ setSpatialUnknownsDofsInteractions()

virtual void AbstractDiscretization::setSpatialUnknownsDofsInteractions ( )
protectedpure virtual

Build dofs neighboring list for spatial unknowns. Different in FE (same element dofs) and VF, FD dofs on neighbor geom elements.

Implemented in FiniteElements.

◆ setupDofsFactory()

void AbstractDiscretization::setupDofsFactory ( )
protected

Prepare the factory that will generate vectors of dofs.

Definition at line 412 of file AbstractDiscretization.cpp.

◆ setZeroDUnknownsDofsInteractions()

void AbstractDiscretization::setZeroDUnknownsDofsInteractions ( )
protectedvirtual

Updates neighboring lists with interaction with 0D unknowns.

Definition at line 379 of file AbstractDiscretization.cpp.

Field Documentation

◆ m_distributedDofs

DistributedInfos<DegreeOfFreedom*>* AbstractDiscretization::m_distributedDofs
protected

[owned] Pointer on the dofs structure

Definition at line 292 of file AbstractDiscretization.hpp.

◆ m_dofsBuilt

CepsBool AbstractDiscretization::m_dofsBuilt
protected

Tells if the class has setup all the dofs.

Definition at line 298 of file AbstractDiscretization.hpp.

◆ m_dofsFactory

DistributedFactory* AbstractDiscretization::m_dofsFactory
protected

[owned] The factory that builds vectors and matrices

Definition at line 296 of file AbstractDiscretization.hpp.

◆ m_dofsForUnknown

CepsMap<CepsUnknownIndex,CepsVector<DegreeOfFreedom*> > AbstractDiscretization::m_dofsForUnknown
protected

Reverse mapping of dofs, with unknown index as key.

Definition at line 300 of file AbstractDiscretization.hpp.

◆ m_dofsTree

DegreeOfFreedomTree AbstractDiscretization::m_dofsTree
protected

[owned] Pointer on the distributed degree of freedom tree

Definition at line 294 of file AbstractDiscretization.hpp.

◆ m_extraAdjacencyToRecv

CepsMap<CepsDofGlobalIndex,CepsSet<CepsDofGlobalIndex> > AbstractDiscretization::m_extraAdjacencyToRecv
protected

A map to correct the missing adjacency that can occur in rare occasions in 3D.

Definition at line 308 of file AbstractDiscretization.hpp.

◆ m_extraAdjacencyToSend

CepsVector<CepsMap<CepsDofGlobalIndex,CepsSet<CepsDofGlobalIndex> > > AbstractDiscretization::m_extraAdjacencyToSend
protected

A map to correct the missing adjacency that can occur in rare occasions in 3D.

Definition at line 306 of file AbstractDiscretization.hpp.

◆ m_geometry

Geometry* AbstractDiscretization::m_geometry
protected

[not owned] Pointer on the geometry

Definition at line 287 of file AbstractDiscretization.hpp.

◆ m_problem

AbstractPdeProblem* AbstractDiscretization::m_problem
protected

[not owned] Pointer on the pde problem

Definition at line 285 of file AbstractDiscretization.hpp.

◆ m_spatialUnknownsLocation

CepsLocationFlag AbstractDiscretization::m_spatialUnknownsLocation
protected

Where unknowns are located. Must be set in child classes.

Definition at line 303 of file AbstractDiscretization.hpp.


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