Flow123d  master-f44eb46
Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | Friends | List of all members
DOFHandlerMultiDim Class Reference

Provides the numbering of the finite element degrees of freedom on the computational mesh. More...

#include <dofhandler.hh>

Inheritance diagram for DOFHandlerMultiDim:
Inheritance graph
[legend]
Collaboration diagram for DOFHandlerMultiDim:
Collaboration graph
[legend]

Public Member Functions

 DOFHandlerMultiDim (MeshBase &_mesh, bool make_elem_part=true)
 Constructor. More...
 
void distribute_dofs (std::shared_ptr< DiscreteSpace > ds)
 Distributes degrees of freedom on the mesh needed for the given discrete space. More...
 
std::shared_ptr< DOFHandlerMultiDimsequential ()
 Returns sequential version of the current dof handler. More...
 
std::shared_ptr< VecScatter > sequential_scatter ()
 Returns scatter context from parallel to sequential vectors. More...
 
virtual VectorMPI create_vector ()
 Allocates PETSc vector according to the dof distribution. More...
 
LongIdx edge_index (int loc_edg) const
 Returns the global index of local edge. More...
 
LongIdx nb_index (int loc_nb) const
 Returns the global index of local neighbour. More...
 
unsigned int n_loc_edges () const
 Returns number of local edges. More...
 
unsigned int n_loc_nb () const
 Returns number of local neighbours. More...
 
void print () const
 Output structure of dof handler. More...
 
std::size_t hash () const override
 
Range< DHCellAccessorown_range () const
 Returns range of DOF handler cells (only range of own without ghost cells) More...
 
Range< DHCellAccessorlocal_range () const
 Returns range over own and ghost cells of DOF handler. More...
 
Range< DHCellAccessorghost_range () const
 Returns range over ghosts DOF handler cells. More...
 
unsigned int n_own_cells () const
 Return size of own range (number of own cells) More...
 
unsigned int n_local_cells () const
 Return size of local range (number of local cells) More...
 
unsigned int n_ghost_cells () const
 Return size of ghost range (number of ghost cells) More...
 
const DHCellAccessor cell_accessor_from_element (unsigned int elm_idx) const
 Return DHCellAccessor appropriate to ElementAccessor of given idx. More...
 
std::shared_ptr< DiscreteSpaceds () const
 Return pointer to discrete space for which the handler distributes dofs. More...
 
const std::vector< LongIdx > & get_local_to_global_map () const
 Get the map between local dof indices and the global ones. More...
 
 ~DOFHandlerMultiDim () override
 Destructor. More...
 
- Public Member Functions inherited from DOFHandlerBase
 DOFHandlerBase (MeshBase &_mesh)
 Constructor. More...
 
unsigned int n_global_dofs () const
 Getter for the number of all mesh dofs required by the given finite element. More...
 
unsigned int lsize () const
 Returns the number of dofs on the current process. More...
 
unsigned int max_elem_dofs () const
 Returns max. number of dofs on one element. More...
 
std::shared_ptr< Distributiondistr () const
 
MeshBasemesh () const
 Returns the mesh. More...
 
virtual ~DOFHandlerBase ()
 Destructor. More...
 

Protected Member Functions

bool el_is_local (int index) const
 
void make_elem_partitioning ()
 Prepare parallel distribution of elements, edges and neighbours. More...
 
void init_cell_starts ()
 Initialize vector of starting indices for elements. More...
 
void init_dof_starts (std::vector< LongIdx > &node_dof_starts, std::vector< LongIdx > &edge_dof_starts)
 Initialize auxiliary vector of starting indices of nodal/edge dofs. More...
 
void init_status (std::vector< short int > &node_status, std::vector< short int > &edge_status)
 Initialize node_status and edge_status. More...
 
void receive_ghost_dofs (unsigned int proc, std::vector< LongIdx > &dofs)
 Obtain dof numbers on ghost elements from other processor. More...
 
