19 #ifndef DOFHANDLER_HH_
20 #define DOFHANDLER_HH_
23 #include <unordered_map>
83 virtual std::size_t
hash()
const =0;
180 std::shared_ptr<DOFHandlerMultiDim>
sequential();
225 std::size_t
hash()
const override;
255 std::shared_ptr<DiscreteSpace>
ds()
const {
return ds_; }
375 return LocDofVec(mem_ptr, ndofs,
false,
false);
398 std::shared_ptr<DiscreteSpace>
ds_;
Cell accessor allow iterate over DOF handler cells.
Side accessor allows to iterate over sides of DOF handler cell.
Class allows to iterate over sides of neighbour.
unsigned int max_elem_dofs() const
Returns max. number of dofs on one element.
DOFHandlerBase(MeshBase &_mesh)
Constructor.
unsigned int max_elem_dofs_
Max. number of dofs per element.
MeshBase * mesh_
Pointer to the mesh to which the dof handler is associated.
unsigned int lsize_
Number of dofs associated to local process.
unsigned int n_global_dofs_
Number of global dofs assigned by the handler.
unsigned int n_global_dofs() const
Getter for the number of all mesh dofs required by the given finite element.
virtual unsigned int get_dof_indices(const DHCellAccessor &cell, std::vector< LongIdx > &indices) const =0
Fill vector of the global indices of dofs associated to the cell.
unsigned int loffset_
Index of the first dof on the local process.
std::shared_ptr< Distribution > dof_ds_
Distribution of dofs associated to local process.
MeshBase * mesh() const
Returns the mesh.
virtual std::size_t hash() const =0
Compute hash value of DOF handler.
virtual ~DOFHandlerBase()
Destructor.
std::shared_ptr< Distribution > distr() const
unsigned int lsize() const
Returns the number of dofs on the current process.
virtual LocDofVec get_loc_dof_indices(unsigned int loc_ele_idx) const =0
Returns a vector of the indices of dofs associated to the cell on the local process.
Provides the numbering of the finite element degrees of freedom on the computational mesh.
static const int ASSIGNED_NFACE
static const int INVALID_NFACE
vector< LongIdx > nb_4_loc
Local neighbour index -> global neighbour index.
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.
unsigned int n_own_cells() const
Return size of own range (number of own cells)
Distribution * el_ds_
Distribution of elements.
LongIdx nb_index(int loc_nb) const
Returns the global index of local neighbour.
std::shared_ptr< DiscreteSpace > ds_
Pointer to the discrete space for which the handler distributes dofs.
std::shared_ptr< VecScatter > scatter_to_seq_
Scatter context for parallel to sequential vectors.
unsigned int get_dof_indices(const DHCellAccessor &cell, std::vector< LongIdx > &indices) const override
Returns the global indices of dofs associated to the cell.
unsigned int n_loc_edges() const
Returns number of local edges.
Range< DHCellAccessor > own_range() const
Returns range of DOF handler cells (only range of own without ghost cells)
void make_elem_partitioning()
Prepare parallel distribution of elements, edges and neighbours.
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.
vector< LongIdx > ghost_4_loc
Indices of ghost cells (neighbouring with local elements).
bool distribute_edge_dofs
Temporary flag which prevents using dof handler on meshes where edges are not allocated (currently BC...
std::vector< LongIdx > cell_starts
Starting indices for local (owned+ghost) element dofs.
unsigned int n_ghost_cells() const
Return size of ghost range (number of ghost cells)
friend class SubDOFHandlerMultiDim
std::size_t hash() const override
LongIdx edge_index(int loc_edg) const
Returns the global index of local edge.
bool el_is_local(int index) const
vector< LongIdx > edg_4_loc
Local edge index -> global edge index.
bool is_parallel_
Indicator for parallel/sequential dof handler.
static const int INVALID_DOF
unsigned int n_local_cells() const
Return size of local range (number of local cells)
set< unsigned int > ghost_proc
Processors of ghost elements.
void distribute_dofs(std::shared_ptr< DiscreteSpace > ds)
Distributes degrees of freedom on the mesh needed for the given discrete space.
std::vector< IntIdx > dof_indices
Dof numbers on local and ghost elements.
std::shared_ptr< DOFHandlerMultiDim > sequential()
Returns sequential version of the current dof handler.
static const int VALID_NFACE
const std::vector< LongIdx > & get_local_to_global_map() const
Get the map between local dof indices and the global ones.
std::shared_ptr< DOFHandlerMultiDim > dh_seq_
Sequential dof handler associated to the current (parallel) one.
void print_cell_dofs(stringstream &s, DHCellAccessor cell) const
Auxiliary function for output of information about dofs in a cell.
virtual VectorMPI create_vector()
Allocates PETSc vector according to the dof distribution.
void init_cell_starts()
Initialize vector of starting indices for elements.
DOFHandlerMultiDim(MeshBase &_mesh, bool make_elem_part=true)
Constructor.
~DOFHandlerMultiDim() override
Destructor.
void create_sequential()
Communicate local dof indices to all processors and create new sequential dof handler.
const DHCellAccessor cell_accessor_from_element(unsigned int elm_idx) const
Return DHCellAccessor appropriate to ElementAccessor of given idx.
void print() const
Output structure of dof handler.
Range< DHCellAccessor > local_range() const
Returns range over own and ghost cells of DOF handler.
Range< DHCellAccessor > ghost_range() const
Returns range over ghosts DOF handler cells.
void init_status(std::vector< short int > &node_status, std::vector< short int > &edge_status)
Initialize node_status and edge_status.
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.
std::vector< LongIdx > local_to_global_dof_idx_
Maps local and ghost dof indices to global ones.
void send_ghost_dofs(unsigned int proc)
Send dof numbers to other processor.
std::shared_ptr< DiscreteSpace > ds() const
Return pointer to discrete space for which the handler distributes dofs.
std::unordered_map< LongIdx, LongIdx > global_to_local_el_idx_
Maps global element index into local/ghost index (obsolete).
unsigned int n_loc_nb() const
Returns number of local neighbours.
std::shared_ptr< VecScatter > sequential_scatter()
Returns scatter context from parallel to sequential vectors.
void receive_ghost_dofs(unsigned int proc, std::vector< LongIdx > &dofs)
Obtain dof numbers on ghost elements from other processor.
map< unsigned int, vector< LongIdx > > ghost_proc_el
Arrays of ghost cells for each neighbouring processor.
unsigned int lsize(int proc) const
get local size
Abstract class for the description of a general finite element on a reference simplex in dim dimensio...
Base class for Mesh and BCMesh.
unsigned int fe_idx_
Index of FE in parent FESystem.
std::vector< LongIdx > parent_dof_idx_
Local indices in the parent handler.
void update_parent_vector(VectorMPI &vec, const VectorMPI &subvec)
Update values in parent vector from values of subvector.
void receive_sub_ghost_dofs(unsigned int proc, vector< unsigned int > &ghost_dof_proc, map< LongIdx, LongIdx > &global_to_local_dof_idx)
Get global dof indices of ghost dofs for sub-handler.
void update_subvector(const VectorMPI &vec, VectorMPI &subvec)
Update values in subvector from parent vector.
void send_sub_ghost_dofs(unsigned int proc, const map< LongIdx, LongIdx > &global_to_local_dof_idx)
Send global indices of dofs that are ghost on other processors.
const std::vector< LongIdx > & parent_indices()
Local indices in the parent handler.
std::shared_ptr< DOFHandlerMultiDim > parent_
Parent dof handler.
Declaration of class which provides the finite element for every mesh cell.
Template Iter serves as general template for internal iterators.
arma::Col< IntIdx > LocDofVec
int LongIdx
Define type that represents indices of large arrays (elements, nodes, dofs etc.)
Implementation of range helper class.