34 #include "rheolef/geo_domain.h"
35 #include "rheolef/space.h"
36 #include "rheolef/space_numbering.h"
37 #include "rheolef/space_mult.h"
38 #include "rheolef/piola_util.h"
39 #include "rheolef/basis_get.h"
47 template <
class T,
class M>
55 if (approx ==
"")
return;
56 const size_type unset = std::numeric_limits<basis_option::size_type>::max();
72 template <
class T,
class M>
77 !get_basis().have_compact_support_inside_element(),
78 "try to [un]block domain `" << act.
get_domain_name() <<
"' on discontinuous or bubble space("
79 << get_basis().
name()<<
"(" <<_omega.name()<<
"))");
80 _acts.push_back (act);
82 template <
class T,
class M>
87 if (! is_hierarchical()) {
89 scalar_constit.
do_act (act);
93 for (
iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
101 template <
class T,
class M>
106 const distributor& start_by_flattened_component)
const
108 check_macro (get_basis().is_initialized(),
"space with undefined finite element method cannot be used");
109 std::vector<size_type> comp_dis_idof_t;
110 distributor dof_ownership = blocked_flag.ownership();
115 for (
const_iterator iter = begin(), last = end(); iter != last; iter++) {
117 const std::string& dom_name = (*iter).get_domain_name();
118 size_type i_comp = (*iter).get_component_index();
120 "invalid blocked domain for the "<<i_comp <<
"-th component in a "
121 <<n_comp<<
"-component non-hierarchical space with basis=\""
122 << get_basis().
name() <<
"\"");
124 if (prev_act != act) {
125 blocked_flag.dis_entry_assembly();
130 distributor ige_ownership = _omega.geo_element_ownership (dom_dim);
134 for (
size_type ioige = 0, noige = dom.size(); ioige < noige; ioige++) {
135 size_type ige = dom.oige (ioige).index();
137 const geo_element& S = _omega.get_geo_element (dom_dim, ige);
144 loc_idof_start =
dim-1;
149 for (
size_type loc_idof = loc_idof_start, loc_ndof = comp_dis_idof_t.size(); loc_idof < loc_ndof; loc_idof += loc_idof_incr) {
150 size_type comp_dis_idof = comp_dis_idof_t [loc_idof];
153 assert_macro (comp_dis_idof >= first_comp_dis_idof,
"unexpected comp_dis_idof");
154 size_type comp_idof = comp_dis_idof - first_comp_dis_idof;
155 size_type comp_start_idof = start_by_flattened_component.
size(iproc);
156 size_type idof = comp_start_idof + comp_idof;
160 blocked_flag.dis_entry (
dis_idof) = blk;
163 basis scalar_basis = ...;
165 for (
size_type loc_scalar_idof = 0, loc_scalar_ndof = comp_dis_idof_t.size()/n_comp; loc_scalar_idof < loc_scalar_ndof; ++loc_scalar_idof) {
166 size_type dis_scalar_idof = scalar_dis_idof_t [loc_idof];
168 normal.dis_entry (dis_scalar_idof) = copy_of_a_previous_computation;
174 template <
class T,
class M>
178 const std::vector<distributor>& start_by_flattened_component,
181 if (! is_hierarchical()) {
183 scalar_constit.
data().build_blocked_flag (blocked_flag, ownership(), start_by_flattened_component [i_flat_comp]);
188 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
190 constit.
data().build_blocked_flag_recursive (blocked_flag, start_by_flattened_component, i_flat_comp);
193 template <
class T,
class M>
199 build_blocked_flag_recursive (blocked_flag, _start_by_flattened_component, i_flat_comp);
200 blocked_flag.dis_entry_assembly();
206 template <
class T,
class M>
210 if (! is_hierarchical()) {
211 return get_scalar().get_geo();
216 check_macro (_hier_constit.size() > 0,
"get_geo: empty space product");
218 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
221 if (ptr_comp_dom->name() == ptr_dom->name())
continue;
222 if (ptr_dom->get_background_geo().name() == ptr_comp_dom->name())
continue;
223 if (ptr_comp_dom->get_background_geo().name() == ptr_dom->name()) {
224 ptr_dom = ptr_comp_dom;
226 error_macro (
"get_geo: incompatible domains: \""<<ptr_dom->name()<<
"\" and \""<<ptr_comp_dom->name() <<
"\"");
231 template <
class T,
class M>
235 if (! is_hierarchical()) {
236 return get_scalar().get_background_geo();
239 check_macro (_hier_constit.size() > 0,
"get_background_geo: empty space product");
244 template <
class T,
class M>
248 check_macro(! is_hierarchical(),
"get_basis: undefined for heterogeneous space products");
249 return get_scalar().get_basis();
254 template <
class T,
class M>
258 if (! is_hierarchical()) {
259 return get_scalar().get_basis().have_compact_support_inside_element();
262 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
268 template <
class T,
class M>
272 if (! is_hierarchical()) {
273 return get_scalar().get_basis().is_discontinuous();
276 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
282 template <
class T,
class M>
286 : _is_initialized(false),
288 _start_by_flattened_component(),
289 _start_by_component(),
291 _valued_tag(space_constant::
scalar),
297 _loc_ndof.fill (std::numeric_limits<size_type>::max());
305 template <
class T,
class M>
307 : _is_initialized(false),
309 _start_by_flattened_component(),
310 _start_by_component(),
312 _valued_tag(space_constant::
mixed),
315 _hier_constit(
expr.size()),
318 _loc_ndof.fill (std::numeric_limits<size_type>::max());
322 *iter_constit = Xi.get_constitution();
326 template <
class T,
class M>
330 if (_is_hier != V2.
_is_hier) {
return false; }
333 if (_hier_constit.size() != V2.
_hier_constit.size())
return false;
334 for (
const_iterator iter1 = _hier_constit.begin(), iter2 = V2.
_hier_constit.begin(), last1 = _hier_constit.end(); iter1 != last1; ++iter1, ++iter2) {
335 if (! (*iter1).data().operator==((*iter2).data())) {
344 template <
class T,
class M>
348 if (! is_hierarchical())
return 1;
350 for (
size_type i_comp = 0, n_comp = size(); i_comp < n_comp; i_comp++) {
352 flattened_size += comp_constit.
data()._init_flattened_size();
354 return flattened_size;
356 template <
class T,
class M>
362 std::vector<distributor>& start_by_flattened_component)
const
364 if (! is_hierarchical()) {
365 if (! get_basis().is_initialized())
return;
366 start_by_flattened_component [i_flat_comp] =
distributor (dis_start_flat_comp_idof, comm(), start_flat_comp_idof);
370 start_flat_comp_idof +=
ndof;
371 dis_start_flat_comp_idof +=
dis_ndof;
375 for (
size_type i_comp = 0, n_comp = size(); i_comp < n_comp; i_comp++) {
377 comp_constit.
data()._init_start_by_flattened_component (i_flat_comp, start_flat_comp_idof, dis_start_flat_comp_idof, start_by_flattened_component);
380 template <
class T,
class M>
384 if (! is_hierarchical()) {
385 if (! get_basis().is_initialized())
return;
386 _start_by_component.resize (1);
387 _start_by_component [0] =
distributor (0, comm(), 0);
391 _start_by_component.resize (size());
394 for (
size_type i_comp = 0, n_comp = size(); i_comp < n_comp; i_comp++) {
396 _start_by_component [i_comp] =
distributor (comp_start_dis_idof, comm(), comp_start_idof);
397 comp_start_idof += comp_constit.
ownership(). size();
401 template <
class T,
class M>
405 if (_is_initialized)
return;
406 _is_initialized =
true;
408 _flattened_size = _init_flattened_size();
409 _start_by_flattened_component.resize (_flattened_size);
413 _init_start_by_flattened_component (i_flat_comp, start_flat_comp_idof, dis_start_flat_comp_idof, _start_by_flattened_component);
414 _ownership =
distributor (dis_start_flat_comp_idof, comm(), start_flat_comp_idof);
415 _init_start_by_component();
422 local_append_external_dis_idof (
423 const std::vector<size_t>& comp_dis_idof_tab,
427 std::set<size_t>& ext_dof_set)
431 for (
size_type loc_idof = 0, loc_ndof = comp_dis_idof_tab.size(); loc_idof < loc_ndof; loc_idof++) {
432 size_type comp_dis_idof = comp_dis_idof_tab [loc_idof];
433 assert_macro (comp_dis_idof < comp_dis_ndof,
"idof " << comp_dis_idof_tab [loc_idof] <<
" out of range[0:"<< comp_dis_ndof <<
"[");
434 if (! comp_ownership.
is_owned (comp_dis_idof)) {
437 assert_macro (comp_dis_idof >= comp_first_dis_idof,
"unexpected comp_dis_idof");
438 size_type comp_idof = comp_dis_idof - comp_first_dis_idof;
439 size_type comp_start_idof = start_by_flattened_component.
size(iproc);
440 size_type idof = comp_start_idof + comp_idof;
444 << first_dis_idof <<
":" << dof_ownership.
last_index() <<
"[");
449 template <
class T,
class M>
453 std::set<size_type>& ext_dof_set,
455 const distributor& start_by_flattened_component)
const
460 std::vector<size_type> comp_dis_idof_tab;
462 for (
size_type ie = 0, ne = omega.size(); ie < ne; ie++) {
464 assembly_dis_idof (omega, K, comp_dis_idof_tab);
465 local_append_external_dis_idof (comp_dis_idof_tab, comp_ownership, dof_ownership, start_by_flattened_component, ext_dof_set);
472 for (
auto iter : omega.get_external_geo_element_map(
variant)) {
474 assembly_dis_idof (omega, K, comp_dis_idof_tab);
475 local_append_external_dis_idof (comp_dis_idof_tab, comp_ownership, dof_ownership, start_by_flattened_component, ext_dof_set);
479 template <
class T,
class M>
482 std::set<size_type>& ext_dof_set,
484 const std::vector<distributor>& start_by_flattened_component,
487 if (! is_hierarchical()) {
489 append_external_dof (scalar_constit.
get_geo(), ext_dof_set, dof_ownership, start_by_flattened_component [i_flat_comp]);
494 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
495 (*iter).data().compute_external_dofs (ext_dof_set, dof_ownership, start_by_flattened_component, i_flat_comp);
498 template <
class T,
class M>
501 std::set<size_type>& ext_dof_set)
const
505 compute_external_dofs (ext_dof_set, ownership(), _start_by_flattened_component, i_flat_comp);
510 template<
class T,
class M>
514 if (! is_hierarchical()) {
516 if (!scalar_constit.
get_basis().is_initialized())
return;
519 if (level > 0)
out <<
"(";
521 for (
size_type i = 0,
n = get_hierarchy().size(); i <
n; i++) {
524 if (i+1 <
n) {
out <<
"*"; }
526 if (level > 0)
out <<
")";
529 template<
class T,
class M>
533 std::ostringstream ostrstr;
535 return ostrstr.str();
540 template <
class T,
class M>
544 if (! is_hierarchical()) {
545 return get_scalar().get_basis().degree();
548 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
550 degree_max = std::max(degree_max, constit.
degree_max());
554 template <
class T,
class M>
558 if (! is_hierarchical()) {
561 get_scalar().neighbour_guard();
565 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
573 template <
class T,
class M>
578 communicator comm =_scalar_constit.get_geo().comm();
582 return ios_ownership.
size();
585 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
591 template <
class T,
class M>
599 template <
class T,
class M>
612 _scalar_constit.set_ios_permutations (comp_idof2ios_dis_idof, comp_ios_idof2dis_idof);
614 size_type comp_ndof = comp_idof2ios_dis_idof.size();
615 size_type comp_dis_ndof = comp_idof2ios_dis_idof.dis_size();
618 for (
size_type comp_idof = 0; comp_idof < comp_ndof; comp_idof++) {
619 size_type comp_ios_dis_idof = comp_idof2ios_dis_idof [comp_idof];
620 size_type ios_dis_idof = comp_start_dis_idof + comp_ios_dis_idof;
621 size_type idof = comp_start_idof + comp_idof;
624 idof2ios_dis_idof [idof] = ios_dis_idof;
626 comp_start_idof += comp_ndof;
627 comp_start_dis_idof += comp_dis_ndof;
631 for (
const_iterator iter = _hier_constit.begin(), last = _hier_constit.end(); iter != last; ++iter) {
636 template <
class T,
class M>
644 _scalar_constit.set_ios_permutations (idof2ios_dis_idof, ios_idof2dis_idof);
651 communicator comm1 = comm();
655 idof2ios_dis_idof.resize (dof_ownership, std::numeric_limits<size_type>::max());
658 set_ios_permutation_recursion (idof2ios_dis_idof, comp_start_idof, comp_start_dis_idof);
663 distributor ios_dof_ownership (dis_ndof1, idof2ios_dis_idof.comm(), ios_ndof1);
664 ios_idof2dis_idof.resize (ios_dof_ownership, std::numeric_limits<size_type>::max());
665 idof2ios_dis_idof.reverse_permutation (ios_idof2dis_idof);
674 #ifdef _RHEOLEF_HAVE_MPI
field::size_type size_type
see the basis page for the full documentation
void set_coordinate_system(coordinate_type sc)
void set_dimension(size_type d)
size_type dimension() const
valued_type valued_tag() const
static std::string standard_naming(std::string family_name, size_t degree, const basis_option &sopt)
see the distributor page for the full documentation
size_type last_index(size_type iproc) const
size_type find_owner(size_type dis_i) const
find iproc associated to a global index dis_i: CPU=log(nproc)
size_type dis_size() const
global and local sizes
size_type size(size_type iproc) const
size_type first_index(size_type iproc) const
global index range and local size owned by ip-th process
bool is_owned(size_type dis_i, size_type iproc) const
true when dis_i in [first_index(iproc):last_index(iproc)[
static const size_type decide
std::allocator< int >::size_type size_type
see the geo_element page for the full documentation
static variant_type last_variant_by_dimension(size_type dim)
static variant_type first_variant_by_dimension(size_type dim)
const std::string & get_domain_name() const
static const size_type unset_index
scalar_type _scalar_constit
void append_external_dof(const geo_basic< T, M > &dom, std::set< size_type > &ext_dof_set, const distributor &dof_ownership, const distributor &start_by_component) const
disarray< size_type, M > build_blocked_flag() const
hierarchy_type _hier_constit
bool have_compact_support_inside_element() const
void set_ios_permutations(disarray< size_type, M > &idof2ios_dis_idof, disarray< size_type, M > &ios_idof2dis_idof) const
const basis_basic< T > & get_basis() const
size_type ios_ndof() const
void _init_start_by_component() const
const geo_basic< T, M > & get_background_geo() const
std::array< size_type, reference_element::max_variant > _loc_ndof
void put(std::ostream &out, size_type level=0) const
size_type _init_flattened_size() const
bool is_discontinuous() const
void build_blocked_flag_recursive(disarray< size_type, M > &blocked_flag, const std::vector< distributor > &start_by_component, size_type &i_comp) const
void _init_start_by_flattened_component(size_type &i_flat_comp, size_type &start_flat_comp_idof, size_type &dis_start_flat_comp_idof, std::vector< distributor > &start_by_flattened_component) const
void set_ios_permutation_recursion(disarray< size_type, M > &idof2ios_dis_idof, size_type &comp_start_idof, size_type &comp_start_dis_idof) const
hierarchy_type::iterator iterator
bool operator==(const space_constitution_rep< T, M > &V2) const
size_type degree_max() const
hierarchy_type::const_iterator const_iterator
hierarchy_type::size_type size_type
space_scalar_constitution< T, M > scalar_type
const geo_basic< T, M > & get_geo() const
void compute_external_dofs(std::set< size_type > &ext_dof_set) const
void neighbour_guard() const
void do_act(const space_act &act)
bool have_compact_support_inside_element() const
size_type degree_max() const
const geo_basic< T, M > & get_geo() const
bool is_discontinuous() const
const geo_basic< T, M > & get_background_geo() const
void set_ios_permutation_recursion(disarray< size_type, M > &idof2ios_dis_idof, size_type &comp_start_idof, size_type &comp_start_dis_idof) const
const distributor & ownership() const
size_type ios_ndof() const
void neighbour_guard() const
void do_act(const space_act &act)
rep::const_iterator const_iterator
container_type::size_type size_type
void set_ios_permutations(disarray< size_type, M > &idof2ios_dis_idof, disarray< size_type, M > &ios_idof2dis_idof) const
container_type::const_iterator const_iterator
void build_blocked_flag(disarray< size_type, M > &blocked_flag, const distributor &comp_ownership, const distributor &start_by_component) const
space_scalar_constitution_rep()
basis_basic< T > _fem_basis
void do_act(const space_act &act)
const geo_basic< T, M > & get_geo() const
void do_act(const space_act &act)
const basis_basic< T > & get_basis() const
#define assert_macro(ok_condition, message)
#define error_macro(message)
void get_geo(istream &in, my_geo &omega)
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
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)
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)
size_type dis_ndof(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
size_type ndof(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
This file is part of Rheolef.
void put(std::ostream &out, std::string name, const tiny_matrix< T > &a)
void basis_parse_from_string(const std::string &str, family_index_option_type &fio)
geo_basic< T, M > compact(const geo_basic< T, M > &gamma)
details::field_expr_v2_nonlinear_terminal_function< details::normal_pseudo_function< Float > > normal()
normal: see the expression page for the full documentation