void send_ghost_dofs (unsigned int proc)
 Send dof numbers to other processor. More...
 
void update_local_dofs (unsigned int proc, const std::vector< bool > &update_cells, const std::vector< LongIdx > &dofs, const std::vector< LongIdx > &node_dof_starts, std::vector< LongIdx > &node_dofs, const std::vector< LongIdx > &edge_dof_starts, std::vector< LongIdx > &edge_dofs)
 Update dofs on local elements from ghost element dofs. More...
 
void create_sequential ()
 Communicate local dof indices to all processors and create new sequential dof handler. More...
 
unsigned int get_dof_indices (const DHCellAccessor &cell, std::vector< LongIdx > &indices) const override
 Returns the global indices of dofs associated to the cell. More...
 
LocDofVec get_loc_dof_indices (unsigned int loc_ele_idx) const override
 Returns the indices of dofs associated to the cell on the local process. More...
 

Protected Attributes

std::shared_ptr< DiscreteSpaceds_
 Pointer to the discrete space for which the handler distributes dofs. More...
 
bool is_parallel_
 Indicator for parallel/sequential dof handler. More...
 
std::shared_ptr< DOFHandlerMultiDimdh_seq_
 Sequential dof handler associated to the current (parallel) one. More...
 
std::shared_ptr< VecScatter > scatter_to_seq_
 Scatter context for parallel to sequential vectors. More...
 
std::vector< LongIdxcell_starts
 Starting indices for local (owned+ghost) element dofs. More...
 
std::vector< IntIdxdof_indices
 Dof numbers on local and ghost elements. More...
 
std::vector< LongIdxlocal_to_global_dof_idx_
 Maps local and ghost dof indices to global ones. More...
 
std::unordered_map< LongIdx, LongIdxglobal_to_local_el_idx_
 Maps global element index into local/ghost index (obsolete). More...
 
Distributionel_ds_
 Distribution of elements. More...
 
vector< LongIdxedg_4_loc
 Local edge index -> global edge index. More...
 
vector< LongIdxnb_4_loc
 Local neighbour index -> global neighbour index. More...
 
vector< LongIdxghost_4_loc
 Indices of ghost cells (neighbouring with local elements). More...
 
set< unsigned int > ghost_proc
 Processors of ghost elements. More...
 
map< unsigned int, vector< LongIdx > > ghost_proc_el
 Arrays of ghost cells for each neighbouring processor. More...
 
bool distribute_edge_dofs
 Temporary flag which prevents using dof handler on meshes where edges are not allocated (currently BCMesh). More...
 
- Protected Attributes inherited from DOFHandlerBase
unsigned int n_global_dofs_
 Number of global dofs assigned by the handler. More...
 
unsigned int lsize_
 Number of dofs associated to local process. More...
 
unsigned int loffset_
 Index of the first dof on the local process. More...
 
unsigned int max_elem_dofs_
 Max. number of dofs per element. More...
 
MeshBasemesh_
 Pointer to the mesh to which the dof handler is associated. More...
 
std::shared_ptr< Distributiondof_ds_
 Distribution of dofs associated to local process. More...
 

Static Protected Attributes

static const int INVALID_NFACE = 1
 
static const int VALID_NFACE = 2
 
static const int ASSIGNED_NFACE = 3
 
static const int INVALID_DOF = -1
 

Private Member Functions

void print_cell_dofs (stringstream &s, DHCellAccessor cell) const
 Auxiliary function for output of information about dofs in a cell. More...
 

Friends

class DHCellAccessor
 
class DHCellSide
 
class DHNeighbSide
 
class SubDOFHandlerMultiDim
 

Detailed Description

Provides the numbering of the finite element degrees of freedom on the computational mesh.

Class DOFHandlerMultiDim distributes the degrees of freedom (dof) for a particular triplet of 1d, 2d and 3d finite elements on the computational mesh and provides mappings between local and global dofs. The template parameter dim denotes the spatial dimension of the reference finite element.

