space_numbering.cc
Go to the documentation of this file.
#include "rheolef/space_numbering.h"
#include "rheolef/geo.h"
#include "rheolef/rheostream.h"
namespace rheolef { namespace space_numbering {
template <class T>
size_type
ndof (
const basis_basic<T>& b,
const geo_size& gs,
size_type map_dim)
{
size_type ndof = 0;
for (size_type variant = 0; variant < reference_element::max_variant; variant++) {
ndof += gs.ownership_by_variant [variant].size() * b.ndof_on_subgeo (map_dim, variant);
}
return ndof;
}
template <class T>
size_type
dis_ndof (
const basis_basic<T>& b,
const geo_size& gs,
size_type map_dim)
{
size_type dis_ndof = 0;
for (size_type variant = 0; variant < reference_element::max_variant; variant++) {
dis_ndof += gs.ownership_by_variant [variant].dis_size() * b.ndof_on_subgeo (map_dim, variant);
}
return dis_ndof;
}
template <class T>
size_type
nnod (
const basis_basic<T>& b,
const geo_size& gs,
size_type map_dim)
{
size_type nnod = 0;
for (size_type variant = 0; variant < reference_element::max_variant; variant++) {
nnod += gs.ownership_by_variant [variant].size() * b.nnod_on_subgeo (map_dim, variant);
}
return nnod;
}
template <class T>
size_type
dis_nnod (
const basis_basic<T>& b,
const geo_size& gs,
size_type map_dim)
{
size_type dis_nnod = 0;
for (size_type variant = 0; variant < reference_element::max_variant; variant++) {
dis_nnod += gs.ownership_by_variant [variant].dis_size() * b.nnod_on_subgeo (map_dim, variant);
}
return dis_nnod;
}
// --------------------------------------------------------------------------
// dis_inod/dis_idof shared internal function
// --------------------------------------------------------------------------
template <class Callback>
static
void
dis_ixxx (
const Callback& nxxx_on_subgeo,
const geo_size& gs,
const geo_element& K,
typename std::vector<size_t>::iterator dis_ixxx_tab)
{
typedef size_t size_type;
reference_element hat_K = K;
size_type K_map_dim = hat_K.dimension();
size_type n_comp = nxxx_on_subgeo.n_component();
size_type first_loc_ixxx = 0;
for (size_type subgeo_dim = 0; subgeo_dim <= K_map_dim; ++subgeo_dim) {
for (size_type subgeo_variant = reference_element::first_variant_by_dimension (subgeo_dim);
subgeo_variant < reference_element:: last_variant_by_dimension (subgeo_dim); ++subgeo_variant) {
if (subgeo_dim == K_map_dim && subgeo_variant != hat_K.variant()) continue;
reference_element hat_S (subgeo_variant);
size_type loc_ngev = hat_K.n_subgeo_by_variant (hat_S.variant());
if (loc_ngev == 0) continue;
size_type loc_nxxx_on_subgeo = nxxx_on_subgeo (hat_S.variant());
size_type loc_comp_nxxx_on_subgeo = loc_nxxx_on_subgeo / n_comp;
// => reduce the basis interface
if (loc_nxxx_on_subgeo == 0) continue;
// loop on all subgeos S like hat_S and that belongs to K
for (size_type loc_ige = 0, loc_nge = hat_K.n_subgeo (hat_S.dimension()); loc_ige < loc_nge; ++ loc_ige) {
reference_element hat_S2 = hat_K.subgeo (hat_S.dimension(), loc_ige);
if (hat_S2.variant() != hat_S.variant()) continue;
// 0) compute loc_igev
size_type loc_igev = hat_K.local_subgeo_index2index_by_variant (hat_S.variant(), loc_ige);
// 1) compute dis_igev
size_type dis_ige = K.subgeo_dis_index (hat_S.dimension(), loc_ige);
if (hat_S.dimension() == 0) { dis_ige = gs.dis_inod2dis_iv (dis_ige); }
size_type dis_igev = gs.dis_ige2dis_igev_by_variant (hat_S.variant(), dis_ige);
// 2) compute dis_ixxx_loc_start
size_type dis_ixxx_loc_start = 0;
size_type iproc = gs.ownership_by_variant [hat_S.variant()].find_owner (dis_igev);
for (size_type prev_subgeo_variant = 0;
prev_subgeo_variant < hat_S.variant();
prev_subgeo_variant++) {
dis_ixxx_loc_start
+= gs.ownership_by_variant [prev_subgeo_variant]. last_index(iproc)
* nxxx_on_subgeo (prev_subgeo_variant);
}
dis_ixxx_loc_start += dis_igev*nxxx_on_subgeo (hat_S.variant());
for (size_type next_subgeo_variant = hat_S.variant()+1;
next_subgeo_variant < reference_element::max_variant;
next_subgeo_variant++) {
dis_ixxx_loc_start
+= gs.ownership_by_variant [next_subgeo_variant].first_index(iproc)
* nxxx_on_subgeo (next_subgeo_variant);
}
// 3) set all ixxxs
size_type loc_ixxx_loc_start = first_loc_ixxx + loc_igev*loc_nxxx_on_subgeo;
for (size_type loc_comp_ixxx_on_subgeo = 0; loc_comp_ixxx_on_subgeo < loc_comp_nxxx_on_subgeo; ++loc_comp_ixxx_on_subgeo) {
size_type loc_comp_ixxx_on_subgeo_fixed = geo_element::fix_indirect (K, hat_S.variant(), loc_ige, loc_comp_ixxx_on_subgeo, nxxx_on_subgeo.degree());
for (size_type i_comp = 0; i_comp < n_comp; ++i_comp) {
size_type loc_ixxx = loc_ixxx_loc_start + loc_comp_ixxx_on_subgeo *n_comp + i_comp;
dis_ixxx_tab [loc_ixxx] = dis_ixxx_loc_start + loc_comp_ixxx_on_subgeo_fixed*n_comp + i_comp;
}
}
}
first_loc_ixxx += loc_ngev*loc_nxxx_on_subgeo;
}
}
}
// --------------------------------------------------------------------------
// dis_inod/dis_idof call-back
// --------------------------------------------------------------------------
// codes for these functions are very similar
// thus create a call-back that allows one to share one main code between them
template<class T>
struct nnod_callback {
nnod_callback (const basis_basic<T>& b, size_t map_d)
: _b(b), _map_d(map_d) {}
size_t operator() (size_t variant) const {
return _b.nnod_on_subgeo (_map_d, variant); }
size_t n_component() const { return 1; }
size_t degree() const { return _b.degree(); }
basis_basic<T> _b;
size_t _map_d;
};
template<class T>
struct ndof_callback {
ndof_callback (const basis_basic<T>& b, size_t map_d)
: _b(b), _map_d(map_d) {}
size_t operator() (size_t variant) const {
return _b.ndof_on_subgeo (_map_d, variant); }
size_t n_component() const { return _b.size(); }
size_t degree() const { return _b.degree(); }
basis_basic<T> _b;
size_t _map_d;
};
template <class T>
void
dis_inod (
const basis_basic<T>& b,
const geo_size& gs,
const geo_element& K,
typename std::vector<size_type>::iterator dis_inod_tab)
{
dis_ixxx (nnod_callback<T> (b, gs.map_dimension()), gs, K, dis_inod_tab);
}
template <class T>
void
dis_idof (
const basis_basic<T>& b,
const geo_size& gs,
const geo_element& K,
typename std::vector<size_type>::iterator dis_idof_tab)
{
dis_ixxx (ndof_callback<T> (b, gs.map_dimension()), gs, K, dis_idof_tab);
}
template <class T>
void
dis_inod (
const basis_basic<T>& b,
const geo_size& gs,
const geo_element& K,
std::vector<size_type>& dis_inod_tab)
{
dis_inod_tab.resize (b.nnod (K.variant()));
#ifdef _RHEOLEF_PARANO
std::fill (dis_inod_tab.begin(), dis_inod_tab.end(), std::numeric_limits<size_type>::max());
#endif // _RHEOLEF_PARANO
dis_inod (b, gs, K, dis_inod_tab.begin());
}
template <class T>
void
dis_idof (
const basis_basic<T>& b,
const geo_size& gs,
const geo_element& K,
std::vector<size_type>& dis_idof_tab)
{
dis_idof_tab.resize (b.ndof (K.variant()));
#ifdef _RHEOLEF_PARANO
std::fill (dis_idof_tab.begin(), dis_idof_tab.end(), std::numeric_limits<size_type>::max());
#endif // _RHEOLEF_PARANO
dis_idof (b, gs, K, dis_idof_tab.begin());
}
#ifdef _RHEOLEF_HAVE_MPI
// -------------------------------------------------------------
// permut to/from ios dof numbering, for distributed environment
// --------------------------------------------------------------------------
// set permutation for idofs: ios_idof2dis_idof & idof2ios_dis_idof
// and redistribute ios_dof[] into _dof[]
// applied for initialiing distributed mesh:
// posibly curved with Pk piola basis.
// --------------------------------------------------------------------------
template <class T>
void
generic_set_ios_permutation (
const basis_basic<T>& b,
size_t map_d,
const geo_size& gs,
const std::array<disarray<size_t,distributed>,reference_element::max_variant>&
igev2ios_dis_igev,
disarray<size_t,distributed>& idof2ios_dis_idof)
{
typedef geo_element::size_type size_type;
// 1) count & resize
size_type ndof = 0;
for (size_type variant = 0,
variant_last = reference_element::last_variant_by_dimension(map_d);
variant < variant_last;
variant++) {
ndof += gs.ownership_by_variant [variant].size() * b.ndof_on_subgeo (map_d, variant);
}
communicator comm = gs.ownership_by_variant[0].comm();
distributor dof_ownership (distributor::decide, comm, ndof);
idof2ios_dis_idof.resize (dof_ownership);
// 2) set idof2ios_dis_idof[]
size_type first_dis_idof= dof_ownership.first_index();
size_type idof = 0;
size_type first_ios_dis_v = 0;
for (size_type variant = 0,
variant_last = reference_element::last_variant_by_dimension(map_d);
variant < variant_last;
variant++) {
size_type loc_ndof = b.ndof_on_subgeo (map_d, variant);
for (size_type igev = 0, ngev = gs.ownership_by_variant [variant].size(); igev < ngev; igev++) {
size_type ios_dis_igev = igev2ios_dis_igev [variant] [igev];
for (size_type loc_idof = 0; loc_idof < loc_ndof; loc_idof++, idof++) {
size_type ios_dis_idof = first_ios_dis_v + ios_dis_igev*loc_ndof + loc_idof;
idof2ios_dis_idof [idof] = ios_dis_idof;
}
}
first_ios_dis_v += gs.ownership_by_variant [variant].dis_size() * loc_ndof;
}
}
template <class T>
void
set_ios_permutations (
const basis_basic<T>& b,
const geo_basic<T,distributed>& omega,
disarray<size_type,distributed>& idof2ios_dis_idof,
disarray<size_type,distributed>& ios_idof2dis_idof)
{
generic_set_ios_permutation (b, omega.map_dimension(), omega.sizes(), omega.get_igev2ios_dis_igev(), idof2ios_dis_idof);
size_type dis_ndof = idof2ios_dis_idof.ownership().dis_size();
distributor ios_dof_ownership (dis_ndof, idof2ios_dis_idof.comm(), distributor::decide);
ios_idof2dis_idof.resize (ios_dof_ownership, std::numeric_limits<size_type>::max());
idof2ios_dis_idof.reverse_permutation (ios_idof2dis_idof);
}
#endif // _RHEOLEF_HAVE_MPI
// ----------------------------------------------------------------------------
// instanciation in library
// ----------------------------------------------------------------------------
#define _RHEOLEF_instanciate(T) \
template size_type ndof ( \
const basis_basic<T>& b, \
const geo_size& gs, \
size_type map_dim); \
template size_type nnod ( \
const basis_basic<T>& b, \
const geo_size& gs, \
size_type map_dim); \
template size_type dis_ndof ( \
const basis_basic<T>& b, \
const geo_size& gs, \
size_type map_dim); \
template size_type dis_nnod ( \
const basis_basic<T>& b, \
const geo_size& gs, \
size_type map_dim); \
template void dis_idof ( \
const basis_basic<T>& b, \
const geo_size& gs, \
const geo_element& K, \
typename std::vector<size_type>::iterator dis_inod_tab); \
template void dis_inod ( \
const basis_basic<T>& b, \
const geo_size& gs, \
const geo_element& K, \
typename std::vector<size_type>::iterator dis_inod_tab); \
template void dis_idof ( \
const basis_basic<T>& b, \
const geo_size& gs, \
const geo_element& K, \
std::vector<size_type>& dis_inod_tab); \
template void dis_inod ( \
const basis_basic<T>& b, \
const geo_size& gs, \
const geo_element& K, \
std::vector<size_type>& dis_inod_tab); \
_RHEOLEF_instanciate(Float)
#ifdef _RHEOLEF_HAVE_MPI
#define _RHEOLEF_instanciate_distributed(T) \
template void generic_set_ios_permutation ( \
const basis_basic<T>& b, \
size_t map_d, \
const geo_size& gs, \
const std::array<disarray<size_t,distributed>,reference_element::max_variant>& \
igev2ios_dis_igev, \
disarray<size_t,distributed>& idof2ios_dis_idof); \
template void set_ios_permutations ( \
const basis_basic<T>& b, \
const geo_basic<T,distributed>& omega, \
disarray<size_type,distributed>& idof2ios_dis_idof, \
disarray<size_type,distributed>& ios_idof2dis_idof); \
_RHEOLEF_instanciate_distributed(Float)
#endif // _RHEOLEF_HAVE_MPI
}} // namespace rheolef::space_numbering
rheolef::space_numbering::ndof
size_type ndof(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
Definition: space_numbering.cc:28
rheolef::reference_element::last_variant_by_dimension
static variant_type last_variant_by_dimension(size_type dim)
Definition: reference_element.h:150
rheolef::space_numbering::generic_set_ios_permutation
void generic_set_ios_permutation(const basis_basic< T > &b, size_t map_d, const geo_size &gs, const std::array< disarray< size_t, distributed >, reference_element::max_variant > &igev2ios_dis_igev, disarray< size_t, distributed > &idof2ios_dis_idof)
Definition: space_numbering.cc:234
rheolef::geo_size::dis_inod2dis_iv
size_type dis_inod2dis_iv(size_type dis_inod) const
Definition: geo_size.cc:110
rheolef::geo_element::fix_indirect
static size_type fix_indirect(const geo_element &K, size_type subgeo_variant, size_type loc_ige, size_type loc_comp_idof_on_subgeo, size_type order)
Definition: geo_element.cc:320
rheolef::space_numbering::dis_ndof
size_type dis_ndof(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
Definition: space_numbering.cc:41
rheolef::geo_size
Definition: geo_size.h:29
rheolef::space_numbering::set_ios_permutations
void set_ios_permutations(const basis_basic< T > &b, const geo_basic< T, distributed > &omega, disarray< size_type, distributed > &idof2ios_dis_idof, disarray< size_type, distributed > &ios_idof2dis_idof)
Definition: space_numbering.cc:275
rheolef::distributor::comm
const communicator_type & comm() const
Definition: distributor.h:145
rheolef::space_constant::n_component
size_type n_component(valued_type valued_tag, size_type d, coordinate_type sys_coord)
Definition: space_constant.cc:136
rheolef::geo_basic< T, distributed >::sizes
const geo_size & sizes() const
Definition: geo.h:1424
rheolef::distributor::find_owner
size_type find_owner(size_type dis_i) const
find iproc associated to a global index dis_i: CPU=log(nproc)
Definition: distributor.cc:106
rheolef::reference_element::n_subgeo
size_type n_subgeo(size_type subgeo_dim) const
Definition: reference_element.h:110
rheolef::geo_element::subgeo_dis_index
size_type subgeo_dis_index(size_type subgeo_dim, size_type i) const
Definition: geo_element.h:214
rheolef::distributor
see the distributor page for the full documentation
Definition: distributor.h:62
rheolef::geo_element
see the geo_element page for the full documentation
Definition: geo_element.h:102
rheolef::distributor::decide
static const size_type decide
Definition: distributor.h:76
rheolef::size_type
size_t size_type
Definition: basis_get.cc:76
rheolef::space_numbering::dis_idof
void dis_idof(const basis_basic< T > &b, const geo_size &gs, const geo_element &K, typename std::vector< size_type >::iterator dis_idof_tab)
Definition: space_numbering.cc:187
rheolef::geo_element::variant
variant_type variant() const
Definition: geo_element.h:161
rheolef::reference_element::first_variant_by_dimension
static variant_type first_variant_by_dimension(size_type dim)
Definition: reference_element.h:148
rheolef::geo_basic< T, distributed >::get_igev2ios_dis_igev
const std::array< disarray< size_type, distributed >, reference_element::max_variant > & get_igev2ios_dis_igev() const
Definition: geo.h:1570
rheolef::reference_element::local_subgeo_index2index_by_variant
size_type local_subgeo_index2index_by_variant(size_type subgeo_variant, size_type i) const
Definition: reference_element.h:124
rheolef::basis_basic
Definition: basis.h:206
rheolef::reference_element
see the reference_element page for the full documentation
Definition: reference_element.h:66
rheolef::reference_element::subgeo
reference_element subgeo(size_type subgeo_dim, size_type loc_isid) const
Definition: reference_element.h:111
mkgeo_ball.variant
variant
Definition: mkgeo_ball.sh:149
rheolef::geo_basic< T, distributed >
distributed mesh with rerefence counting
Definition: geo.h:1367
rheolef::geo_size::map_dimension
size_type map_dimension() const
Definition: geo_size.h:41
rheolef::space_numbering::size_type
geo_element::size_type size_type
Definition: space_numbering.h:50
rheolef::geo_element::size_type
reference_element::size_type size_type
Definition: geo_element.h:125
rheolef::space_numbering::dis_nnod
size_type dis_nnod(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
Definition: space_numbering.cc:67
_RHEOLEF_instanciate_distributed
#define _RHEOLEF_instanciate_distributed(T)
rheolef::reference_element::variant
variant_type variant() const
Definition: reference_element.h:99
rheolef::geo_basic< T, distributed >::map_dimension
size_type map_dimension() const
Definition: geo.h:1410
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::distributor::dis_size
size_type dis_size() const
global and local sizes
Definition: distributor.h:207
rheolef::space_numbering::dis_inod
void dis_inod(const basis_basic< T > &b, const geo_size &gs, const geo_element &K, typename std::vector< size_type >::iterator dis_inod_tab)
Definition: space_numbering.cc:177
Float
see the Float page for the full documentation
rheolef::geo_size::ownership_by_variant
distributor ownership_by_variant[reference_element::max_variant]
Definition: geo_size.h:64
rheolef::reference_element::max_variant
static const variant_type max_variant
Definition: reference_element.h:82
rheolef::disarray
see the disarray page for the full documentation
Definition: disarray.h:459
rheolef::reference_element::dimension
size_type dimension() const
Definition: reference_element.h:101
mkgeo_ball.b
b
Definition: mkgeo_ball.sh:152
rheolef::space_numbering::nnod
size_type nnod(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
Definition: space_numbering.cc:54
size_type
field::size_type size_type
Definition: branch.cc:425
rheolef::geo_size::dis_ige2dis_igev_by_variant
size_type dis_ige2dis_igev_by_variant(size_type variant, size_type dis_ige) const
Definition: geo_size.cc:98
rheolef::distributor::first_index
size_type first_index(size_type iproc) const
global index range and local size owned by ip-th process
Definition: distributor.h:151
rheolef::_RHEOLEF_instanciate
_RHEOLEF_instanciate(Float, sequential) _RHEOLEF_instanciate(Float
mkgeo_ball.map_dim
map_dim
Definition: mkgeo_ball.sh:337
rheolef::reference_element::n_subgeo_by_variant
size_type n_subgeo_by_variant(size_type subgeo_variant) const
Definition: reference_element.h:107
rheolef::distributor::size
size_type size(size_type iproc) const
Definition: distributor.h:163