14#include "arcane/utils/StringBuilder.h"
15#include "arcane/utils/Iostream.h"
16#include "arcane/utils/PlatformUtils.h"
17#include "arcane/utils/Convert.h"
18#include "arcane/utils/NotImplementedException.h"
19#include "arcane/utils/ArgumentException.h"
20#include "arcane/utils/FloatingPointExceptionSentry.h"
21#include "arcane/utils/ValueConvert.h"
23#include "arcane/ISubDomain.h"
24#include "arcane/IParallelMng.h"
25#include "arcane/ItemEnumerator.h"
26#include "arcane/IMesh.h"
27#include "arcane/IMeshSubMeshTransition.h"
28#include "arcane/ItemGroup.h"
29#include "arcane/Service.h"
30#include "arcane/Timer.h"
31#include "arcane/FactoryService.h"
32#include "arcane/ItemPrinter.h"
33#include "arcane/IItemFamily.h"
34#include "arcane/MeshVariable.h"
35#include "arcane/VariableBuildInfo.h"
37#include "arcane/std/MeshPartitionerBase.h"
38#include "arcane/std/MetisMeshPartitioner_axl.h"
41#define MPICH_SKIP_MPICXX
42#define OMPI_SKIP_MPICXX
46#if (PARMETIS_MAJOR_VERSION < 4)
47#error "Your version of Parmetis is too old .Version 4.0+ is required"
50typedef real_t realtype;
53#include "arcane/std/PartitionConverter.h"
54#include "arcane/std/GraphDistributor.h"
55#include "arcane/std/internal/MetisWrapper.h"
65using MetisCallStrategy = TypesMetisMeshPartitioner::MetisCallStrategy;
66using MetisEmptyPartitionStrategy = TypesMetisMeshPartitioner::MetisEmptyPartitionStrategy;
73class MetisMeshPartitioner
95 bool m_disable_floatingexception;
106 const String& Name)
const;
73class MetisMeshPartitioner {
…};
112MetisMeshPartitioner::
113MetisMeshPartitioner(
const ServiceBuildInfo& sbi)
114: ArcaneMetisMeshPartitionerObject(sbi)
117, m_disable_floatingexception(false)
119 m_parallel_mng = mesh()->parallelMng();
121 if (s==
"1" || s==
"TRUE")
122 m_disable_floatingexception =
true;
131 Int32 nb_part = m_parallel_mng->commSize();
144 if (m_disable_floatingexception){
163 bool force_partition =
false;
170 info() <<
"INFO: Unable to test load balancing on a single sub-domain";
178 bool force_full_repartition =
false;
180 force_full_repartition |= (m_nb_refine < 0);
182 if (nb_part != nb_rank) {
183 force_full_repartition =
true;
184 force_partition =
true;
187 info() <<
"WARNING: compensating the potential lack of 'Metis' options in case of manual instanciation";
188 Integer max_diffusion_count = 10;
189 Real imbalance_relative_tolerance = 4;
190 float tolerance_target = 1.05f;
191 bool dump_graph =
false;
192 MetisCallStrategy call_strategy = MetisCallStrategy::one_processor_per_node;
193 bool in_out_digest =
false;
198 if (call_strategy_env ==
"all-processors"){
199 call_strategy = MetisCallStrategy::all_processors;
200 }
else if (call_strategy_env ==
"one-processor-per-node") {
201 call_strategy = MetisCallStrategy::one_processor_per_node;
202 }
else if (call_strategy_env ==
"two-processors-two-nodes") {
203 call_strategy = MetisCallStrategy::two_processors_two_nodes;
204 }
else if (call_strategy_env ==
"two-gathered-processors") {
205 call_strategy = MetisCallStrategy::two_gathered_processors;
206 }
else if (call_strategy_env ==
"two-scattered-processors") {
207 call_strategy = MetisCallStrategy::two_scattered_processors;
208 }
else if (!call_strategy_env.
null()) {
209 ARCANE_FATAL(
"Invalid value '{0}' for ARCANE_METIS_CALL_STRATEGY environment variable",call_strategy_env);
213 in_out_digest = (v.value()!=0);
215 dump_graph = (v.value()!=0);
218 max_diffusion_count =
options()->maxDiffusiveCount();
219 imbalance_relative_tolerance =
options()->imbalanceRelativeTolerance();
220 tolerance_target = (float)(
options()->toleranceTarget());
221 dump_graph =
options()->dumpGraph();
222 call_strategy =
options()->metisCallStrategy();
223 in_out_digest =
options()->inputOutputDigest();
226 if (max_diffusion_count > 0)
227 force_full_repartition |= (m_nb_refine>=max_diffusion_count);
229 force_full_repartition |= (
imbalance()>1.0);
232 initConstraints(
false);
238 if (is_shared_memory)
239 call_strategy = MetisCallStrategy::one_processor_per_node;
241 idxtype nb_weight = nbCellWeight();
243 info() <<
"Load balancing with Metis nb_weight=" << nb_weight <<
" initial=" << initial_partition
244 <<
" call_strategy=" << (int)call_strategy
245 <<
" is_shared_memory?=" << is_shared_memory
246 <<
" disabling_fpe?=" << m_disable_floatingexception
247 <<
" sizeof(idxtype)==" <<
sizeof(idxtype);
250 initial_partition =
true;
261 Integer nb_own_cell = nbOwnCellsWithConstraints();
263 Int32 nb_empty_part = 0;
266 metis_vtkdist[0] = 0;
267 for(
Integer i=0; i<nb_rank; ++i ){
269 Int32 n = global_nb_own_cell[i];
273 metis_vtkdist[i+1] =
static_cast<idxtype
>(total_nb_cell);
278 info() <<
"Total nb_cell=" << total_nb_cell <<
" nb_empty_partition=" << nb_empty_part;
279 if (total_nb_cell==0){
280 info() <<
"INFO: mesh '" <<
mesh->name() <<
" has no cell. No partitioning is needed";
297 Integer nb_max_face_neighbour_cell = 0;
303 const Cell& item = *i_item;
304 if (cellUsedWithConstraints(item)){
305 cell_metis_uid[item] = mid;
308 nb_max_face_neighbour_cell += item.
nbFace();
310 cell_metis_uid.synchronize();
313 _initUidRef(cell_metis_uid);
316 cell_metis_uid.setUsed(
false);
319 metis_xadj.
reserve(nb_own_cell+1);
322 metis_adjncy.
reserve(nb_max_face_neighbour_cell);
333 if (!cellUsedWithConstraints(item))
336 metis_xadj.
add(metis_adjncy.
size());
338 getNeighbourCellsUidWithConstraints(item, neighbour_cells, &edge_weights);
340 metis_adjncy.
add(
static_cast<idxtype
>(neighbour_cells[z]));
342 metis_xadj.
add(metis_adjncy.
size());
346 idxtype numflags = 0;
347 idxtype nparts =
static_cast<int>(nb_part);
354 if (!s.
null() && (s.
upper() ==
"PARTICLES"))
360 bool cond = (nb_weight == 1);
376 metis_ubvec.
fill(tolerance_target);
381 cond = converter.isBalancable<realtype>(cells_weights, metis_ubvec, nb_part);
383 info() <<
"We have to increase imbalance :";
384 info() << metis_ubvec;
385 for (
auto& tol: metis_ubvec ) {
395 const bool do_print_weight =
false;
396 if (do_print_weight){
399 Integer iteration =
mesh->subDomain()->commonVariables().globalIteration();
406 for (
int i = 0; i < metis_xadj.
size()-1 ; ++i) {
407 dumpy <<
" Weight uid=" << i;
408 for(
int j=0 ; j < nb_weight ; ++j ){
409 dumpy <<
" w=" << *(metis_vwgt.
begin()+i*nb_weight+j)
410 <<
"(" << cells_weights[i*nb_weight+j] <<
")";
419 converter.computeContrib(edge_weights);
426 if (!initial_partition) {
427 cells_size = cellsSizeWithConstraints();
428 converter.computeContrib(cells_size, (
Real)itr);
433 idxtype metis_options[4];
434 metis_options[0] = 0;
442 metis_options[0] = 1;
443 metis_options[1] = 0;
444 metis_options[2] = m_random_seed;
445 metis_options[3] = 1;
451 idxtype metis_edgecut = 0;
456 std::chrono::high_resolution_clock clock;
457 auto start_time = clock.now();
469 bool redistribute =
true;
470 bool redistribute_by_wrapper =
false;
472 if (call_strategy == MetisCallStrategy::all_processors || call_strategy == MetisCallStrategy::two_scattered_processors) {
473 redistribute =
false;
476 if (call_strategy == MetisCallStrategy::two_processors_two_nodes || call_strategy == MetisCallStrategy::two_scattered_processors) {
477 redistribute_by_wrapper =
true;
494 bool gd_allow_only_one_rank =
false;
495 if (is_shared_memory || (nb_empty_part!=0 && initial_partition))
496 gd_allow_only_one_rank =
true;
501 if ((nb_empty_part+1)==nb_rank){
502 info() <<
"Initialize GraphDistributor with max rank=1";
503 gd.initWithMaxRank(1);
505 else if (call_strategy == MetisCallStrategy::two_gathered_processors && (nb_rank > 2)) {
507 info() <<
"Initialize GraphDistributor with max rank=2";
508 gd.initWithMaxRank(2);
511 info() <<
"Initialize GraphDistributor with one rank per node"
512 <<
" (allow_one?=" << gd_allow_only_one_rank <<
")";
518 info() <<
"Using GraphDistributor?=" << redistribute;
520 metis_xadj = gd.convert<idxtype>(metis_xadj, &metis_part,
true);
521 metis_vwgt = gd.convert<idxtype>(metis_vwgt);
522 metis_adjncy = gd.convert<idxtype>(metis_adjncy);
523 metis_ewgt = gd.convert<idxtype>(metis_ewgt);
524 if (!initial_partition)
525 metis_vsize = gd.convert<idxtype>(metis_vsize);
526 metis_pm = gd.subParallelMng();
528 metis_vtkdist.
resize(gd.size()+1);
529 if (gd.contribute()) {
531 Integer nbVertices = metis_part.size();
533 metis_vtkdist[0] = 0;
534 for (
int i = 1 ; i < metis_vtkdist.
size() ; ++i) {
535 metis_vtkdist[i] = metis_vtkdist[i-1] + buffer[i-1];
538 metis_options[3] = PARMETIS_PSR_UNCOUPLED ;
544 MPI_Comm metis_mpicomm =
static_cast<MPI_Comm
>(metis_pm->
communicator());
546 bool do_call_metis =
true;
548 do_call_metis = gd.contribute();
555 _writeGraph(metis_pm, metis_vtkdist, metis_xadj, metis_adjncy, metis_vwgt, metis_ewgt, name.
toString());
561 float fill_value = (float)(1.0/(
double)nparts);
562 tpwgt.
fill(fill_value);
564 int retval = METIS_OK;
568 force_partition |= initial_partition;
569 force_full_repartition |= force_partition;
570 if (force_full_repartition){
572 if (force_partition) {
573 info() <<
"Metis: use a complete partitioning.";
575 redistribute_by_wrapper,
576 metis_vtkdist.
data(),
581 &wgtflag,&numflags,&nb_weight,
582 &nparts, tpwgt.
data(),
584 metis_options,&metis_edgecut,
588 info() <<
"Metis: use a complete REpartitioning.";
590 redistribute_by_wrapper,
591 metis_vtkdist.
data(),
597 &wgtflag,&numflags,&nb_weight,
598 &nparts,tpwgt.
data(),
600 &itr, metis_options,&metis_edgecut,
608 info() <<
"Metis: use a diffusive REpartitioning";
609 retval = ParMETIS_V3_RefineKway(metis_vtkdist.
data(),
614 &wgtflag,&numflags,&nb_weight,
615 &nparts,tpwgt.
data(),
617 metis_options,&metis_edgecut,
621 if (retval != METIS_OK) {
622 ARCANE_FATAL(
"Error while computing ParMetis partitioning r='{0}'",retval);
626 metis_part = gd.convertBack<idxtype>(metis_part, nb_own_cell);
629 auto end_time = clock.now();
630 std::chrono::microseconds duration = std::chrono::duration_cast<std::chrono::microseconds>(end_time - start_time);
631 info() <<
"Time to partition using parmetis = " << duration.count() <<
" us ";
637 switch(
options()->emptyPartitionStrategy()){
638 case MetisEmptyPartitionStrategy::DoNothing:
640 case MetisEmptyPartitionStrategy::TakeFromBiggestPartitionV1:
643 case MetisEmptyPartitionStrategy::TakeFromBiggestPartitionV2:
656 if (!cellUsedWithConstraints(item))
661 changeCellOwner(item, cells_new_owner, new_owner);
668 cells_new_owner.synchronize();
701 for (
int i =0 ; i < nb_own_cell ; i++) {
704 pm->
computeMinMaxSum(elem_by_part, min_part, max_part, sum_part, min_roc, max_proc);
707 Int32 max_part_id = -1;
708 Int32 max_part_nbr = -1;
710 for(
int i = 0; i < nb_part ; i++) {
711 if (sum_part[i] == 0) {
714 if(max_part_nbr < sum_part[i]) {
715 max_part_nbr = sum_part[i];
719 info() <<
"Parmetis: number empty parts " << nb_hole;
723 if(my_rank == max_proc[max_part_id]) {
725 for(
int i = 0; i < nb_part ; i++) {
727 if (sum_part[i] == 0 && offset < nb_own_cell) {
728 while(offset < nb_own_cell && metis_part[offset] != max_part_id) {
732 if(offset == nb_own_cell)
735 metis_part[offset] = i;
756 const Int32 nb_own_cell = metis_part.
size();
764 for (
int i =0 ; i < nb_own_cell ; i++) {
765 elem_by_part[metis_part[i]]++;
767 pm->
computeMinMaxSum(elem_by_part, min_part, max_part, sum_part, min_rank, max_rank);
772 Int32 max_part_id = -1;
773 Int64 max_part_nbr = -1;
775 Int64 total_nb_cell = 0;
776 for(
int i = 0; i < nb_part ; i++) {
777 Int64 current_sum_part = sum_part[i];
778 total_nb_cell += current_sum_part;
779 if (current_sum_part == 0) {
780 empty_part_ranks.
add(i);
783 if (max_part_nbr < current_sum_part) {
784 max_part_nbr = current_sum_part;
789 ARCANE_FATAL(
"Bad value max_part_id ({0})",max_part_id);
790 info() <<
"Parmetis: check empty parts: (" << algo_iteration <<
") nb_empty_parts=" << nb_hole
791 <<
" nb_max_part=" << max_part_nbr
792 <<
" max_part_rank=" << max_part_id
793 <<
" max_proc_max_part_id=" << max_rank[max_part_id]
794 <<
" empty_part_ranks=" << empty_part_ranks
795 <<
" total_nb_cell=" << total_nb_cell;
802 if (my_rank == max_rank[max_part_id]) {
804 Int32 max_remove_cell = nb_own_cell - 1;
806 for(
int i = 0; i < nb_part ; i++) {
808 if (sum_part[i] == 0 && offset < nb_own_cell) {
809 while (offset < max_remove_cell && metis_part[offset] != max_part_id) {
813 if (offset == max_remove_cell)
816 metis_part[offset] = i;
836 bool do_continue =
true;
837 Int32 last_nb_hole = 0;
840 info() <<
"Parmetis: nb_empty_partition=" << nb_hole <<
" last_nb_partition=" << last_nb_hole;
847 if (last_nb_hole>0 && last_nb_hole<=nb_hole){
848 pwarning() <<
"Can not remove all empty partitions. This is probably because you try"
849 <<
" to cut in too many partitions";
852 last_nb_hole = nb_hole;
875 bool have_ewgt =
false;
879 info() <<
"COMM_SIZE=" << commsize <<
" RANK=" << commrank;
887#define METIS_ERROR ARCANE_FATAL("_writeGraph")
891 filename += commrank;
894 if (metis_vtkdist.
size() != commsize + 1)
897 nvtx = metis_vtkdist[commrank+1] - metis_vtkdist[commrank];
898 if (metis_xadj.
size() != nvtx + 1) {
899 std::cerr <<
"_writeGraph : nvtx+1 = " << nvtx <<
" != " << metis_xadj.
size() << std::endl;
904 nwgt = metis_vwgt.
size()/nvtx;
905 if (nwgt != 0 && metis_vwgt.
size() % nwgt != 0)
908 have_ewgt = (metis_ewgt.
size() != 0) ;
909 if (have_ewgt && metis_ewgt.
size() != metis_adjncy.
size())
915 Int64 localEdges = metis_xadj[nvtx];
918 file <<
"2" << std::endl;
919 file << commsize <<
"\t" << commrank << std::endl;
920 file << metis_vtkdist[commsize] <<
"\t" << globalEdges << std::endl;
921 file << nvtx <<
"\t" << localEdges << std::endl;
922 file <<
"0\t" <<
"0" << have_ewgt << nwgt << std::endl;
933 for (idxtype vertnum = 0 ; vertnum < nvtx ; vertnum++) {
935 for (
int dim = 0 ; dim < nwgt ; ++ dim)
936 file << metis_vwgt[vertnum*nwgt+dim] <<
'\t';
937 file << metis_xadj[vertnum + 1] - metis_xadj[vertnum];
938 for (idxtype edgenum = metis_xadj[vertnum] ;
939 edgenum < metis_xadj[vertnum + 1] ; ++edgenum) {
941 file <<
'\t' << metis_ewgt[edgenum];
942 file <<
'\t' << metis_adjncy[edgenum];
963#if ARCANE_DEFAULT_PARTITIONER == METIS_DEFAULT_PARTITIONER
#define ARCANE_THROW(exception_class,...)
Macro pour envoyer une exception avec formattage.
#define ARCANE_FATAL(...)
Macro envoyant une exception FatalErrorException.
#define ARCANE_SERVICE_INTERFACE(ainterface)
Macro pour déclarer une interface lors de l'enregistrement d'un service.
Integer size() const
Nombre d'éléments du vecteur.
ArcaneMetisMeshPartitionerObject(const Arcane::ServiceBuildInfo &sbi)
Constructeur.
CaseOptionsMetisMeshPartitioner * options() const
Options du jeu de données du service.
Exception lorsqu'un argument est invalide.
Conversion d'un tableau d'un type vers un autre type.
Vue modifiable d'un tableau d'un type T.
constexpr Integer size() const noexcept
Retourne la taille du tableau.
void fill(const DataType &data)
Remplissage du tableau.
ConstArrayView< T > constView() const
Vue constante sur ce tableau.
void resize(Int64 s)
Change le nombre d'éléments du tableau à s.
void reserve(Int64 new_capacity)
Réserve le mémoire pour new_capacity éléments.
const T * data() const
Accès à la racine du tableau hors toute protection.
iterator begin()
Itérateur sur le premier élément du tableau.
void add(ConstReferenceType val)
Ajoute l'élément val à la fin du tableau.
Int32 nbFace() const
Nombre de faces de la maille.
Int32 globalIteration() const
Numéro de l'itération courante.
Vue constante d'un tableau de type T.
Classe permettant d'activer/désactiver temporairement les exceptions flottantes du processeur.
Redistribute graph data to another "communicator".
void initWithOneRankPerNode(bool allow_only_one_rank)
Automatic distribution : one partitioning process per node.
Interface d'un partitionneur de maillage.
Interface d'un partitionneur de maillage.
Interface du gestionnaire de parallélisme pour un sous-domaine.
virtual bool isThreadImplementation() const =0
Indique si l'implémentation utilise les threads.
virtual void computeMinMaxSum(char val, char &min_val, char &max_val, char &sum_val, Int32 &min_rank, Int32 &max_rank)=0
Calcule en une opération la somme, le min, le max d'une valeur.
virtual Int32 commRank() const =0
Rang de cette instance dans le communicateur.
virtual Int32 commSize() const =0
Nombre d'instance dans le communicateur.
virtual void allGather(ConstArrayView< char > send_buf, ArrayView< char > recv_buf)=0
Effectue un regroupement sur tous les processeurs. Il s'agit d'une opération collective....
virtual Parallel::Communicator communicator() const =0
Communicateur MPI associé à ce gestionnaire.
virtual char reduce(eReduceType rt, char v)=0
Effectue la réduction de type rt sur le réel v et retourne la valeur.
virtual void barrier()=0
Effectue une barière.
Interface du gestionnaire d'un sous-domaine.
virtual const CommonVariables & commonVariables() const =0
Informations sur les variables standards.
virtual void flush()=0
Flush tous les flots.
@ PNoDump
Indique que la variable ne doit pas être sauvegardée.
Real imbalance() const override
Déséquilibre de temps de calcul.
Real maxImbalance() const override
Déséquilibre maximal autorisé
virtual void changeOwnersFromCells()
Positionne les nouveaux propriétaires des noeuds, arêtes et faces à partir des mailles.
Partitioneur de maillage utilisant la bibliothèque PARMetis.
void partitionMesh(bool initial_partition) override
void _removeEmptyPartsV2(Int32 nb_part, ArrayView< idxtype > metis_part)
Comble les partitions vides (version 2).
Int32 _removeEmptyPartsV2Helper(Int32 nb_part, ArrayView< idxtype > metis_part, Int32 algo_iteration)
Applique une itération de l'algorithme de suppression des partitions vides.
void build() override
Construction de niveau build du service.
void _partitionMesh(bool initial_partition, Int32 nb_part)
int _writeGraph(IParallelMng *pm, Span< const idxtype > metis_vtkdist, Span< const idxtype > metis_xadj, Span< const idxtype > metis_adjncy, Span< const idxtype > metis_vwgt, Span< const idxtype > metis_ewgt, const String &Name) const
void _removeEmptyPartsV1(const Int32 nb_part, const Int32 nb_own_cell, ArrayView< idxtype > metis_part)
Comble les partitions vides (version 1).
Wrapper autour des appels de Parmetis.
int callPartKway(const bool print_digest, const bool gather, idx_t *vtxdist, idx_t *xadj, idx_t *adjncy, idx_t *vwgt, idx_t *adjwgt, idx_t *wgtflag, idx_t *numflag, idx_t *ncon, idx_t *nparts, real_t *tpwgts, real_t *ubvec, idx_t *options, idx_t *edgecut, idx_t *part)
Simple wrapper autour de la routine ParMetis "ParMETIS_V3_PartKway".
int callAdaptiveRepart(const bool print_digest, const bool gather, idx_t *vtxdist, idx_t *xadj, idx_t *adjncy, idx_t *vwgt, idx_t *vsize, idx_t *adjwgt, idx_t *wgtflag, idx_t *numflag, idx_t *ncon, idx_t *nparts, real_t *tpwgts, real_t *ubvec, real_t *ipc2redist, idx_t *options, idx_t *edgecut, idx_t *part)
Simple wrapper autour de la routine ParMetis "ParMETIS_V3_AdaptiveRepart".
Conversion d'un tableau de flottants vers un tableau d'entiers/longs. \abstract Cette classe gere le ...
Structure contenant les informations pour créer un service.
Propriétés de création d'un service.
Vecteur 1D de données avec sémantique par référence.
constexpr __host__ __device__ SizeType size() const noexcept
Retourne la taille du tableau.
Vue d'un tableau d'éléments de type T.
Constructeur de chaîne de caractère unicode.
String toString() const
Retourne la chaîne de caractères construite.
Chaîne de caractères unicode.
bool null() const
Retourne true si la chaîne est nulle.
String upper() const
Transforme tous les caractères de la chaîne en majuscules.
const char * localstr() const
Retourne la conversion de l'instance dans l'encodage UTF-8.
Postionne le nom de l'action en cours d'exécution.
TraceMessage info() const
Flot pour un message d'information.
ITraceMng * traceMng() const
Gestionnaire de trace.
TraceMessage pwarning() const
Vecteur 1D de données avec sémantique par valeur (style STL).
Paramètres nécessaires à la construction d'une variable.
ItemGroupT< Cell > CellGroup
Groupe de mailles.
#define ARCANE_REGISTER_SERVICE(aclass, a_service_property,...)
Macro pour enregistrer un service.
MeshVariableScalarRefT< Cell, Integer > VariableCellInteger
Grandeur au centre des mailles de type entier.
ItemVariableScalarRefT< Int32 > VariableItemInt32
Grandeur de type entier 32 bits.
Integer toInteger(Real r)
Converti un Int64 en un Integer.
Int32 toInt32(Int64 v)
Converti un Int64 en un Int32.
Integer toInteger(Real r)
Converti un Real en Integer.
@ ReduceSum
Somme des valeurs.
@ ReduceMax
Maximum des valeurs.
-*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
UniqueArray< Int64 > Int64UniqueArray
Tableau dynamique à une dimension d'entiers 64 bits.
std::int64_t Int64
Type entier signé sur 64 bits.
Int32 Integer
Type représentant un entier.
@ ST_SubDomain
Le service s'utilise au niveau du sous-domaine.
@ IK_Cell
Entité de maillage de genre maille.
double Real
Type représentant un réel.
std::int32_t Int32
Type entier signé sur 32 bits.