Currently the functionality is restricted to finite elements with internal and nodal dofs, i.e. the neighboring elements can share only dofs on nodes.

Definition at line 151 of file dofhandler.hh.

Constructor & Destructor Documentation

◆ DOFHandlerMultiDim()

DOFHandlerMultiDim::DOFHandlerMultiDim ( MeshBase _mesh,
bool  make_elem_part = true 
)

Constructor.

Parameters
_meshThe mesh.
make_elem_partAllow switch off make_element_partitioning, necessary for boundary DOF handler.

Definition at line 49 of file dofhandler.cc.

◆ ~DOFHandlerMultiDim()

DOFHandlerMultiDim::~DOFHandlerMultiDim ( )
override

Destructor.

Definition at line 579 of file dofhandler.cc.

Member Function Documentation

◆ cell_accessor_from_element()

const DHCellAccessor DOFHandlerMultiDim::cell_accessor_from_element ( unsigned int  elm_idx) const

Return DHCellAccessor appropriate to ElementAccessor of given idx.

Definition at line 701 of file dofhandler.cc.

Here is the caller graph for this function:

◆ create_sequential()

void DOFHandlerMultiDim::create_sequential ( )
protected

Communicate local dof indices to all processors and create new sequential dof handler.

Collective on all processors.

Definition at line 464 of file dofhandler.cc.

Here is the caller graph for this function:

◆ create_vector()

VectorMPI DOFHandlerMultiDim::create_vector ( )
virtual

Allocates PETSc vector according to the dof distribution.

Definition at line 553 of file dofhandler.cc.

◆ distribute_dofs()

void DOFHandlerMultiDim::distribute_dofs ( std::shared_ptr< DiscreteSpace ds)

Distributes degrees of freedom on the mesh needed for the given discrete space.

By default, the dof handler is parallel, meaning that each processor has access to dofs on the local elements and on one layer of ghost elements (owned by neighbouring elements). This can be changed by setting sequential to true.

Parameters
dsThe discrete space consisting of finite elements for each mesh element.
sequentialIf true then each processor will have information about all dofs.

Definition at line 324 of file dofhandler.cc.

◆ ds()

std::shared_ptr<DiscreteSpace> DOFHandlerMultiDim::ds ( ) const
inline

Return pointer to discrete space for which the handler distributes dofs.

Definition at line 255 of file dofhandler.hh.

Here is the caller graph for this function:

◆ edge_index()

LongIdx DOFHandlerMultiDim::edge_index ( int  loc_edg) const
inline

Returns the global index of local edge.

Parameters
loc_edgLocal index of edge.

Definition at line 200 of file dofhandler.hh.

◆ el_is_local()

bool DOFHandlerMultiDim::el_is_local ( int  index) const
protected

Returns true if element is on local process.

Parameters
indexGlobal element index.

Definition at line 669 of file dofhandler.cc.

Here is the caller graph for this function:

◆ get_dof_indices()

unsigned int DOFHandlerMultiDim::get_dof_indices ( const DHCellAccessor cell,
std::vector< LongIdx > &  indices 
) const
overrideprotectedvirtual

Returns the global indices of dofs associated to the cell.

Parameters
cellThe cell.
indicesArray of dof indices on the cell.

Implements DOFHandlerBase.

Definition at line 567 of file dofhandler.cc.

Here is the caller graph for this function:

◆ get_loc_dof_indices()

LocDofVec DOFHandlerMultiDim::get_loc_dof_indices ( unsigned int  loc_ele_idx) const
inlineoverrideprotectedvirtual

Returns the indices of dofs associated to the cell on the local process.

Parameters
loc_ele_idxlocal element index.

Implements DOFHandlerBase.

Definition at line 369 of file dofhandler.hh.

Here is the caller graph for this function:

◆ get_local_to_global_map()

const std::vector<LongIdx>& DOFHandlerMultiDim::get_local_to_global_map ( ) const
inline

Get the map between local dof indices and the global ones.

