21
21
#include "rheolef/config.h"
22
22
#ifdef _RHEOLEF_HAVE_MPI
23
23
#include "rheolef/geo.h"
24
#include "rheolef/geo_domain.h"
25
#include "rheolef/dis_macros.h"
26
#include "rheolef/rheostream.h"
27
#include "rheolef/index_set.h"
29
25
namespace rheolef {
32
array<polymorphic_array<geo_element>::size_type>
33
geo_partition (const polymorphic_array<geo_element>& elts, size_t dis_nv, size_t map_dim);
35
// --------------------------------------------------------------------------
36
// edges & faces renumbering subroutine
37
// --------------------------------------------------------------------------
38
template <class PolymorphicArray>
40
geo_element_renumbering (
41
const std::vector<geo_element::size_type>& new_global_vertex_owner,
42
const std::vector<geo_element::size_type>& global_new_num_vert,
43
const PolymorphicArray& ios_ge,
46
array<geo_element::size_type>& ios_ige2dis_ige
50
typedef typename geo_element::size_type size_type;
52
communicator comm = ios_ge.ownership().comm();
53
size_type dis_nge = ios_ge.dis_size();
54
vector<size_type> new_local_ge_owner (dis_nge, 0);
57
// 1) global all_reduce
58
// TODO: not balanced: ge nums
59
// TODO: also, not optimal: O(N) in communication, memory & CPU, instead of O(N/nproc)
61
for (size_type ige = 0, nge = ios_ge.size(); ige < nge; ige++) {
62
const geo_element& S = ios_ge [ige];
63
size_type ios_dis_ige = S.ios_dis_ie();
65
for (size_type iloc = 0; iloc < S.size(); iloc++) {
66
size_type ios_dis_iv = S [iloc];
67
assert_macro (ios_dis_iv < dis_nv, "vertex index "<<ios_dis_iv<<" is out of range [0:"<<dis_nv<<"[");
68
size_type iproc = new_global_vertex_owner [ios_dis_iv];
69
owner = std::max (owner,iproc);
71
new_local_ge_owner [ios_dis_ige] = owner;
73
vector<size_type> new_global_ge_owner (dis_nge, 0);
76
new_local_ge_owner.begin().operator->(),
78
new_global_ge_owner.begin().operator->(),
79
mpi::maximum<size_type>());
81
distributor ios_ge_ownership = ios_ge.ownership();
83
// 2) redistribute the vertex partition
85
array<size_type> ge_partition (ios_ge_ownership);
86
for (size_type ios_dis_ige = ge_partition.ownership().first_index(); ios_dis_ige < ge_partition.ownership().last_index(); ios_dis_ige++) {
87
ge_partition.dis_entry (ios_dis_ige) = new_global_ge_owner[ios_dis_ige];
89
ge_partition.dis_entry_assembly();
91
array<size_type> ige2ios_dis_ige;
99
// 3) vertices S[iloc] of new numbered element table still have ios numbering: fix it
101
size_type first_dis_ige = ge.ownership().first_index();
102
for (size_type ige = 0, nge = ge.size(); ige < nge; ige++) {
103
geo_element& S = ge [ige];
104
S.set_dis_ie (first_dis_ige + ige);
105
for (size_type iloc = 0, nloc = S.size(); iloc < nloc; iloc++) {
106
assert_macro (S[iloc] < dis_nv, "vertex index "<<S[iloc]<<" out of range [0:"<<dis_nv<<"[");
107
assert_macro (global_new_num_vert[S[iloc]] < dis_nv, "new vertex index "<<global_new_num_vert[S[iloc]] <<" out of range [0:"<<dis_nv<<"[");
108
S[iloc] = global_new_num_vert [S[iloc]];
112
// --------------------------------------------------------------------------
114
// --------------------------------------------------------------------------
117
geo_rep<T,distributed>::get (idiststream& ips)
120
check_macro (ips.good(), "bad input stream for geo.");
121
communicator comm = base::_geo_element[0].ownership().comm();
123
size_type io_proc = idiststream::io_proc();
124
size_type my_proc = comm.rank();
125
size_type nproc = comm.size();
126
if (my_proc == io_proc && !scatch(ips.is(),"\nmesh"))
127
error_macro("input stream does not contains a geo.");
128
// ------------------------------------------------------------------------
130
// ------------------------------------------------------------------------
136
size_type dis_nedg = 0;
137
size_type dis_n_fac = 0;
139
base::_name = "unnamed";
141
ips >> base::_version
146
if (base::_version < 2) {
147
warning_macro ("mesh version < 2 no more supported in distributed version");
149
if (base::_dimension >= 3) {
152
if (base::_dimension >= 2) {
156
std::fill (_ios_size_by_variant, _ios_size_by_variant+reference_element::max_size, 0);
158
// 1.2) get coordinates
160
array<vertex_type> ios_vertex (dis_nv);
161
ios_vertex.get_values (ips, _point_get<T>(base::_dimension));
162
check_macro (ips.good(), "bad input stream for geo.");
163
// set ios_dis_iv index as fisrt field of the idx_vertex pair:
164
distributor ios_vertex_ownership = ios_vertex.ownership();
165
polymorphic_array<geo_element> ios_geo_element_p (ios_vertex_ownership);
166
size_type first_ios_dis_iv = ios_vertex.ownership().first_index();
167
for (size_type ios_iv = 0, ios_nv = ios_vertex.size(); ios_iv < ios_nv; ios_iv++) {
168
size_type ios_dis_iv = first_ios_dis_iv + ios_iv;
169
ios_geo_element_p[ios_iv] = geo_element_p(ios_dis_iv);
170
geo_element& P = ios_geo_element_p[ios_iv];
171
P.set_ios_dis_ie (ios_dis_iv);
172
_ios_size_by_variant [P.variant()]++;
177
polymorphic_array<geo_element> ios_elts (dis_ne);
178
ios_elts.get_values (ips);
179
size_type ios_dis_ie_start = ios_elts.ownership().first_index();
180
base::_map_dimension = 0;
181
for (size_type ie = 0, ne = ios_elts.size(); ie < ne; ie++) {
182
geo_element& K = ios_elts [ie];
183
base::_map_dimension = std::max (K.dimension(), base::_map_dimension);
184
K.set_ios_dis_ie (ios_dis_ie_start + ie);
185
if (K.dimension() > 0) { // point already counted; used when dim>0 and map_dim=0
186
_ios_size_by_variant [K.variant()]++;
189
// merge: when np > n_element, map_dim=0 on some procs without any elts...
190
base::_map_dimension = mpi::all_reduce (comm, base::_map_dimension, mpi::maximum<size_type>());
192
// 1.4) get faces & edges
194
polymorphic_array<geo_element,distributed> ios_faces;
195
if (base::_version >= 2 && base::_dimension >= 3) {
196
ios_faces.resize (dis_n_fac);
197
ios_faces.get_values (ips);
198
size_type ios_dis_ifac_start = ios_faces.ownership().first_index();
199
for (size_type ifac = 0, nfac = ios_faces.size(); ifac < nfac; ifac++) {
200
geo_element& F = ios_faces [ifac];
201
F.set_ios_dis_ie (ios_dis_ifac_start + ifac);
202
_ios_size_by_variant [F.variant()]++;
205
polymorphic_array<geo_element,distributed> ios_edges;
206
if (base::_version >= 2 && base::_dimension >= 2) {
207
ios_edges.resize (dis_nedg);
208
ios_edges.get_values (ips);
209
size_type ios_dis_iedg_start = ios_edges.ownership().first_index();
210
for (size_type iedg = 0, nedg = ios_edges.size(); iedg < nedg; iedg++) {
211
geo_element& E = ios_edges [iedg];
212
E.set_ios_dis_ie (ios_dis_iedg_start + iedg);
213
_ios_size_by_variant [E.variant()]++;
216
// ------------------------------------------------------------------------
217
// 2) mesh partition & element renumbering
218
// ------------------------------------------------------------------------
219
array<size_type> partition = geo_partition (ios_elts, dis_nv, map_dimension());
221
// 2.1) elements renumbering
223
array<size_type> ie2ios_dis_ie; // no more used
224
ios_elts.repartition (
226
base::_geo_element[base::_map_dimension],
228
_ios_ige2dis_ige[base::_map_dimension]);
230
size_type first_dis_ie = base::_geo_element[base::_map_dimension].ownership().first_index();
231
for (size_type ie = 0, ne = size(); ie < ne; ie++) {
232
geo_element& K = operator[] (ie);
233
K.set_dis_ie (first_dis_ie + ie);
235
// ------------------------------------------------------------------------
236
// 3) vertices renumbering
237
// ------------------------------------------------------------------------
238
// 3.1) global all_reduce
239
// TODO: not balanced: vertex nums
240
// TODO: also, not optimal: O(N) in communication, memory & CPU, instead of O(N/nproc)
242
vector<size_type> new_local_vertex_owner (dis_nv, 0);
243
for (size_type ie = 0, ne = size(); ie < ne; ie++) {
244
const geo_element& K = operator[] (ie);
245
for (size_type iloc = 0; iloc < K.size(); iloc++) {
246
new_local_vertex_owner [K[iloc]] = my_proc;
249
vector<size_type> new_global_vertex_owner (dis_nv, 0);
252
new_local_vertex_owner.begin().operator->(),
254
new_global_vertex_owner.begin().operator->(),
255
mpi::maximum<size_type>());
257
// 3.2) redistribute the vertex partition
258
array<size_type> vertex_partition (ios_vertex.ownership());
259
for (size_type ios_dis_iv = vertex_partition.ownership().first_index(); ios_dis_iv < vertex_partition.ownership().last_index(); ios_dis_iv++) {
260
vertex_partition.dis_entry (ios_dis_iv) = new_global_vertex_owner[ios_dis_iv];
262
vertex_partition.dis_entry_assembly();
264
array<size_type> iv2ios_dis_iv;
265
ios_geo_element_p.repartition (
267
base::_geo_element[0],
269
_ios_ige2dis_ige[0]);
271
distributor vertex_ownership = base::_geo_element[0].ownership();
272
base::_vertex.resize (vertex_ownership);
273
ios_vertex.permutation_apply (
277
// 3.3) set the element[0] array
278
base::_geo_element[0].resize (vertex_ownership);
279
size_type first_dis_iv = vertex_ownership.first_index();
280
size_type last_dis_iv = vertex_ownership.last_index();
281
for (size_type iv = 0, nv = base::_geo_element[0].size(); iv < nv; iv++) {
282
geo_element& P = base::_geo_element [0][iv];
283
P[0] = first_dis_iv + iv;
284
P.set_dis_ie (first_dis_iv + iv);
287
// 3.3) vertices K[iloc] of new numbered element table K still have ios numbering: fix it
289
vector<size_type> local_new_num_vert (dis_nv, 0);
290
for (size_type dis_iv = _ios_ige2dis_ige[0].ownership().first_index(),
291
dis_nv = _ios_ige2dis_ige[0].ownership().last_index();
292
dis_iv < dis_nv; dis_iv++) {
293
size_type iv = dis_iv - _ios_ige2dis_ige[0].ownership().first_index();
294
local_new_num_vert [dis_iv] = _ios_ige2dis_ige[0][iv];
296
vector<size_type> global_new_num_vert (dis_nv, 0);
299
local_new_num_vert.begin().operator->(),
301
global_new_num_vert.begin().operator->(),
302
mpi::maximum<size_type>());
304
for (size_type ie = 0; ie < size(); ie++) {
305
geo_element& K = operator[] (ie);
306
for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
307
assert_macro (K[iloc] < dis_nv, "vertex index "<<K[iloc]<<" out of range [0:"<<dis_nv<<"[");
308
assert_macro (global_new_num_vert[K[iloc]] < dis_nv, "new vertex index "<<global_new_num_vert[K[iloc]] <<" out of range [0:"<<dis_nv<<"[");
309
K[iloc] = global_new_num_vert [K[iloc]];
312
// ------------------------------------------------------------------------
313
// 4) edge & face renumbering
314
// ------------------------------------------------------------------------
315
if (base::_version >= 2 && base::_map_dimension >= 2) {
316
geo_element_renumbering (
317
new_global_vertex_owner,
321
base::_geo_element[1],
322
_ios_ige2dis_ige[1]);
324
if (base::_version >= 2 && base::_map_dimension >= 3) {
325
geo_element_renumbering (
326
new_global_vertex_owner,
330
base::_geo_element[2],
331
_ios_ige2dis_ige[2]);
333
// ------------------------------------------------------------------------
334
// 6) get domain, until end-of-file (requires ios_ige2dis_ige renumbering)
335
// ------------------------------------------------------------------------
337
domain_indirect_basic<distributed> dom;
338
bool status = dom.get (ips, *this);
340
base::_domains.push_back (dom);
342
// ------------------------------------------------------------------------
343
// 7) external entities
344
// ------------------------------------------------------------------------
346
build_external_entities ();
347
set_element_edge_index();
348
set_element_face_index();
350
/** TODO: autres champs a initialiser :
352
_geo_element [dim] : pour dim=1,..,map_dim-1
355
_ios_ige2dis_ige[dim] : pour dim=1,..,map_dim
362
geo_rep<T,distributed>::reset_size_by ()
364
// _ios_size_by_variant[] may also be set !
365
geo_base_rep<T,distributed>::reset_size_by();
366
std::fill (_ios_size_by_dimension, _ios_size_by_dimension+4, 0);
367
std::fill (_dis_size_by_dimension, _dis_size_by_dimension+4, 0);
368
std::fill (_dis_size_by_variant, _dis_size_by_variant+reference_element::max_size, 0);
369
for (size_type dim = 0; dim <= base::_map_dimension; dim++) {
370
_ios_size_by_dimension [dim] = geo_element_ios_ownership(dim).size();
371
_dis_size_by_dimension [dim] = base::_geo_element [dim].ownership().dis_size();
373
mpi::all_reduce (base::comm(), base::_size_by_variant, reference_element::max_size, _dis_size_by_variant, std::plus<size_type>());
375
/** ------------------------------------------------------------------------
376
* loop on geo_element (edges, faces, etc):
377
* identify some vertices, that are referenced
378
* by locally-managed geo_elements, but these vertices are managed
379
* by another processor: e.g. vertices at a partition boundary.
380
* ------------------------------------------------------------------------
384
geo_rep<T,distributed>::build_external_entities ()
386
distributor vertex_ownership = base::_geo_element[0].ownership();
387
size_type first_dis_iv = vertex_ownership.first_index();
388
size_type last_dis_iv = vertex_ownership.last_index();
389
size_type dis_nv = vertex_ownership.dis_size();
391
// 1) list external vertex indexes:
392
std::set<size_type> ext_vertex_set;
393
for (size_type dim = 1; dim <= base::_map_dimension; dim++) {
394
for (size_type ige = 0, nge = base::_geo_element[dim].size(); ige < nge; ige++) {
395
const geo_element& K = base::_geo_element[dim][ige];
396
for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
397
assert_macro (K[iloc] < dis_nv, "vertex index "<<K[iloc]<<" out of range [0:"<<dis_nv<<"[");
398
size_type dis_iv = K [iloc];
399
if (dis_iv >= first_dis_iv && dis_iv < last_dis_iv) continue;
400
ext_vertex_set.insert (dis_iv);
404
warning_macro ("external: ext_iv_set.size="<<ext_vertex_set.size());
405
// 2) get external vertices:
406
base::_vertex.get_dis_entry (
410
base::_geo_element[0].get_dis_entry (
415
// ----------------------------------------------------------
416
// this paradigm appears very often: x[dis_i]
417
// TODO: embed it directly in array<T,distributed>
418
// ----------------------------------------------------------
419
template<class T, class M>
423
class ext_array<T,distributed> : public array<T,distributed> {
426
typedef array<T,distributed> base;
427
typedef typename base::size_type size_type;
428
typedef std::map<size_type,T> map_type;
430
ext_array(const distributor& ownership = distributor(), const T init_val = T())
431
: base(ownership,init_val), _ext_x() {}
432
void resize (const distributor& ownership = distributor(), const T init_val = T()) {
433
base::resize (ownership,init_val);
436
void set_external_indexes (const std::set<size_type>& ext_ix_set) {
438
base::get_dis_entry (ext_ix_set, _ext_x);
441
const T& dis_at (const size_type dis_i) const {
442
if (dis_i >= base::ownership().first_index() && dis_i < base::ownership().last_index()) {
443
size_type i = dis_i - base::ownership().first_index();
444
return base::operator[](i);
446
typename map_type::const_iterator iter = _ext_x.find (dis_i);
447
check_macro (iter != _ext_x.end(), "unexpected external index="<<dis_i);
448
return (*iter).second;
455
/** ------------------------------------------------------------------------
456
* on any 2d or 3d geo_element K, set K.dis_iedge(iloc) number
457
* ------------------------------------------------------------------------
461
geo_rep<T,distributed>::set_element_edge_index()
463
if (map_dimension() < 2) return;
464
warning_macro ("set edge...");
465
size_type nv = geo_element_ownership(0).size();
466
size_type first_dis_iv = geo_element_ownership(0).first_index();
467
size_type last_dis_iv = geo_element_ownership(0). last_index();
468
size_type first_dis_iedg = geo_element_ownership(1).first_index();
470
// ------------------------------------------------------------------------
471
// 1) ball(X) := { E; X is a vertex of E }
472
// ------------------------------------------------------------------------
473
warning_macro ("set edge [1]");
474
index_set empty_set; // TODO: add a global allocator to empty_set
475
array<index_set,distributed> ball (geo_element_ownership(0), empty_set);
476
std::set<size_t> ext_iv_set; // size=O((N/nproc)^((d-1)/d))
477
for (size_type iedg = 0, nedg = geo_element_ownership(1).size(); iedg < nedg; iedg++) {
478
const geo_element& E = get_geo_element (1, iedg);
479
size_type dis_iedg = first_dis_iedg + iedg;
480
for (size_type iloc = 0; iloc < 2; iloc++) {
481
size_type dis_iv = E[iloc];
482
if (dis_iv >= first_dis_iv && dis_iv < last_dis_iv) {
483
size_type iv = dis_iv - first_dis_iv;
484
ball [iv] += dis_iedg;
486
index_set dis_iedg_set;
487
dis_iedg_set += dis_iedg;
488
ball.dis_entry (dis_iv) += dis_iedg_set; // not so efficient
489
ext_iv_set.insert (dis_iv);
493
ball.dis_entry_assembly (index_set_add_op<index_set,index_set>());
494
// ------------------------------------------------------------------------
495
// 2) for all the dis_iv that are not handled by the current process
496
// but are referenced by at least an edge, get ext_ball[dis_iv]
497
// ------------------------------------------------------------------------
498
// Question: faudra-t'il inclure d'autres sommets dans ext_iv_set ?
499
// -> ceux issus des triangles, tetra, ect, et externes a la partition ?
500
// Question : est-ce que ca en ajouterait ? reponse : OUI (teste')
501
// mais c'est pas utile pour l'instant : ball sert juste aux aretes
502
warning_macro ("set edge [2] ext_iv_set.size="<<ext_iv_set.size());
504
ball.set_dis_indexes (ext_iv_set);
507
warning_macro ("set edge done");
509
/** ------------------------------------------------------------------------
510
* on any 3d geo_element K, set K.dis_iface(iloc) number
511
* ------------------------------------------------------------------------
515
geo_rep<T,distributed>::set_element_face_index()
517
if (map_dimension() < 3) return;
519
// ----------------------------------------------------------------------------
521
// ----------------------------------------------------------------------------
524
geo_rep<T,distributed>::load (
525
std::string filename,
526
const communicator& comm)
529
ips.open (filename, "geo", comm);
530
check_macro(ips.good(), "\"" << filename << "[.geo[.gz]]\" not found.");
532
std::string root_name = delete_suffix (delete_suffix(filename, "gz"), "geo");
533
std::string name = get_basename (root_name);
536
// ----------------------------------------------------------------------------
538
// ----------------------------------------------------------------------------
539
/// @brief helper permutation class for geo i/o
540
template <class V = typename polymorphic_traits<geo_element>::derived_type>
541
struct geo_element_perm {
542
typedef geo_element::size_type size_type;
543
geo_element_perm (const polymorphic_array<geo_element,distributed,V>& elts)
545
size_type operator[] (size_type ie) const {
546
const geo_element& K = _elts [ie];
547
return K.ios_dis_ie();
549
const polymorphic_array<geo_element,distributed,V>& _elts;
553
geo_rep<T,distributed>::put (odiststream& ops) const
556
communicator comm = base::_geo_element[0].ownership().comm();
557
size_type io_proc = odiststream::io_proc();
558
size_type my_proc = comm.rank();
559
size_type nproc = comm.size();
560
size_type dis_nv = base::_vertex.dis_size ();
561
size_type dis_ne = base::_geo_element[base::_map_dimension].dis_size ();
562
ops << "#!geo" << endl
566
<< " " << base::_dimension
569
if (base::_version >= 2) {
570
if (base::_dimension >= 3) {
571
ops << " " << base::_geo_element[2].dis_size();
573
if (base::_dimension >= 2) {
574
ops << " " << base::_geo_element[1].dis_size();
579
// build a permutationh array (could be avoided, but requires a template Permutation arg in array::permuted_put
580
array<size_type> iv2ios_dis_iv (base::_vertex.ownership());
581
for (size_type iv = 0, niv = iv2ios_dis_iv.size(); iv < niv; iv++) {
582
const geo_element& P = base::_geo_element[0][iv];
583
iv2ios_dis_iv [iv] = P.ios_dis_ie();
585
base::_vertex.permuted_put_values (ops, iv2ios_dis_iv, _point_put<T>(base::_dimension));
588
// elements are permuted back to original order and may
589
// refer to original vertices numbering
590
std::vector<size_type> vertex_perm ((my_proc == io_proc) ? dis_nv : 0);
591
size_type tag_gather = distributor::get_new_tag();
592
if (my_proc == io_proc) {
593
size_type i_start = iv2ios_dis_iv.ownership().first_index(my_proc);
594
size_type i_size = iv2ios_dis_iv.ownership().size (my_proc);
595
for (size_type i = 0; i < i_size; i++) {
596
vertex_perm [i_start + i] = iv2ios_dis_iv [i];
598
for (size_type iproc = 0; iproc < nproc; iproc++) {
599
if (iproc == my_proc) continue;
600
size_type i_start = iv2ios_dis_iv.ownership().first_index(iproc);
601
size_type i_size = iv2ios_dis_iv.ownership().size (iproc);
602
comm.recv (iproc, tag_gather, vertex_perm.begin().operator->() + i_start, i_size);
605
comm.send (0, tag_gather, iv2ios_dis_iv.begin().operator->(), iv2ios_dis_iv.size());
607
geo_element_permuted_put put_element (vertex_perm);
609
base::_geo_element[base::_map_dimension].permuted_put_values (
611
geo_element_perm<> (base::_geo_element[base::_map_dimension]),
617
if (base::_version >= 2 && base::_dimension >= 3) {
618
base::_geo_element[2].permuted_put_values (
620
geo_element_perm<> (base::_geo_element[2]),
623
if (base::_version >= 2 && base::_dimension >= 2) {
624
base::_geo_element[1].permuted_put_values (
626
geo_element_perm<> (base::_geo_element[1]),
632
for (typename std::vector<domain_indirect_basic<distributed> >::const_iterator
633
iter = base::_domains.begin(),
634
last = base::_domains.end();
635
iter != last; ++iter) {
637
(*iter).put (ops, *this);
643
geo_rep<T,distributed>::dump (std::string name) const {
644
base::_vertex.dump (name + "-vert");
645
base::_geo_element[base::_map_dimension].dump(name + "-elem");
647
27
// --------------------------------------------------------------------------
648
28
// accessors to distributed data
649
29
// --------------------------------------------------------------------------
650
30
template <class T>
31
typename geo_rep<T,distributed>::const_reference
32
geo_rep<T,distributed>::dis_get_geo_element (size_type dim, size_type dis_ige) const
34
if (base::_gs.ownership_by_dimension[dim].is_owned (dis_ige)) {
35
size_type first_dis_ige = base::_gs.ownership_by_dimension[dim].first_index();
36
size_type ige = dis_ige - first_dis_ige;
37
return get_geo_element (dim, ige);
39
// element is owned by another proc ; get its variant
40
size_type iproc = base::_gs.ownership_by_dimension[dim].find_owner(dis_ige);
41
size_type first_dis_ige = base::_gs.ownership_by_dimension[dim].first_index(iproc);
42
size_type first_dis_v = first_dis_ige;
43
size_type last_dis_v = first_dis_v;
45
for (size_type variant = reference_element::first_variant_by_dimension(dim);
46
variant < reference_element:: last_variant_by_dimension(dim); variant++) {
47
last_dis_v += base::_geo_element [variant].ownership().size (iproc);
48
if (dis_ige < last_dis_v) {
49
assert_macro (dis_ige >= shift, "unexpected index computation");
50
size_type dis_igev = dis_ige - shift;
51
return base::_geo_element [variant].dis_at (dis_igev);
53
shift += base::_geo_element [variant].ownership().first_index (iproc);
54
first_dis_v = last_dis_v;
56
error_macro ("geo_element dis_index " << dis_ige
57
<< " cause problem; its range is [0:"<< base::_gs.ownership_by_dimension[dim].dis_size() << "[");
58
return base::_geo_element [0].dis_at(0); // not reached
651
61
typename geo_rep<T,distributed>::size_type
652
62
geo_rep<T,distributed>::dis_ige2ios_dis_ige (size_type dim, size_type dis_ige) const
654
check_macro (dim == 0, "dimension="<<dim<<": not yet !");
655
if (dis_ige >= vertex_ownership().first_index()
656
&& dis_ige < vertex_ownership().last_index()) {
657
size_type ige = dis_ige - vertex_ownership().first_index();
658
const geo_element& P = base::_geo_element[0][ige];
659
return P.ios_dis_ie();
661
typename polymorphic_map<geo_element>::const_iterator iter = _ext_geo_element_0.find (dis_ige);
662
check_macro (iter != _ext_geo_element_0.end(), "unexpected vertex index:"<<dis_ige);
663
const geo_element& P = *iter;
664
return P.ios_dis_ie();
667
const typename geo_rep<T,distributed>::vertex_type&
668
geo_rep<T,distributed>::dis_vertex (size_type dis_iv) const
670
if (dis_iv >= vertex_ownership().first_index()
671
&& dis_iv < vertex_ownership().last_index()) {
672
size_type iv = dis_iv - vertex_ownership().first_index();
673
return base::_vertex [iv];
675
// here, dis_iv is not managed by current proc
676
// try on external associative table
677
typename vertex_map_type::const_iterator iter = _ext_vertex.find (dis_iv);
678
check_macro (iter != _ext_vertex.end(), "unexpected vertex index:"<<dis_iv);
679
return (*iter).second;
64
const geo_element& K = dis_get_geo_element(dim,dis_ige);
65
return K.ios_dis_ie();
681
67
// --------------------------------------------------------------------------
682
68
// access by geo_element(dim,idx)
added
removed
nfem/plib/geo_mpi.cc