Definition at line 258 of file dofhandler.hh.

Here is the caller graph for this function:

◆ ghost_range()

Range< DHCellAccessor > DOFHandlerMultiDim::ghost_range ( ) const

Returns range over ghosts DOF handler cells.

Definition at line 694 of file dofhandler.cc.

Here is the caller graph for this function:

◆ hash()

std::size_t DOFHandlerMultiDim::hash ( ) const
overridevirtual

Implements DOFHandlerBase::hash.

Implements DOFHandlerBase.

Definition at line 675 of file dofhandler.cc.

◆ init_cell_starts()

void DOFHandlerMultiDim::init_cell_starts ( )
protected

Initialize vector of starting indices for elements.

Definition at line 83 of file dofhandler.cc.

Here is the caller graph for this function:

◆ init_dof_starts()

void DOFHandlerMultiDim::init_dof_starts ( std::vector< LongIdx > &  node_dof_starts,
std::vector< LongIdx > &  edge_dof_starts 
)
protected

Initialize auxiliary vector of starting indices of nodal/edge dofs.

Parameters
node_dof_startsVector of nodal starting indices (output).
edge_dof_startsVector of edge starting indices (output).

Definition at line 97 of file dofhandler.cc.

Here is the caller graph for this function:

◆ init_status()

void DOFHandlerMultiDim::init_status ( std::vector< short int > &  node_status,
std::vector< short int > &  edge_status 
)
protected

Initialize node_status and edge_status.

Set VALID_NFACE for nodes/edges owned by local elements and INVALID_NFACE for nodes/edges owned by ghost elements.

Parameters
node_statusVector of nodal status (output).
edge_statusVector of edge status (output).

Definition at line 123 of file dofhandler.cc.

Here is the caller graph for this function:

◆ local_range()

Range< DHCellAccessor > DOFHandlerMultiDim::local_range ( ) const

Returns range over own and ghost cells of DOF handler.

Definition at line 687 of file dofhandler.cc.

Here is the caller graph for this function:

◆ make_elem_partitioning()

void DOFHandlerMultiDim::make_elem_partitioning ( )
protected

Prepare parallel distribution of elements, edges and neighbours.

Definition at line 584 of file dofhandler.cc.

Here is the caller graph for this function:

◆ n_ghost_cells()

unsigned int DOFHandlerMultiDim::n_ghost_cells ( ) const
inline

Return size of ghost range (number of ghost cells)

Definition at line 247 of file dofhandler.hh.

◆ n_loc_edges()

unsigned int DOFHandlerMultiDim::n_loc_edges ( ) const
inline

Returns number of local edges.

Definition at line 212 of file dofhandler.hh.

◆ n_loc_nb()

unsigned int DOFHandlerMultiDim::n_loc_nb ( ) const
inline

Returns number of local neighbours.

Definition at line 217 of file dofhandler.hh.

◆ n_local_cells()

unsigned int DOFHandlerMultiDim::n_local_cells ( ) const
inline

Return size of local range (number of local cells)

Definition at line 242 of file dofhandler.hh.

◆ n_own_cells()

unsigned int DOFHandlerMultiDim::n_own_cells ( ) const
inline

Return size of own range (number of own cells)

Definition at line 237 of file dofhandler.hh.

◆ nb_index()

LongIdx DOFHandlerMultiDim::nb_index ( int  loc_nb) const
inline

Returns the global index of local neighbour.

Parameters
loc_nbLocal index of neighbour.

Definition at line 207 of file dofhandler.hh.

◆ own_range()

Range< DHCellAccessor > DOFHandlerMultiDim::own_range ( ) const

Returns range of DOF handler cells (only range of own without ghost cells)

Definition at line 680 of file dofhandler.cc.

Here is the caller graph for this function:

◆ print()

void DOFHandlerMultiDim::print ( ) const

Output structure of dof handler.

Definition at line 722 of file dofhandler.cc.

◆ print_cell_dofs()

void DOFHandlerMultiDim::print_cell_dofs ( stringstream &  s,
DHCellAccessor  cell 
) const
private

Auxiliary function for output of information about dofs in a cell.

Definition at line 707 of file dofhandler.cc.

Here is the caller graph for this function:

◆ receive_ghost_dofs()

void DOFHandlerMultiDim::receive_ghost_dofs ( unsigned int  proc,
std::vector< LongIdx > &  dofs 
)
protected

Obtain dof numbers on ghost elements from other processor.

Parameters
procNeighbouring processor.
dofsArray where dofs are stored (output).

Definition at line 178 of file dofhandler.cc.

Here is the caller graph for this function:

◆ send_ghost_dofs()

void DOFHandlerMultiDim::send_ghost_dofs ( unsigned int  proc)
protected

Send dof numbers to other processor.

Parameters
procNeighbouring processor.

Definition at line 199 of file dofhandler.cc.

Here is the caller graph for this function:

◆ sequential()

std::shared_ptr< DOFHandlerMultiDim > DOFHandlerMultiDim::sequential ( )

Returns sequential version of the current dof handler.

Collective on all processors.

Definition at line 69 of file dofhandler.cc.

◆ sequential_scatter()

std::shared_ptr< VecScatter > DOFHandlerMultiDim::sequential_scatter ( )

Returns scatter context from parallel to sequential vectors.

For sequential dof handler it returns null pointer. Collective on all processors.

Definition at line 76 of file dofhandler.cc.

◆ update_local_dofs()

void DOFHandlerMultiDim::update_local_dofs ( unsigned int  proc,
const std::vector< bool > &  update_cells,
const std::vector< LongIdx > &  dofs,
const std::vector< LongIdx > &  node_dof_starts,
std::vector< LongIdx > &  node_dofs,
const std::vector< LongIdx > &  edge_dof_starts,
std::vector< LongIdx > &  edge_dofs 
)
protected

Update dofs on local elements from ghost element dofs.

Parameters
procNeighbouring processor.
update_cellsVector of global indices of elements which need to be updated from ghost elements.
dofsVector of dof indices on ghost elements from processor proc.
node_dof_startsVector of starting indices of nodal dofs.
node_dofsVector of nodal dof indices (output).
edge_dof_startsVector of starting indices of edge dofs.
edge_dofsVector of edge dof indices (output).

Definition at line 230 of file dofhandler.cc.

Here is the caller graph for this function:

Friends And Related Function Documentation

◆ DHCellAccessor

friend class DHCellAccessor
friend

Definition at line 265 of file dofhandler.hh.

◆ DHCellSide

friend class DHCellSide
friend

Definition at line 266 of file dofhandler.hh.

◆ DHNeighbSide

friend class DHNeighbSide
friend

Definition at line 267 of file dofhandler.hh.

◆ SubDOFHandlerMultiDim

friend class SubDOFHandlerMultiDim
friend

Definition at line 268 of file dofhandler.hh.

Member Data Documentation

◆ ASSIGNED_NFACE

const int DOFHandlerMultiDim::ASSIGNED_NFACE = 3
staticprotected

Definition at line 393 of file dofhandler.hh.

◆ cell_starts

std::vector<LongIdx> DOFHandlerMultiDim::cell_starts
protected

Starting indices for local (owned+ghost) element dofs.

E.g. dof_indices[cell_starts[idx]] = dof number for first dof on the cell with local index idx within the parallel structure. To use with DHCellAccessor use the following:

DHCellAccessor<3> cell; ... // i-th dof number on the cell dof_indices[cell_starts[cell.local_idx()]+i] = ...

For sequential dof handler, dofs from all elements are stored and elements are ordered as in mesh_->get_row_4_el().

Definition at line 424 of file dofhandler.hh.

◆ dh_seq_

std::shared_ptr<DOFHandlerMultiDim> DOFHandlerMultiDim::dh_seq_
protected

Sequential dof handler associated to the current (parallel) one.

Definition at line 404 of file dofhandler.hh.

◆ distribute_edge_dofs

bool DOFHandlerMultiDim::distribute_edge_dofs
protected

Temporary flag which prevents using dof handler on meshes where edges are not allocated (currently BCMesh).

Definition at line 466 of file dofhandler.hh.

◆ dof_indices

std::vector<IntIdx> DOFHandlerMultiDim::dof_indices
protected

Dof numbers on local and ghost elements.

Dofs are ordered accordingly with cell_starts and local dof order given by the finite element. See cell_starts for more description.

Definition at line 432 of file dofhandler.hh.

◆ ds_

std::shared_ptr<DiscreteSpace> DOFHandlerMultiDim::ds_
protected

Pointer to the discrete space for which the handler distributes dofs.

Definition at line 398 of file dofhandler.hh.

◆ edg_4_loc

vector<LongIdx> DOFHandlerMultiDim::edg_4_loc
protected

Local edge index -> global edge index.

Definition at line 451 of file dofhandler.hh.

◆ el_ds_

Distribution* DOFHandlerMultiDim::el_ds_
protected

Distribution of elements.

Definition at line 448 of file dofhandler.hh.

◆ ghost_4_loc

vector<LongIdx> DOFHandlerMultiDim::ghost_4_loc
protected

Indices of ghost cells (neighbouring with local elements).

Definition at line 457 of file dofhandler.hh.

◆ ghost_proc

set<unsigned int> DOFHandlerMultiDim::ghost_proc
protected

Processors of ghost elements.

Definition at line 460 of file dofhandler.hh.

◆ ghost_proc_el

map<unsigned int, vector<LongIdx> > DOFHandlerMultiDim::ghost_proc_el
protected

Arrays of ghost cells for each neighbouring processor.

Definition at line 463 of file dofhandler.hh.

◆ global_to_local_el_idx_

std::unordered_map<LongIdx,LongIdx> DOFHandlerMultiDim::global_to_local_el_idx_
protected

Maps global element index into local/ghost index (obsolete).

Definition at line 445 of file dofhandler.hh.

◆ INVALID_DOF

const int DOFHandlerMultiDim::INVALID_DOF = -1
staticprotected

Definition at line 394 of file dofhandler.hh.

◆ INVALID_NFACE

const int DOFHandlerMultiDim::INVALID_NFACE = 1
staticprotected

Flags used during distribution of dofs to mark n-face and dof status.

INVALID_NFACE means that on this node/edge/side/cell the current processor will not distribute dofs. VALID_NFACE means that on this node/edge/side/cell the current processor will distribute dofs. ASSIGNED_NFACE means that dofs on this n-face are already distributed.

INVALID_DOF marks dofs whose global number has to be set by neighbouring processor.

Definition at line 391 of file dofhandler.hh.

◆ is_parallel_

bool DOFHandlerMultiDim::is_parallel_
protected

Indicator for parallel/sequential dof handler.

Definition at line 401 of file dofhandler.hh.

◆ local_to_global_dof_idx_

std::vector<LongIdx> DOFHandlerMultiDim::local_to_global_dof_idx_
protected

Maps local and ghost dof indices to global ones.

First lsize_ entries correspond to dofs owned by local processor, the remaining entries are ghost dofs sorted by neighbouring processor id.

Definition at line 440 of file dofhandler.hh.

◆ nb_4_loc

vector<LongIdx> DOFHandlerMultiDim::nb_4_loc
protected

Local neighbour index -> global neighbour index.

Definition at line 454 of file dofhandler.hh.

◆ scatter_to_seq_

std::shared_ptr<VecScatter> DOFHandlerMultiDim::scatter_to_seq_
protected

Scatter context for parallel to sequential vectors.

Definition at line 407 of file dofhandler.hh.

◆ VALID_NFACE

const int DOFHandlerMultiDim::VALID_NFACE = 2
staticprotected

Definition at line 392 of file dofhandler.hh.


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