Arcane  v3.16.2.0
Documentation développeur
Chargement...
Recherche...
Aucune correspondance
VtkMeshIOService.cc
1// -*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
2//-----------------------------------------------------------------------------
3// Copyright 2000-2025 CEA (www.cea.fr) IFPEN (www.ifpenergiesnouvelles.com)
4// See the top-level COPYRIGHT file for details.
5// SPDX-License-Identifier: Apache-2.0
6//-----------------------------------------------------------------------------
7/*---------------------------------------------------------------------------*/
8/* VtkMeshIOService.cc (C) 2000-2025 */
9/* */
10/* Lecture/Ecriture d'un maillage au format Vtk historique (legacy). */
11/*---------------------------------------------------------------------------*/
12/*---------------------------------------------------------------------------*/
13
14#include "arcane/utils/Collection.h"
15#include "arcane/utils/Enumerator.h"
16#include "arcane/utils/HashTableMap.h"
17#include "arcane/utils/IOException.h"
18#include "arcane/utils/ITraceMng.h"
19#include "arcane/utils/Iostream.h"
20#include "arcane/utils/OStringStream.h"
21#include "arcane/utils/ScopedPtr.h"
22#include "arcane/utils/StdHeader.h"
23#include "arcane/utils/String.h"
24#include "arcane/utils/ValueConvert.h"
25#include "arcane/utils/Real3.h"
26
27#include "arcane/core/BasicService.h"
28#include "arcane/core/FactoryService.h"
29#include "arcane/core/ICaseMeshReader.h"
30#include "arcane/core/IItemFamily.h"
31#include "arcane/core/IPrimaryMesh.h"
32#include "arcane/core/IMeshBuilder.h"
33#include "arcane/core/IMeshReader.h"
34#include "arcane/core/IMeshUtilities.h"
35#include "arcane/core/IMeshWriter.h"
36#include "arcane/core/IParallelMng.h"
37#include "arcane/core/IVariableAccessor.h"
38#include "arcane/core/IXmlDocumentHolder.h"
39#include "arcane/core/Item.h"
41#include "arcane/core/IVariableMng.h"
42#include "arcane/core/VariableTypes.h"
43#include "arcane/core/XmlNode.h"
44#include "arcane/core/XmlNodeList.h"
45#include "arcane/core/UnstructuredMeshAllocateBuildInfo.h"
46
47#include "arcane/core/internal/IVariableMngInternal.h"
48#include "arcane/core/internal/VtkCellTypes.h"
49
50/*---------------------------------------------------------------------------*/
51/*---------------------------------------------------------------------------*/
52
53namespace Arcane
54{
55
56/*---------------------------------------------------------------------------*/
57/*---------------------------------------------------------------------------*/
58
59class VtkFile;
60using namespace Arcane::VtkUtils;
61
62/*---------------------------------------------------------------------------*/
63/*---------------------------------------------------------------------------*/
101class VtkMeshIOService
102: public TraceAccessor
103{
104 public:
105
106 explicit VtkMeshIOService(ITraceMng* tm)
107 : TraceAccessor(tm)
108 {}
109
110 public:
111
112 void build() {}
113
114 public:
115
116 enum eMeshType
117 {
118 VTK_MT_Unknown,
119 VTK_MT_StructuredGrid,
120 VTK_MT_UnstructuredGrid
121 };
122
124 {
125 public:
126 };
127
129 : public VtkMesh
130 {
131 public:
132
133 int m_nb_x;
134 int m_nb_y;
135 int m_nb_z;
136 };
137
138 public:
139
140 bool readMesh(IPrimaryMesh* mesh, const String& file_name, const String& dir_name, bool use_internal_partition);
141
142 private:
143
144 bool _readStructuredGrid(IPrimaryMesh* mesh, VtkFile&, bool use_internal_partition);
145 bool _readUnstructuredGrid(IPrimaryMesh* mesh, VtkFile& vtk_file, bool use_internal_partition);
146 void _readCellVariable(IMesh* mesh, VtkFile& vtk_file, const String& name_str, Integer nb_cell);
147 void _readItemGroup(IMesh* mesh, VtkFile& vtk_file, const String& name_str, Integer nb_item,
148 eItemKind ik, ConstArrayView<Int32> local_id);
149 void _readNodeGroup(IMesh* mesh, VtkFile& vtk_file, const String& name, Integer nb_item);
150 void _createFaceGroup(IMesh* mesh, const String& name, Int32ConstArrayView faces_lid);
151 bool _readData(IMesh* mesh, VtkFile& vtk_file, bool use_internal_partition, eItemKind cell_kind,
152 Int32ConstArrayView local_id, Integer nb_node);
153 void _readNodesUnstructuredGrid(IMesh* mesh, VtkFile& vtk_file, Array<Real3>& node_coords);
155 Array<Integer>& cells_nb_node,
156 Array<ItemTypeId>& cells_type,
157 Array<Int64>& cells_connectivity);
158 void _readFacesMesh(IMesh* mesh, const String& file_name,
159 const String& dir_name, bool use_internal_partition);
160 bool _readMetadata(IMesh* mesh, VtkFile& vtk_file);
161
162 private:
163};
164
165/*---------------------------------------------------------------------------*/
166/*---------------------------------------------------------------------------*/
167
168class VtkFile
169{
170 public:
171
172 static const int BUFSIZE = 10000;
173
174 public:
175
176 explicit VtkFile(std::istream* stream)
177 : m_stream(stream)
178 , m_is_init(false)
180 , m_is_eof(false)
181 , m_is_binary_file(false)
182 , m_buf{}
183 {}
184
185 const char* getCurrentLine();
186 bool isEmptyNextLine();
187 const char* getNextLine();
188
189 void checkString(const String& current_value, const String& expected_value);
190 void checkString(const String& current_value,
191 const String& expected_value1,
192 const String& expected_value2);
193
194 static bool isEqualString(const String& current_value, const String& expected_value);
195
196 void reReadSameLine() { m_need_reread_current_line = true; }
197
198 bool isEof() { return m_is_eof; }
199
200 template <class T>
201 void getBinary(T& type);
202 float getFloat();
203 double getDouble();
204 int getInt();
205
206 void setIsBinaryFile(bool new_val) { m_is_binary_file = new_val; }
207
208 private:
209
211 std::istream* m_stream = nullptr;
212
215
218
221
224
226 char m_buf[BUFSIZE];
227};
228
229/*---------------------------------------------------------------------------*/
230/*---------------------------------------------------------------------------*/
231
237const char* VtkFile::
239{
240 if (!m_is_init)
241 getNextLine();
242 return m_buf;
243}
244
245/*---------------------------------------------------------------------------*/
246/*---------------------------------------------------------------------------*/
247
259{
260 m_is_init = true;
261
262 // On veut que getNextLine lise une nouvelle ligne.
263 // (on met à false dans le cas où cette méthode serai
264 // appelée plusieurs fois à la suite).
266
267 // Si l'on est arrivé à la fin du fichier lors du précédent appel de cette méthode ou
268 // de getNextLine, on envoie une erreur.
269 if (m_is_eof) {
270 throw IOException("VtkFile::isEmptyNextLine()", "Unexpected EndOfFile");
271 }
272
273 if (m_stream->good()) {
274 // Le getline s'arrete (par défaut) au char '\n' et ne l'inclus pas dans le buf
275 // mais le remplace par '\0'.
276 m_stream->getline(m_buf, sizeof(m_buf) - 1);
277
278 // Si on arrive au bout du fichier, on return true (pour dire oui, il y a une ligne vide,
279 // à l'appelant de gérer ça).
280 if (m_stream->eof()) {
281 m_is_eof = true;
282 return true;
283 }
284
285 // Sous Windows, une ligne vide commence par \r.
286 // getline remplace \n par \0, que ce soit sous Windows ou Linux.
287 if (m_buf[0] == '\r' || m_buf[0] == '\0') {
288 getNextLine();
289
290 // On demande à ce que le prochain appel à getNextLine renvoie la ligne
291 // qui vient tout juste d'être bufferisée.
293 return true;
294 }
295 else {
296 bool is_comment = true;
297
298 // On retire le commentaire, s'il y en a un, en remplaçant '#' par '\0'.
299 for (int i = 0; i < BUFSIZE && m_buf[i] != '\0'; ++i) {
300 if (!isspace(m_buf[i]) && m_buf[i] != '#' && is_comment) {
301 is_comment = false;
302 }
303 if (m_buf[i] == '#') {
304 m_buf[i] = '\0';
305 break;
306 }
307 }
308
309 // Si ce n'est pas un commentaire, on supprime juste le '\r' final (si windows).
310 if (!is_comment) {
311 // Supprime le '\r' final
312 for (int i = 0; i < BUFSIZE && m_buf[i] != '\0'; ++i) {
313 if (m_buf[i] == '\r') {
314 m_buf[i] = '\0';
315 break;
316 }
317 }
318 }
319
320 // Si c'était un commentaire, on recherche la prochaine ligne "valide"
321 // en appelant getNextLine.
322 else {
323 getNextLine();
324 }
325 }
327 return false;
328 }
329 throw IOException("VtkFile::isEmptyNextLine()", "Not Good");
330}
331
332/*---------------------------------------------------------------------------*/
333/*---------------------------------------------------------------------------*/
334
340const char* VtkFile::
342{
343 m_is_init = true;
344
345 // On return le buffer actuel, si celui-ci n'a pas été utilisé.
348 return getCurrentLine();
349 }
350
351 // Si l'on est arrivé à la fin du fichier lors du précédent appel de cette méthode ou
352 // de isEmptyNextLine, on envoie une erreur.
353 if (m_is_eof) {
354 throw IOException("VtkFile::isEmptyNextLine()", "Unexpected EndOfFile");
355 }
356
357 while (m_stream->good()) {
358 // Le getline s'arrete (par défaut) au char '\n' et ne l'inclus pas dans le buf mais le remplace par '\0'.
359 m_stream->getline(m_buf, sizeof(m_buf) - 1);
360
361 // Si on arrive au bout du fichier, on return le buffer avec \0 au début (c'est à l'appelant d'appeler
362 // isEof() pour savoir si le fichier est fini ou non).
363 if (m_stream->eof()) {
364 m_is_eof = true;
365 m_buf[0] = '\0';
366 return m_buf;
367 }
368
369 bool is_comment = true;
370
371 // Sous Windows, une ligne vide commence par \r.
372 // getline remplace \n par \0, que ce soit sous Windows ou Linux.
373 if (m_buf[0] == '\0' || m_buf[0] == '\r')
374 continue;
375
376 // On retire le commentaire, s'il y en a un, en remplaçant '#' par '\0'.
377 for (int i = 0; i < BUFSIZE && m_buf[i] != '\0'; ++i) {
378 if (!isspace(m_buf[i]) && m_buf[i] != '#' && is_comment) {
379 is_comment = false;
380 }
381 if (m_buf[i] == '#') {
382 m_buf[i] = '\0';
383 break;
384 }
385 }
386
387 // Si ce n'est pas un commentaire, on supprime juste le '\r' final (si windows).
388 if (!is_comment) {
389 // Supprime le '\r' final
390 for (int i = 0; i < BUFSIZE && m_buf[i] != '\0'; ++i) {
391 if (m_buf[i] == '\r') {
392 m_buf[i] = '\0';
393 break;
394 }
395 }
396 return m_buf;
397 }
398 }
399 throw IOException("VtkFile::getNextLine()", "Not good");
400}
401
402/*---------------------------------------------------------------------------*/
403/*---------------------------------------------------------------------------*/
404
411getFloat()
412{
413 float v = 0.;
414 if (m_is_binary_file) {
415 getBinary(v);
416 return v;
417 }
418 (*m_stream) >> ws >> v;
419
420 if (m_stream->good())
421 return v;
422
423 throw IOException("VtkFile::getFloat()", "Bad float");
424}
425
426/*---------------------------------------------------------------------------*/
427/*---------------------------------------------------------------------------*/
428
435getDouble()
436{
437 double v = 0.;
438 if (m_is_binary_file) {
439 getBinary(v);
440 return v;
441 }
442 (*m_stream) >> ws >> v;
443
444 if (m_stream->good())
445 return v;
446
447 throw IOException("VtkFile::getDouble()", "Bad double");
448}
449
450/*---------------------------------------------------------------------------*/
451/*---------------------------------------------------------------------------*/
452
459getInt()
460{
461 int v = 0;
462 if (m_is_binary_file) {
463 getBinary(v);
464 return v;
465 }
466 (*m_stream) >> ws >> v;
467
468 if (m_stream->good())
469 return v;
470
471 throw IOException("VtkFile::getInt()", "Bad int");
472}
473
474/*---------------------------------------------------------------------------*/
475/*---------------------------------------------------------------------------*/
476
482template <class T>
484getBinary(T& value)
485{
486 constexpr size_t sizeofT = sizeof(T);
487
488 // Le fichier VTK est en big endian et les CPU actuels sont en little endian.
489 Byte big_endian[sizeofT];
490 Byte little_endian[sizeofT];
491
492 // On lit les 'sizeofT' prochains octets que l'on met dans big_endian.
493 m_stream->read((char*)big_endian, sizeofT);
494
495 // On transforme le big_endian en little_endian.
496 for (size_t i = 0; i < sizeofT; i++) {
497 little_endian[sizeofT - 1 - i] = big_endian[i];
498 }
499
500 // On 'cast' la liste d'octet en type 'T'.
501 T* conv = new (little_endian) T;
502 value = *conv;
503}
504
505/*---------------------------------------------------------------------------*/
506/*---------------------------------------------------------------------------*/
507
518checkString(const String& current_value, const String& expected_value)
519{
520 String current_value_low = current_value.lower();
521 String expected_value_low = expected_value.lower();
522
523 if (current_value_low != expected_value_low) {
524 String s = "Expecting chain '" + expected_value + "', found '" + current_value + "'";
525 throw IOException("VtkFile::checkString()", s);
526 }
527}
528
529/*---------------------------------------------------------------------------*/
530/*---------------------------------------------------------------------------*/
531
543checkString(const String& current_value, const String& expected_value1, const String& expected_value2)
544{
545 String current_value_low = current_value.lower();
546 String expected_value1_low = expected_value1.lower();
547 String expected_value2_low = expected_value2.lower();
548
549 if (current_value_low != expected_value1_low && current_value_low != expected_value2_low) {
550 String s = "Expecting chain '" + expected_value1 + "' or '" + expected_value2 + "', found '" + current_value + "'";
551 throw IOException("VtkFile::checkString()", s);
552 }
553}
554
555/*---------------------------------------------------------------------------*/
556/*---------------------------------------------------------------------------*/
557
568isEqualString(const String& current_value, const String& expected_value)
569{
570 String current_value_low = current_value.lower();
571 String expected_value_low = expected_value.lower();
572 return (current_value_low == expected_value_low);
573}
574
575/*---------------------------------------------------------------------------*/
576/*---------------------------------------------------------------------------*/
577
578/*---------------------------------------------------------------------------*/
579/*---------------------------------------------------------------------------*/
580
591readMesh(IPrimaryMesh* mesh, const String& file_name, const String& dir_name, bool use_internal_partition)
592{
593 std::ifstream ifile(file_name.localstr(), std::ifstream::binary);
594
595 if (!ifile) {
596 error() << "Unable to read file '" << file_name << "'";
597 return true;
598 }
599
600 debug() << "Fichier ouvert : " << file_name.localstr();
601
602 VtkFile vtk_file(&ifile);
603 const char* buf = 0;
604
605 // Lecture de la description
606 // Lecture title.
607 String title = vtk_file.getNextLine();
608
609 info() << "Titre du fichier VTK : " << title.localstr();
610
611 // Lecture format.
612 String format = vtk_file.getNextLine();
613
614 debug() << "Format du fichier VTK : " << format.localstr();
615
616 if (VtkFile::isEqualString(format, "BINARY")) {
617 vtk_file.setIsBinaryFile(true);
618 }
619
620 eMeshType mesh_type = VTK_MT_Unknown;
621
622 // Lecture du type de maillage
623 // TODO: en parallèle, avec use_internal_partition vrai, seul le processeur 0
624 // lit les données. Dans ce cas, inutile que les autres ouvrent le fichier.
625 {
626 buf = vtk_file.getNextLine();
627
628 std::istringstream mesh_type_line(buf);
629 std::string dataset_str;
630 std::string mesh_type_str;
631
632 mesh_type_line >> ws >> dataset_str >> ws >> mesh_type_str;
633
634 vtk_file.checkString(dataset_str, "DATASET");
635
636 if (VtkFile::isEqualString(mesh_type_str, "STRUCTURED_GRID")) {
637 mesh_type = VTK_MT_StructuredGrid;
638 }
639
640 if (VtkFile::isEqualString(mesh_type_str, "UNSTRUCTURED_GRID")) {
641 mesh_type = VTK_MT_UnstructuredGrid;
642 }
643
644 if (mesh_type == VTK_MT_Unknown) {
645 error() << "Support exists only for 'STRUCTURED_GRID' and 'UNSTRUCTURED_GRID' formats (format=" << mesh_type_str << "')";
646 return true;
647 }
648 }
649 debug() << "Lecture en-tête OK";
650
651 bool ret = true;
652 switch (mesh_type) {
653 case VTK_MT_StructuredGrid:
654 ret = _readStructuredGrid(mesh, vtk_file, use_internal_partition);
655 break;
656
657 case VTK_MT_UnstructuredGrid:
658 ret = _readUnstructuredGrid(mesh, vtk_file, use_internal_partition);
659 debug() << "Lecture _readUnstructuredGrid OK";
660 if (!ret) {
661 // Tente de lire le fichier des faces s'il existe
662 _readFacesMesh(mesh, file_name + "faces.vtk", dir_name, use_internal_partition);
663 debug() << "Lecture _readFacesMesh OK";
664 }
665 break;
666
667 case VTK_MT_Unknown:
668 break;
669 }
670 /*while ( (buf=vtk_file.getNextLine()) != 0 ){
671 info() << " STR " << buf;
672 }*/
673
674 ifile.close();
675 return ret;
676}
677
678/*---------------------------------------------------------------------------*/
679/*---------------------------------------------------------------------------*/
680
690_readStructuredGrid(IPrimaryMesh* mesh, VtkFile& vtk_file, bool use_internal_partition)
691{
692 // Lecture du nombre de points: DIMENSIONS nx ny nz
693 const char* buf = nullptr;
694 Integer nb_node_x = 0;
695 Integer nb_node_y = 0;
696 Integer nb_node_z = 0;
697 {
698 buf = vtk_file.getNextLine();
699 std::istringstream iline(buf);
700 std::string dimension_str;
701 iline >> ws >> dimension_str >> ws >> nb_node_x >> ws >> nb_node_y >> ws >> nb_node_z;
702
703 if (!iline) {
704 error() << "Syntax error while reading grid dimensions";
705 return true;
706 }
707
708 vtk_file.checkString(dimension_str, "DIMENSIONS");
709 if (nb_node_x <= 1 || nb_node_y <= 1 || nb_node_z <= 1) {
710 error() << "Invalid dimensions: x=" << nb_node_x << " y=" << nb_node_y << " z=" << nb_node_z;
711 return true;
712 }
713 }
714 info() << " Infos: " << nb_node_x << " " << nb_node_y << " " << nb_node_z;
715 Integer nb_node = nb_node_x * nb_node_y * nb_node_z;
716
717 // Lecture du nombre de points: POINTS nb float
718 std::string float_str;
719 {
720 buf = vtk_file.getNextLine();
721 std::istringstream iline(buf);
722 std::string points_str;
723 Integer nb_node_read = 0;
724 iline >> ws >> points_str >> ws >> nb_node_read >> ws >> float_str;
725 if (!iline) {
726 error() << "Syntax error while reading grid dimensions";
727 return true;
728 }
729 vtk_file.checkString(points_str, "POINTS");
730 if (nb_node_read != nb_node) {
731 error() << "Number of invalid nodes: expected=" << nb_node << " found=" << nb_node_read;
732 return true;
733 }
734 }
735
736 Int32 rank = mesh->parallelMng()->commRank();
737
738 Integer nb_cell_x = nb_node_x - 1;
739 Integer nb_cell_y = nb_node_y - 1;
740 Integer nb_cell_z = nb_node_z - 1;
741
742 if (use_internal_partition && rank != 0) {
743 nb_node_x = 0;
744 nb_node_y = 0;
745 nb_node_z = 0;
746 nb_cell_x = 0;
747 nb_cell_y = 0;
748 nb_cell_z = 0;
749 }
750
751 const Integer nb_node_yz = nb_node_y * nb_node_z;
752 const Integer nb_node_xy = nb_node_x * nb_node_y;
753
754 Integer nb_cell = nb_cell_x * nb_cell_y * nb_cell_z;
755 UniqueArray<Int32> cells_local_id(nb_cell);
756
757 // Creation du maillage
758 {
759 UniqueArray<Integer> nodes_unique_id(nb_node);
760
761 info() << " NODE YZ = " << nb_node_yz;
762 // Création des noeuds
763 //Integer nb_node_local_id = 0;
764 {
765 Integer node_local_id = 0;
766 for (Integer x = 0; x < nb_node_x; ++x) {
767 for (Integer z = 0; z < nb_node_z; ++z) {
768 for (Integer y = 0; y < nb_node_y; ++y) {
769
770 Integer node_unique_id = y + (z)*nb_node_y + x * nb_node_y * nb_node_z;
771
772 nodes_unique_id[node_local_id] = node_unique_id;
773
774 ++node_local_id;
775 }
776 }
777 }
778 //nb_node_local_id = node_local_id;
779 //warning() << " NB NODE LOCAL ID=" << node_local_id;
780 }
781
782 // Création des mailles
783
784 // Infos pour la création des mailles
785 // par maille: 1 pour son unique id,
786 // 1 pour son type,
787 // 8 pour chaque noeud
788 UniqueArray<Int64> cells_infos(nb_cell * 10);
789
790 {
791 Integer cell_local_id = 0;
792 Integer cells_infos_index = 0;
793
794 // Normalement ne doit pas arriver car les valeurs de nb_node_x et
795 // nb_node_y sont testées lors de la lecture.
796 if (nb_node_xy == 0)
797 ARCANE_FATAL("Null value for nb_node_xy");
798
799 //Integer index = 0;
800 for (Integer z = 0; z < nb_cell_z; ++z) {
801 for (Integer y = 0; y < nb_cell_y; ++y) {
802 for (Integer x = 0; x < nb_cell_x; ++x) {
803 Integer current_cell_nb_node = 8;
804
805 //Integer cell_unique_id = y + (z)*nb_cell_y + x*nb_cell_y*nb_cell_z;
806 Int64 cell_unique_id = x + y * nb_cell_x + z * nb_cell_x * nb_cell_y;
807
808 cells_infos[cells_infos_index] = IT_Hexaedron8;
809 ++cells_infos_index;
810
811 cells_infos[cells_infos_index] = cell_unique_id;
812 ++cells_infos_index;
813
814 //Integer base_id = y + z*nb_node_y + x*nb_node_yz;
815 Integer base_id = x + y * nb_node_x + z * nb_node_xy;
816 cells_infos[cells_infos_index + 0] = nodes_unique_id[base_id];
817 cells_infos[cells_infos_index + 1] = nodes_unique_id[base_id + 1];
818 cells_infos[cells_infos_index + 2] = nodes_unique_id[base_id + nb_node_x + 1];
819 cells_infos[cells_infos_index + 3] = nodes_unique_id[base_id + nb_node_x + 0];
820 cells_infos[cells_infos_index + 4] = nodes_unique_id[base_id + nb_node_xy];
821 cells_infos[cells_infos_index + 5] = nodes_unique_id[base_id + nb_node_xy + 1];
822 cells_infos[cells_infos_index + 6] = nodes_unique_id[base_id + nb_node_xy + nb_node_x + 1];
823 cells_infos[cells_infos_index + 7] = nodes_unique_id[base_id + nb_node_xy + nb_node_x + 0];
824 cells_infos_index += current_cell_nb_node;
825 cells_local_id[cell_local_id] = cell_local_id;
826 ++cell_local_id;
827 }
828 }
829 }
830 }
831
832 mesh->setDimension(3);
833 mesh->allocateCells(nb_cell, cells_infos, false);
834 mesh->endAllocate();
835
836 // Positionne les coordonnées
837 {
838 UniqueArray<Real3> coords(nb_node);
839
840 if (vtk_file.isEqualString(float_str, "int")) {
841 for (Integer z = 0; z < nb_node_z; ++z) {
842 for (Integer y = 0; y < nb_node_y; ++y) {
843 for (Integer x = 0; x < nb_node_x; ++x) {
844 Real nx = vtk_file.getInt();
845 Real ny = vtk_file.getInt();
846 Real nz = vtk_file.getInt();
847 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
848 coords[node_unique_id] = Real3(nx, ny, nz);
849 }
850 }
851 }
852 }
853 else if (vtk_file.isEqualString(float_str, "float")) {
854 for (Integer z = 0; z < nb_node_z; ++z) {
855 for (Integer y = 0; y < nb_node_y; ++y) {
856 for (Integer x = 0; x < nb_node_x; ++x) {
857 Real nx = vtk_file.getFloat();
858 Real ny = vtk_file.getFloat();
859 Real nz = vtk_file.getFloat();
860 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
861 coords[node_unique_id] = Real3(nx, ny, nz);
862 }
863 }
864 }
865 }
866 else if (vtk_file.isEqualString(float_str, "double")) {
867 for (Integer z = 0; z < nb_node_z; ++z) {
868 for (Integer y = 0; y < nb_node_y; ++y) {
869 for (Integer x = 0; x < nb_node_x; ++x) {
870 Real nx = vtk_file.getDouble();
871 Real ny = vtk_file.getDouble();
872 Real nz = vtk_file.getDouble();
873 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
874 coords[node_unique_id] = Real3(nx, ny, nz);
875 }
876 }
877 }
878 }
879 else {
880 throw IOException("_readStructuredGrid", "Invalid type name");
881 }
882
883 VariableNodeReal3& nodes_coord_var(mesh->nodesCoordinates());
884 ENUMERATE_NODE (inode, mesh->allNodes()) {
885 Node node = *inode;
886 nodes_coord_var[inode] = coords[node.uniqueId().asInt32()];
887 }
888 }
889 }
890
891 // Créé les groupes de faces des côtés du parallélépipède
892 {
893 Int32UniqueArray xmin_surface_lid;
894 Int32UniqueArray xmax_surface_lid;
895 Int32UniqueArray ymin_surface_lid;
896 Int32UniqueArray ymax_surface_lid;
897 Int32UniqueArray zmin_surface_lid;
898 Int32UniqueArray zmax_surface_lid;
899
900 ENUMERATE_FACE (iface, mesh->allFaces()) {
901 Face face = *iface;
902 Integer face_local_id = face.localId();
903 bool is_xmin = true;
904 bool is_xmax = true;
905 bool is_ymin = true;
906 bool is_ymax = true;
907 bool is_zmin = true;
908 bool is_zmax = true;
909 for (Node node : face.nodes() ) {
910 Int64 node_unique_id = node.uniqueId().asInt64();
911 Int64 node_z = node_unique_id / nb_node_xy;
912 Int64 node_y = (node_unique_id - node_z * nb_node_xy) / nb_node_x;
913 Int64 node_x = node_unique_id - node_z * nb_node_xy - node_y * nb_node_x;
914 if (node_x != 0)
915 is_xmin = false;
916 if (node_x != (nb_node_x - 1))
917 is_xmax = false;
918 if (node_y != 0)
919 is_ymin = false;
920 if (node_y != (nb_node_y - 1))
921 is_ymax = false;
922 if (node_z != 0)
923 is_zmin = false;
924 if (node_z != (nb_node_z - 1))
925 is_zmax = false;
926 }
927 if (is_xmin)
928 xmin_surface_lid.add(face_local_id);
929 if (is_xmax)
930 xmax_surface_lid.add(face_local_id);
931 if (is_ymin)
932 ymin_surface_lid.add(face_local_id);
933 if (is_ymax)
934 ymax_surface_lid.add(face_local_id);
935 if (is_zmin)
936 zmin_surface_lid.add(face_local_id);
937 if (is_zmax)
938 zmax_surface_lid.add(face_local_id);
939 }
940 _createFaceGroup(mesh, "XMIN", xmin_surface_lid);
941 _createFaceGroup(mesh, "XMAX", xmax_surface_lid);
942 _createFaceGroup(mesh, "YMIN", ymin_surface_lid);
943 _createFaceGroup(mesh, "YMAX", ymax_surface_lid);
944 _createFaceGroup(mesh, "ZMIN", zmin_surface_lid);
945 _createFaceGroup(mesh, "ZMAX", zmax_surface_lid);
946 }
947
948 _readMetadata(mesh, vtk_file);
949
950 // Maintenant, regarde s'il existe des données associées aux fichier
951 bool r = _readData(mesh, vtk_file, use_internal_partition, IK_Cell, cells_local_id, nb_node);
952 if (r)
953 return r;
954
955 return false;
956}
957
958/*---------------------------------------------------------------------------*/
959/*---------------------------------------------------------------------------*/
960
970{
971 ARCANE_UNUSED(mesh);
972
973 const char* func_name = "VtkMeshIOService::_readNodesUnstructuredGrid()";
974 const char* buf = vtk_file.getNextLine();
975 std::istringstream iline(buf);
976 std::string points_str;
977 std::string data_type_str;
978 Integer nb_node = 0;
979
980 iline >> ws >> points_str >> ws >> nb_node >> ws >> data_type_str;
981
982 if (!iline)
983 throw IOException(func_name, "Syntax error while reading number of nodes");
984
985 vtk_file.checkString(points_str, "POINTS");
986
987 if (nb_node < 0)
988 throw IOException(A_FUNCINFO, String::format("Invalid number of nodes: n={0}", nb_node));
989
990 info() << "VTK file : number of nodes = " << nb_node;
991
992 // Lecture les coordonnées
993 node_coords.resize(nb_node);
994 {
995 if (vtk_file.isEqualString(data_type_str, "int")) {
996 for (Integer i = 0; i < nb_node; ++i) {
997 Real nx = vtk_file.getInt();
998 Real ny = vtk_file.getInt();
999 Real nz = vtk_file.getInt();
1000 node_coords[i] = Real3(nx, ny, nz);
1001 }
1002 }
1003 else if (vtk_file.isEqualString(data_type_str, "float")) {
1004 for (Integer i = 0; i < nb_node; ++i) {
1005 Real nx = vtk_file.getFloat();
1006 Real ny = vtk_file.getFloat();
1007 Real nz = vtk_file.getFloat();
1008 node_coords[i] = Real3(nx, ny, nz);
1009 }
1010 }
1011 else if (vtk_file.isEqualString(data_type_str, "double")) {
1012 for (Integer i = 0; i < nb_node; ++i) {
1013 Real nx = vtk_file.getDouble();
1014 Real ny = vtk_file.getDouble();
1015 Real nz = vtk_file.getDouble();
1016 node_coords[i] = Real3(nx, ny, nz);
1017 }
1018 }
1019 else {
1020 throw IOException(func_name, "Invalid type name");
1021 }
1022 }
1023 _readMetadata(mesh, vtk_file);
1024}
1025
1026/*---------------------------------------------------------------------------*/
1027/*---------------------------------------------------------------------------*/
1028
1042 Array<Integer>& cells_nb_node,
1043 Array<ItemTypeId>& cells_type,
1044 Array<Int64>& cells_connectivity)
1045{
1046 ARCANE_UNUSED(mesh);
1047
1048 const char* func_name = "VtkMeshIOService::_readCellsUnstructuredGrid()";
1049 const char* buf = vtk_file.getNextLine();
1050
1051 //String buftest = vtk_file.getCurrentLine(); // DEBUG
1052 //pinfo() << "Ligne lu : " << buftest.localstr();
1053
1054 std::istringstream iline(buf);
1055 std::string cells_str;
1056 Integer nb_cell = 0;
1057 Integer nb_cell_node = 0;
1058
1059 iline >> ws >> cells_str >> ws >> nb_cell >> ws >> nb_cell_node;
1060
1061 if (!iline)
1062 throw IOException(func_name, "Syntax error while reading cells");
1063
1064 vtk_file.checkString(cells_str, "CELLS");
1065
1066 info() << "VTK file : number of cells = " << nb_cell;
1067
1068 if (nb_cell < 0 || nb_cell_node < 0) {
1069 throw IOException(A_FUNCINFO,
1070 String::format("Invalid dimensions: nb_cell={0} nb_cell_node={1}",
1071 nb_cell, nb_cell_node));
1072 }
1073
1074 cells_nb_node.resize(nb_cell);
1075 cells_type.resize(nb_cell);
1076 cells_connectivity.resize(nb_cell_node);
1077
1078 {
1079 Integer connectivity_index = 0;
1080 for (Integer i = 0; i < nb_cell; ++i) {
1081 Integer n = vtk_file.getInt();
1082 cells_nb_node[i] = n;
1083 for (Integer j = 0; j < n; ++j) {
1084 Integer id = vtk_file.getInt();
1085 cells_connectivity[connectivity_index] = id;
1086 ++connectivity_index;
1087 }
1088 }
1089 }
1090
1091 _readMetadata(mesh, vtk_file);
1092
1093 // Lecture du type des mailles
1094 {
1095 buf = vtk_file.getNextLine();
1096 std::istringstream iline(buf);
1097 std::string cell_types_str;
1098 Integer nb_cell_type;
1099 iline >> ws >> cell_types_str >> ws >> nb_cell_type;
1100
1101 if (!iline) {
1102 throw IOException(func_name, "Syntax error while reading cell types");
1103 }
1104
1105 vtk_file.checkString(cell_types_str, "CELL_TYPES");
1106 if (nb_cell_type != nb_cell) {
1107 ARCANE_THROW(IOException, "Inconsistency in number of CELL_TYPES: v={0} nb_cell={1}",
1108 nb_cell_type, nb_cell);
1109 }
1110 }
1111
1112 for (Integer i = 0; i < nb_cell; ++i) {
1113 Integer vtk_ct = vtk_file.getInt();
1114 Int16 it = vtkToArcaneCellType(vtk_ct, cells_nb_node[i]);
1115 cells_type[i] = ItemTypeId{ it };
1116 }
1117 _readMetadata(mesh, vtk_file);
1118}
1119
1120/*---------------------------------------------------------------------------*/
1121/*---------------------------------------------------------------------------*/
1122
1131_readMetadata(IMesh* mesh, VtkFile& vtk_file)
1132{
1133 ARCANE_UNUSED(mesh);
1134
1135 //const char* func_name = "VtkMeshIOService::_readMetadata()";
1136
1137 if (vtk_file.isEof())
1138 return false;
1139 String meta = vtk_file.getNextLine();
1140
1141 // METADATA ?
1142 if (!vtk_file.isEqualString(meta, "METADATA")) {
1143 // S'il n'y a pas de METADATA, on demande à ce que la ligne lue soit relue la prochaine fois.
1144 vtk_file.reReadSameLine();
1145 return false;
1146 }
1147
1148 // Tant qu'il n'y a pas de ligne vide, on lit.
1149 while (!vtk_file.isEmptyNextLine() && !vtk_file.isEof()) {
1150 }
1151 return false;
1152
1153 // // Si l'on a besoin de faire quelque chose avec les METADATA un jour, voilà un code non testé.
1154 // std::string trash;
1155 // Real trash_real;
1156 // const char* buf = vtk_file.getNextLine();
1157
1158 // // INFORMATION ou COMPONENT_NAMES
1159 // std::istringstream iline(buf);
1160
1161 // String name_str;
1162 // String data_type_str;
1163 // Integer nb_info = 0;
1164 // Integer size_vector = 0;
1165
1166 // iline >> ws >> name_str;
1167
1168 // if(vtk_file.isEqualString(name_str, "INFORMATION"))
1169 // {
1170 // iline >> ws >> nb_info;
1171 // for( Integer i = 0; i < nb_info; i++)
1172 // {
1173 // buf = vtk_file.getNextLine();
1174 // std::istringstream iline(buf);
1175
1176 // // NAME [key name] LOCATION [key location (e.g. class name)]
1177 // iline >> ws >> trash >> ws >> data_type_str >> ws >> trash >> ws >> name_str;
1178
1179 // buf = vtk_file.getNextLine();
1180 // std::istringstream iline(buf);
1181
1182 // if(vtk_file.isEqualString(data_type_str, "StringVector"))
1183 // {
1184 // iline >> ws >> trash >> ws >> size_vector;
1185 // for(Integer j = 0; j < size_vector; j++)
1186 // {
1187 // trash = vtk_file.getNextLine();
1188 // }
1189 // }
1190
1191 // // else if(vtk_file.isEqualString(data_type_str, "Double")
1192 // // || vtk_file.isEqualString(data_type_str, "IdType")
1193 // // || vtk_file.isEqualString(data_type_str, "Integer")
1194 // // || vtk_file.isEqualString(data_type_str, "UnsignedLong")
1195 // // || vtk_file.isEqualString(data_type_str, "String"))
1196 // // {
1197 // // iline >> ws >> trash >> ws >> trash_real;
1198 // // }
1199
1200 // // else if(vtk_file.isEqualString(data_type_str, "DoubleVector")
1201 // // || vtk_file.isEqualString(data_type_str, "IntegerVector")
1202 // // || vtk_file.isEqualString(data_type_str, "StringVector"))
1203 // // {
1204 // // iline >> ws >> trash >> ws >> size_vector;
1205 // // for(Integer j = 0; j < size_vector; j++)
1206 // // {
1207 // // iline >> ws >> trash_real;
1208 // // }
1209 // // }
1210 // }
1211 // }
1212
1213 // else if(vtk_file.isEqualString(name_str, "COMPONENT_NAMES"))
1214 // {
1215 // while(!vtk_file.isEmptyNextLine())
1216 // {
1217 // trash = vtk_file.getCurrentLine();
1218 // }
1219 // }
1220
1221 // else
1222 // {
1223 // throw IOException(func_name,"Syntax error after METADATA tag");
1224 // }
1225}
1226
1227/*---------------------------------------------------------------------------*/
1228/*---------------------------------------------------------------------------*/
1229
1239_readUnstructuredGrid(IPrimaryMesh* mesh, VtkFile& vtk_file, bool use_internal_partition)
1240{
1241 Integer nb_node = 0;
1242 Integer nb_cell = 0;
1243 Integer nb_cell_node = 0;
1244 Int32 sid = mesh->parallelMng()->commRank();
1245 UniqueArray<Real3> node_coords;
1246 UniqueArray<Int32> cells_local_id;
1247
1248 // Si on utilise le partitionneur interne, seul le sous-domaine lit le maillage
1249 bool need_read = true;
1250
1251 if (use_internal_partition)
1252 need_read = (sid == 0);
1253
1254 std::array<Int64, 4> nb_cell_by_dimension = {};
1255 Int32 mesh_dimension = -1;
1256 ItemTypeMng* itm = mesh->itemTypeMng();
1257 UnstructuredMeshAllocateBuildInfo mesh_build_info(mesh);
1258
1259 if (need_read) {
1260 // Lecture première partie du fichier (après header).
1261 _readNodesUnstructuredGrid(mesh, vtk_file, node_coords);
1262 debug() << "Lecture _readNodesUnstructuredGrid OK";
1263 nb_node = node_coords.size();
1264
1265 // Lecture des infos des mailles
1266 // Lecture de la connectivité
1267 UniqueArray<Integer> cells_nb_node;
1268 UniqueArray<Int64> cells_connectivity;
1269 UniqueArray<ItemTypeId> cells_type;
1270 _readCellsUnstructuredGrid(mesh, vtk_file, cells_nb_node, cells_type, cells_connectivity);
1271 debug() << "Lecture _readCellsUnstructuredGrid OK";
1272
1273 nb_cell = cells_nb_node.size();
1274 nb_cell_node = cells_connectivity.size();
1275 cells_local_id.resize(nb_cell);
1276
1277 // Création des mailles
1278 mesh_build_info.preAllocate(nb_cell, nb_cell_node);
1279
1280 {
1281 Int32 connectivity_index = 0;
1282 for (Integer i = 0; i < nb_cell; ++i) {
1283 Int32 current_cell_nb_node = cells_nb_node[i];
1284 Int64 cell_unique_id = i;
1285
1286 cells_local_id[i] = i;
1287
1288 Int16 cell_dim = itm->typeFromId(cells_type[i])->dimension();
1289 if (cell_dim >= 0 && cell_dim <= 3)
1290 ++nb_cell_by_dimension[cell_dim];
1291
1292 auto cell_nodes = cells_connectivity.subView(connectivity_index,current_cell_nb_node);
1293 mesh_build_info.addCell(cells_type[i],cell_unique_id, cell_nodes);
1294 connectivity_index += current_cell_nb_node;
1295 }
1296 }
1297 // Vérifie qu'on n'a pas de mailles de différentes dimensions
1298 Int32 nb_different_dim = 0;
1299 for (Int32 i = 0; i < 4; ++i)
1300 if (nb_cell_by_dimension[i] != 0) {
1301 ++nb_different_dim;
1302 mesh_dimension = i;
1303 }
1304 if (nb_different_dim > 1)
1305 ARCANE_FATAL("The mesh contains cells of different dimension. nb0={0} nb1={1} nb2={2} nb3={3}",
1306 nb_cell_by_dimension[0], nb_cell_by_dimension[1], nb_cell_by_dimension[2], nb_cell_by_dimension[3]);
1307 }
1308
1309 // Positionne la dimension du maillage.
1310 {
1311 Integer wanted_dimension = mesh_dimension;
1312 wanted_dimension = mesh->parallelMng()->reduce(Parallel::ReduceMax, wanted_dimension);
1313 mesh->setDimension(wanted_dimension);
1314 mesh_build_info.setMeshDimension(wanted_dimension);
1315 }
1316
1317 mesh_build_info.allocateMesh();
1318
1319 // Positionne les coordonnées
1320 {
1321 VariableNodeReal3& nodes_coord_var(mesh->nodesCoordinates());
1322 ENUMERATE_NODE (inode, mesh->allNodes()) {
1323 Node node = *inode;
1324 nodes_coord_var[inode] = node_coords[node.uniqueId().asInt32()];
1325 }
1326 }
1327
1328 // Maintenant, regarde s'il existe des données associées aux fichier
1329 bool r = _readData(mesh, vtk_file, use_internal_partition, IK_Cell, cells_local_id, nb_node);
1330 debug() << "Lecture _readData OK";
1331
1332 return r;
1333}
1334
1335/*---------------------------------------------------------------------------*/
1336/*---------------------------------------------------------------------------*/
1337
1347_readFacesMesh(IMesh* mesh, const String& file_name, const String& dir_name,
1348 bool use_internal_partition)
1349{
1350 ARCANE_UNUSED(dir_name);
1351
1352 std::ifstream ifile(file_name.localstr(), std::ifstream::binary);
1353 if (!ifile) {
1354 info() << "No face descriptor file found '" << file_name << "'";
1355 return;
1356 }
1357
1358 VtkFile vtk_file(&ifile);
1359 const char* buf = 0;
1360
1361 // Lecture de la description
1362 String title = vtk_file.getNextLine();
1363 info() << "Reading VTK file '" << file_name << "'";
1364 info() << "Title of VTK file: " << title;
1365
1366 String format = vtk_file.getNextLine();
1367 if (VtkFile::isEqualString(format, "BINARY")) {
1368 vtk_file.setIsBinaryFile(true);
1369 }
1370
1371 eMeshType mesh_type = VTK_MT_Unknown;
1372 // Lecture du type de maillage
1373 // TODO: en parallèle, avec use_internal_partition vrai, seul le processeur 0
1374 // lit les données. Dans ce cas, inutile que les autres ouvre le fichier.
1375 {
1376 buf = vtk_file.getNextLine();
1377 std::istringstream mesh_type_line(buf);
1378 std::string dataset_str;
1379 std::string mesh_type_str;
1380 mesh_type_line >> ws >> dataset_str >> ws >> mesh_type_str;
1381 vtk_file.checkString(dataset_str, "DATASET");
1382
1383 if (VtkFile::isEqualString(mesh_type_str, "UNSTRUCTURED_GRID")) {
1384 mesh_type = VTK_MT_UnstructuredGrid;
1385 }
1386
1387 if (mesh_type == VTK_MT_Unknown) {
1388 error() << "Face descriptor file type must be 'UNSTRUCTURED_GRID' (format=" << mesh_type_str << "')";
1389 return;
1390 }
1391 }
1392
1393 {
1394 IParallelMng* pm = mesh->parallelMng();
1395 Integer nb_face = 0;
1396 Int32 sid = pm->commRank();
1397
1398 UniqueArray<Int32> faces_local_id;
1399
1400 // Si on utilise le partitionneur interne, seul le sous-domaine lit le maillage
1401 bool need_read = true;
1402 if (use_internal_partition)
1403 need_read = (sid == 0);
1404
1405 if (need_read) {
1406 {
1407 // Lit des noeuds, mais ne conserve pas leur coordonnées car cela n'est
1408 // pas nécessaire.
1409 UniqueArray<Real3> node_coords;
1410 _readNodesUnstructuredGrid(mesh, vtk_file, node_coords);
1411 //nb_node = node_coords.size();
1412 }
1413
1414 // Lecture des infos des faces
1415 // Lecture de la connectivité
1416 UniqueArray<Integer> faces_nb_node;
1417 UniqueArray<Int64> faces_connectivity;
1418 UniqueArray<ItemTypeId> faces_type;
1419 _readCellsUnstructuredGrid(mesh, vtk_file, faces_nb_node, faces_type, faces_connectivity);
1420 nb_face = faces_nb_node.size();
1421 //nb_face_node = faces_connectivity.size();
1422
1423 // Il faut à partir de la connectivité retrouver les localId() des faces
1424 faces_local_id.resize(nb_face);
1425 {
1426 IMeshUtilities* mu = mesh->utilities();
1427 mu->getFacesLocalIdFromConnectivity(faces_type, faces_connectivity, faces_local_id);
1428 }
1429 }
1430
1431 // Maintenant, regarde s'il existe des données associées aux fichiers
1432 _readData(mesh, vtk_file, use_internal_partition, IK_Face, faces_local_id, 0);
1433 }
1434}
1435
1436/*---------------------------------------------------------------------------*/
1437/*---------------------------------------------------------------------------*/
1438
1451_readData(IMesh* mesh, VtkFile& vtk_file, bool use_internal_partition,
1452 eItemKind cell_kind, Int32ConstArrayView local_id, Integer nb_node)
1453{
1454 // Seul le sous-domain maitre lit les valeurs. Par contre, les autres
1455 // sous-domaines doivent connaitre la liste des variables et groupes créées.
1456 // Si une donnée porte le nom 'GROUP_*', on considère qu'il s'agit d'un
1457 // groupe
1458
1459 IParallelMng* pm = mesh->parallelMng();
1460 IVariableMng* variable_mng = mesh->variableMng();
1461
1462 Int32 sid = pm->commRank();
1463
1464 // Si pas de données, retourne immédiatement.
1465 {
1466 Byte has_data = 1;
1467 if ((sid == 0) && vtk_file.isEof())
1468 has_data = 0;
1469
1470 ByteArrayView bb(1, &has_data);
1471 pm->broadcast(bb, 0);
1472 // Pas de data.
1473 if (!has_data)
1474 return false;
1475 }
1476
1477 OStringStream created_infos_str;
1478 created_infos_str() << "<?xml version='1.0' ?>\n";
1479 created_infos_str() << "<infos>";
1480
1481 Integer nb_cell_kind = mesh->nbItem(cell_kind);
1482 const char* buf = 0;
1483
1484 if (sid == 0) {
1485 bool reading_node = false;
1486 bool reading_cell = false;
1487 while (((buf = vtk_file.getNextLine()) != 0) && !vtk_file.isEof()) {
1488 debug() << "Read line";
1489 std::istringstream iline(buf);
1490 std::string data_str;
1491 iline >> data_str;
1492
1493 // Si l'on a un bloc "CELL_DATA".
1494 if (VtkFile::isEqualString(data_str, "CELL_DATA")) {
1495 Integer nb_item = 0;
1496 iline >> ws >> nb_item;
1497 reading_node = false;
1498 reading_cell = true;
1499 if (nb_item != nb_cell_kind)
1500 error() << "Size expected = " << nb_cell_kind << " found = " << nb_item;
1501 }
1502
1503 // Si l'on a un bloc "POINT_DATA".
1504 else if (VtkFile::isEqualString(data_str, "POINT_DATA")) {
1505 Integer nb_item = 0;
1506 iline >> ws >> nb_item;
1507 reading_node = true;
1508 reading_cell = false;
1509 if (nb_item != nb_node)
1510 error() << "Size expected = " << nb_node << " found = " << nb_item;
1511 }
1512
1513 // Si l'on a un bloc "FIELD".
1514 else if (VtkFile::isEqualString(data_str, "FIELD")) {
1515 std::string name_str;
1516 int nb_fields;
1517
1518 iline >> ws >> name_str >> ws >> nb_fields;
1519
1520 Integer nb_item = 0;
1521 std::string type_str;
1522 std::string s_name_str;
1523 int nb_component = 1;
1524 bool is_group = false;
1525
1526 for (Integer i = 0; i < nb_fields; i++) {
1527 buf = vtk_file.getNextLine();
1528 std::istringstream iline(buf);
1529 iline >> ws >> s_name_str >> ws >> nb_component >> ws >> nb_item >> ws >> type_str;
1530
1531 if (nb_item != nb_cell_kind && reading_cell && !reading_node)
1532 error() << "Size expected = " << nb_cell_kind << " found = " << nb_item;
1533
1534 if (nb_item != nb_node && !reading_cell && reading_node)
1535 error() << "Size expected = " << nb_node << " found = " << nb_item;
1536
1537 String name_str = s_name_str;
1538 String cstr = name_str.substring(0, 6);
1539
1540 if (cstr == "GROUP_") {
1541 is_group = true;
1542 String new_name = name_str.substring(6);
1543 debug() << "** ** ** GROUP ! name=" << new_name;
1544 name_str = new_name;
1545 }
1546
1547 if (is_group) {
1548 if (!VtkFile::isEqualString(type_str, "int")) {
1549 error() << "Group type must be 'int', found=" << type_str;
1550 return true;
1551 }
1552
1553 if (reading_node) {
1554 created_infos_str() << "<node-group name='" << name_str << "'/>";
1555 _readNodeGroup(mesh, vtk_file, name_str, nb_node);
1556 }
1557
1558 if (reading_cell) {
1559 created_infos_str() << "<cell-group name='" << name_str << "'/>";
1560 _readItemGroup(mesh, vtk_file, name_str, nb_cell_kind, cell_kind, local_id);
1561 }
1562 }
1563
1564 // TODO : Voir un exemple si possible.
1565 else {
1566 if (!VtkFile::isEqualString(type_str, "float") && !VtkFile::isEqualString(type_str, "double")) {
1567 error() << "Expecting 'float' or 'double' data type, found=" << type_str;
1568 return true;
1569 }
1570
1571 if (reading_node) {
1572 fatal() << "Unable to read POINT_DATA: feature not implemented";
1573 }
1574 if (reading_cell) {
1575 created_infos_str() << "<cell-variable name='" << name_str << "'/>";
1576
1577 if (cell_kind != IK_Cell)
1578 throw IOException("Unable to read face variables: feature not supported");
1579
1580 _readCellVariable(mesh, vtk_file, name_str, nb_cell_kind);
1581 }
1582 }
1583 }
1584 }
1585
1586 else {
1587 // Lecture des valeurs (bloc "CELL_DATA" ou "POINT_DATA")
1588 if (reading_node || reading_cell) {
1589 std::string type_str;
1590 std::string s_name_str;
1591 //String name_str;
1592 bool is_group = false;
1593 int nb_component = 1;
1594
1595 iline >> ws >> s_name_str >> ws >> type_str >> ws >> nb_component;
1596 debug() << "** ** ** READNAME: name=" << s_name_str << " type=" << type_str;
1597
1598 String name_str = s_name_str;
1599 String cstr = name_str.substring(0, 6);
1600
1601 if (cstr == "GROUP_") {
1602 is_group = true;
1603 String new_name = name_str.substring(6);
1604 info() << "** ** ** GROUP ! name=" << new_name;
1605 name_str = new_name;
1606 }
1607
1608 if (!VtkFile::isEqualString(data_str, "SCALARS")) {
1609 error() << "Expecting 'SCALARS' data type, found=" << data_str;
1610 return true;
1611 }
1612
1613 if (is_group) {
1614 if (!VtkFile::isEqualString(type_str, "int")) {
1615 error() << "Group type must be 'int', found=" << type_str;
1616 return true;
1617 }
1618
1619 // Pour lire LOOKUP_TABLE
1620 buf = vtk_file.getNextLine();
1621
1622 if (reading_node) {
1623 created_infos_str() << "<node-group name='" << name_str << "'/>";
1624 _readNodeGroup(mesh, vtk_file, name_str, nb_node);
1625 }
1626
1627 if (reading_cell) {
1628 created_infos_str() << "<cell-group name='" << name_str << "'/>";
1629 _readItemGroup(mesh, vtk_file, name_str, nb_cell_kind, cell_kind, local_id);
1630 }
1631 }
1632 else {
1633 if (!VtkFile::isEqualString(type_str, "float") && !VtkFile::isEqualString(type_str, "double")) {
1634 error() << "Expecting 'float' or 'double' data type, found=" << type_str;
1635 return true;
1636 }
1637
1638 // Pour lire LOOKUP_TABLE
1639 /*buf = */ vtk_file.getNextLine();
1640 if (reading_node) {
1641 fatal() << "Unable to read POINT_DATA: feature not implemented";
1642 }
1643 if (reading_cell) {
1644 created_infos_str() << "<cell-variable name='" << name_str << "'/>";
1645
1646 if (cell_kind != IK_Cell)
1647 throw IOException("Unable to read face variables: feature not supported");
1648
1649 _readCellVariable(mesh, vtk_file, name_str, nb_cell_kind);
1650 }
1651 }
1652 }
1653 else {
1654 error() << "Expecting value CELL_DATA or POINT_DATA, found='" << data_str << "'";
1655 return true;
1656 }
1657 }
1658 }
1659 }
1660 created_infos_str() << "</infos>";
1661 if (use_internal_partition) {
1662 ByteUniqueArray bytes;
1663 if (sid == 0) {
1664 String str = created_infos_str.str();
1665 ByteConstArrayView bv = str.utf8();
1666 Integer len = bv.size();
1667 bytes.resize(len + 1);
1668 bytes.copy(bv);
1669 }
1670
1671 pm->broadcastMemoryBuffer(bytes, 0);
1672
1673 if (sid != 0) {
1674 String str = String::fromUtf8(bytes);
1675 info() << "FOUND STR=" << bytes.size() << " " << str;
1677 XmlNode doc_node = doc->documentNode();
1678
1679 // Lecture des variables
1680 {
1681 XmlNodeList vars = doc_node.documentElement().children("cell-variable");
1682 for (XmlNode xnode : vars) {
1683 String name = xnode.attrValue("name");
1684 info() << "Building variable: " << name;
1686 variable_mng->_internalApi()->addAutoDestroyVariable(var);
1687 }
1688 }
1689
1690 // Lecture des groupes de mailles
1691 {
1692 XmlNodeList vars = doc_node.documentElement().children("cell-group");
1693 IItemFamily* cell_family = mesh->itemFamily(cell_kind);
1694 for (XmlNode xnode : vars) {
1695 String name = xnode.attrValue("name");
1696 info() << "Building group: " << name;
1697 cell_family->createGroup(name);
1698 }
1699 }
1700
1701 // Lecture des groupes de noeuds
1702 {
1703 XmlNodeList vars = doc_node.documentElement().children("node-group");
1704 IItemFamily* node_family = mesh->nodeFamily();
1705 for (XmlNode xnode : vars) {
1706 String name = xnode.attrValue("name");
1707 info() << "Create node group: " << name;
1708 node_family->createGroup(name);
1709 }
1710 }
1711 }
1712 }
1713 return false;
1714}
1715
1716/*---------------------------------------------------------------------------*/
1717/*---------------------------------------------------------------------------*/
1718
1727_createFaceGroup(IMesh* mesh, const String& name, Int32ConstArrayView faces_lid)
1728{
1729 info() << "Building face group '" << name << "'"
1730 << " size=" << faces_lid.size();
1731
1732 mesh->faceFamily()->createGroup(name, faces_lid);
1733}
1734
1735/*---------------------------------------------------------------------------*/
1736/*---------------------------------------------------------------------------*/
1737
1747_readCellVariable(IMesh* mesh, VtkFile& vtk_file, const String& var_name, Integer nb_cell)
1748{
1749 //TODO Faire la conversion uniqueId() vers localId() correcte
1750 info() << "Reading values for variable: " << var_name << " n=" << nb_cell;
1751 auto* var = new VariableCellReal(VariableBuildInfo(mesh, var_name));
1752 mesh->variableMng()->_internalApi()->addAutoDestroyVariable(var);
1753 RealArrayView values(var->asArray());
1754 for (Integer i = 0; i < nb_cell; ++i) {
1755 Real v = vtk_file.getDouble();
1756 values[i] = v;
1757 }
1758 _readMetadata(mesh, vtk_file);
1759 info() << "Variable build finished: " << vtk_file.isEof();
1760}
1761
1762/*---------------------------------------------------------------------------*/
1763/*---------------------------------------------------------------------------*/
1764
1776_readItemGroup(IMesh* mesh, VtkFile& vtk_file, const String& name, Integer nb_item,
1777 eItemKind ik, Int32ConstArrayView local_id)
1778{
1779 IItemFamily* item_family = mesh->itemFamily(ik);
1780 info() << "Reading group info for group: " << name;
1781
1782 Int32UniqueArray ids;
1783 for (Integer i = 0; i < nb_item; ++i) {
1784 Integer v = vtk_file.getInt();
1785 if (v != 0)
1786 ids.add(local_id[i]);
1787 }
1788 info() << "Building group: " << name << " nb_element=" << ids.size();
1789
1790 item_family->createGroup(name, ids);
1791
1792 _readMetadata(mesh, vtk_file);
1793}
1794
1795/*---------------------------------------------------------------------------*/
1796/*---------------------------------------------------------------------------*/
1797
1807_readNodeGroup(IMesh* mesh, VtkFile& vtk_file, const String& name, Integer nb_item)
1808{
1809 IItemFamily* item_family = mesh->itemFamily(IK_Node);
1810 info() << "Lecture infos groupes de noeuds pour le groupe: " << name;
1811
1812 Int32UniqueArray ids;
1813 for (Integer i = 0; i < nb_item; ++i) {
1814 Integer v = vtk_file.getInt();
1815 if (v != 0)
1816 ids.add(i);
1817 }
1818 info() << "Création groupe: " << name << " nb_element=" << ids.size();
1819
1820 item_family->createGroup(name, ids);
1821
1822 _readMetadata(mesh, vtk_file);
1823}
1824
1825/*---------------------------------------------------------------------------*/
1826/*---------------------------------------------------------------------------*/
1827
1828/*---------------------------------------------------------------------------*/
1829/*---------------------------------------------------------------------------*/
1830
1831class VtkLegacyMeshWriter
1832: public BasicService
1833, public IMeshWriter
1834{
1835 public:
1836
1837 explicit VtkLegacyMeshWriter(const ServiceBuildInfo& sbi)
1838 : BasicService(sbi)
1839 {}
1840
1841 public:
1842
1843 void build() override {}
1844
1845 public:
1846
1847 bool writeMeshToFile(IMesh* mesh, const String& file_name) override;
1848
1849 private:
1850
1851 void _writeMeshToFile(IMesh* mesh, const String& file_name, eItemKind cell_kind);
1852 void _saveGroups(IItemFamily* family, std::ostream& ofile);
1853};
1854
1855/*---------------------------------------------------------------------------*/
1856/*---------------------------------------------------------------------------*/
1857
1859 ServiceProperty("VtkLegacyMeshWriter", ST_SubDomain),
1861
1862/*---------------------------------------------------------------------------*/
1863/*---------------------------------------------------------------------------*/
1879writeMeshToFile(IMesh* mesh, const String& file_name)
1880{
1881 String fname = file_name;
1882 // Ajoute l'extension '.vtk' si elle n'y est pas.
1883 if (!fname.endsWith(".vtk"))
1884 fname = fname + ".vtk";
1885 _writeMeshToFile(mesh, fname, IK_Cell);
1886 _writeMeshToFile(mesh, fname + "faces.vtk", IK_Face);
1887 return false;
1888}
1889
1890/*---------------------------------------------------------------------------*/
1891/*---------------------------------------------------------------------------*/
1899_writeMeshToFile(IMesh* mesh, const String& file_name, eItemKind cell_kind)
1900{
1901 std::ofstream ofile(file_name.localstr());
1902 ofile.precision(FloatInfo<Real>::maxDigit());
1903 if (!ofile)
1904 throw IOException("VtkMeshIOService::writeMeshToFile(): Unable to open file");
1905 ofile << "# vtk DataFile Version 2.0\n";
1906 ofile << "Maillage Arcane\n";
1907 ofile << "ASCII\n";
1908 ofile << "DATASET UNSTRUCTURED_GRID\n";
1909 UniqueArray<Integer> nodes_local_id_to_current(mesh->itemFamily(IK_Node)->maxLocalId());
1910 nodes_local_id_to_current.fill(NULL_ITEM_ID);
1911
1912 Integer nb_node = mesh->nbNode();
1913 IItemFamily* cell_kind_family = mesh->itemFamily(cell_kind);
1914 Integer nb_cell_kind = cell_kind_family->nbItem();
1915
1916 // Sauve les nœuds
1917 {
1918 ofile << "POINTS " << nb_node << " double\n";
1919 VariableNodeReal3& coords(mesh->toPrimaryMesh()->nodesCoordinates());
1920 Integer node_index = 0;
1921 ENUMERATE_NODE (inode, mesh->allNodes()) {
1922 Node node = *inode;
1923 nodes_local_id_to_current[node.localId()] = node_index;
1924 Real3 xyz = coords[inode];
1925 ofile << xyz.x << ' ' << xyz.y << ' ' << xyz.z << '\n';
1926 ++node_index;
1927 }
1928 }
1929
1930 // Sauve les mailles ou faces
1931 {
1932 Integer nb_node_cell_kind = nb_cell_kind;
1933 ENUMERATE_ITEMWITHNODES(iitem, cell_kind_family->allItems())
1934 {
1935 nb_node_cell_kind += (*iitem).nbNode();
1936 }
1937 ofile << "CELLS " << nb_cell_kind << ' ' << nb_node_cell_kind << "\n";
1938 ENUMERATE_ITEMWITHNODES(iitem, cell_kind_family->allItems())
1939 {
1940 ItemWithNodes item = *iitem;
1941 Integer item_nb_node = item.nbNode();
1942 ofile << item_nb_node;
1943 for (NodeLocalId node_id : item.nodes() ) {
1944 ofile << ' ' << nodes_local_id_to_current[node_id];
1945 }
1946 ofile << '\n';
1947 }
1948 // Le type doit être coherent avec celui de vtkCellType.h
1949 ofile << "CELL_TYPES " << nb_cell_kind << "\n";
1950 ENUMERATE_ (ItemWithNodes, iitem, cell_kind_family->allItems()) {
1951 Int16 arcane_type = (*iitem).type();
1952 int type = arcaneToVtkCellType(arcane_type);
1953 ofile << type << '\n';
1954 }
1955 }
1956
1957 // Si on est dans le maillage des mailles, sauve les groupes de noeuds.
1958 if (cell_kind == IK_Cell) {
1959 ofile << "POINT_DATA " << nb_node << "\n";
1960 _saveGroups(mesh->itemFamily(IK_Node), ofile);
1961 }
1962
1963 // Sauve les groupes de mailles
1964 ofile << "CELL_DATA " << nb_cell_kind << "\n";
1965 _saveGroups(mesh->itemFamily(cell_kind), ofile);
1966}
1967
1968/*---------------------------------------------------------------------------*/
1969/*---------------------------------------------------------------------------*/
1970
1971void VtkLegacyMeshWriter::
1972_saveGroups(IItemFamily* family, std::ostream& ofile)
1973{
1974 info() << "Saving groups for family name=" << family->name();
1975 UniqueArray<char> in_group_list(family->maxLocalId());
1976 for (ItemGroupCollection::Enumerator igroup(family->groups()); ++igroup;) {
1977 ItemGroup group = *igroup;
1978 // Inutile de sauver le groupe de toutes les entités
1979 if (group == family->allItems())
1980 continue;
1981 //HACK: a supprimer
1982 if (group.name() == "OuterFaces")
1983 continue;
1984 ofile << "SCALARS GROUP_" << group.name() << " int 1\n";
1985 ofile << "LOOKUP_TABLE default\n";
1986 in_group_list.fill('0');
1987 ENUMERATE_ITEM (iitem, group) {
1988 in_group_list[(*iitem).localId()] = '1';
1989 }
1990 ENUMERATE_ITEM (iitem, family->allItems()) {
1991 ofile << in_group_list[(*iitem).localId()] << '\n';
1992 }
1993 }
1994}
1995
1996/*---------------------------------------------------------------------------*/
1997/*---------------------------------------------------------------------------*/
1998
1999/*---------------------------------------------------------------------------*/
2000/*---------------------------------------------------------------------------*/
2001
2002class VtkMeshReader
2003: public AbstractService
2004, public IMeshReader
2005{
2006 public:
2007
2008 explicit VtkMeshReader(const ServiceBuildInfo& sbi)
2009 : AbstractService(sbi)
2010 {}
2011
2012 public:
2013
2014 bool allowExtension(const String& str) override { return str == "vtk"; }
2015 eReturnType readMeshFromFile(IPrimaryMesh* mesh, const XmlNode& mesh_node, const String& file_name,
2016 const String& dir_name, bool use_internal_partition) override
2017
2018 {
2019 ARCANE_UNUSED(mesh_node);
2020 VtkMeshIOService vtk_service(traceMng());
2021 bool ret = vtk_service.readMesh(mesh, file_name, dir_name, use_internal_partition);
2022 if (ret)
2023 return RTError;
2024
2025 return RTOk;
2026 }
2027};
2028
2029/*---------------------------------------------------------------------------*/
2030/*---------------------------------------------------------------------------*/
2031
2033
2035 ServiceProperty("VtkLegacyMeshReader", ST_SubDomain),
2037
2038/*---------------------------------------------------------------------------*/
2039/*---------------------------------------------------------------------------*/
2040
2041/*---------------------------------------------------------------------------*/
2042/*---------------------------------------------------------------------------*/
2043
2044class VtkLegacyCaseMeshReader
2045: public AbstractService
2046, public ICaseMeshReader
2047{
2048 public:
2049
2050 class Builder
2051 : public IMeshBuilder
2052 {
2053 public:
2054
2055 explicit Builder(ITraceMng* tm, const CaseMeshReaderReadInfo& read_info)
2056 : m_trace_mng(tm)
2057 , m_read_info(read_info)
2058 {}
2059
2060 public:
2061
2062 void fillMeshBuildInfo(MeshBuildInfo& build_info) override
2063 {
2064 ARCANE_UNUSED(build_info);
2065 }
2067 {
2068 VtkMeshIOService vtk_service(m_trace_mng);
2069 String fname = m_read_info.fileName();
2070 m_trace_mng->info() << "VtkLegacy Reader (ICaseMeshReader) file_name=" << fname;
2071 bool ret = vtk_service.readMesh(pm, fname, m_read_info.directoryName(), m_read_info.isParallelRead());
2072 if (ret)
2073 ARCANE_FATAL("Can not read VTK File");
2074 }
2075
2076 private:
2077
2078 ITraceMng* m_trace_mng;
2079 CaseMeshReaderReadInfo m_read_info;
2080 };
2081
2082 public:
2083
2084 explicit VtkLegacyCaseMeshReader(const ServiceBuildInfo& sbi)
2085 : AbstractService(sbi)
2086 {}
2087
2088 public:
2089
2091 {
2092 IMeshBuilder* builder = nullptr;
2093 if (read_info.format() == "vtk")
2094 builder = new Builder(traceMng(), read_info);
2095 return makeRef(builder);
2096 }
2097};
2098
2099/*---------------------------------------------------------------------------*/
2100/*---------------------------------------------------------------------------*/
2101
2103 ServiceProperty("VtkLegacyCaseMeshReader", ST_SubDomain),
2105
2106/*---------------------------------------------------------------------------*/
2107/*---------------------------------------------------------------------------*/
2108
2109} // End namespace Arcane
2110
2111/*---------------------------------------------------------------------------*/
2112/*---------------------------------------------------------------------------*/
#define ARCANE_THROW(exception_class,...)
Macro pour envoyer une exception avec formattage.
#define ARCANE_FATAL(...)
Macro envoyant une exception FatalErrorException.
Types et macros pour itérer sur les entités du maillage.
#define ENUMERATE_FACE(name, group)
Enumérateur générique d'un groupe de faces.
#define ENUMERATE_(type, name, group)
Enumérateur générique d'un groupe d'entité
#define ENUMERATE_ITEM(name, group)
Enumérateur générique d'un groupe de noeuds.
#define ENUMERATE_NODE(name, group)
Enumérateur générique d'un groupe de noeuds.
#define ARCANE_REGISTER_SUB_DOMAIN_FACTORY(aclass, ainterface, aname)
Enregistre un service de fabrique pour la classe aclass.
#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.
Classe de base d'un service.
AbstractService(const ServiceBuildInfo &)
Constructeur à partir d'un ServiceBuildInfo.
Tableau d'items de types quelconques.
void fill(const DataType &data)
Remplissage du tableau.
void resize(Int64 s)
Change le nombre d'éléments du tableau à s.
ArrayView< T > subView(Int64 abegin, Integer asize)
Sous-vue à partir de l'élément abegin et contenant asize éléments.
void copy(Span< const T > rhs)
Copie les valeurs de rhs dans l'instance.
void add(ConstReferenceType val)
Ajoute l'élément val à la fin du tableau.
Informations nécessaires pour la lecture d'un fichier de maillage.
EnumeratorT< ItemGroup > Enumerator
Definition Collection.h:129
Vue constante d'un tableau de type T.
constexpr Integer size() const noexcept
Nombre d'éléments du tableau.
Face d'une maille.
Definition Item.h:952
Informations sur le type flottant.
Definition Limits.h:48
Interface du service de lecture du maillage à partir du jeu de données.
Interface d'une famille d'entités.
Definition IItemFamily.h:84
virtual ItemGroupCollection groups() const =0
Liste des groupes de cette famille.
virtual ItemGroup allItems() const =0
Groupe de toutes les entités.
virtual ItemGroup createGroup(const String &name, Int32ConstArrayView local_ids, bool do_override=false)=0
Créé un groupe d'entités de nom name contenant les entités local_ids.
virtual Int32 maxLocalId() const =0
virtual String name() const =0
Nom de la famille.
virtual Integer nbItem() const =0
Nombre d'entités.
Interface d'un service de création/lecture du maillage.
Interface du service gérant la lecture d'un maillage.
Definition IMeshReader.h:37
eReturnType
Types des codes de retour d'une lecture ou écriture.
Definition IMeshReader.h:42
@ RTError
Erreur lors de l'opération.
Definition IMeshReader.h:44
@ RTOk
Opération effectuée avec succès.
Definition IMeshReader.h:43
Interface d'une classe proposant des fonctions utilitaires sur maillage.
virtual void getFacesLocalIdFromConnectivity(ConstArrayView< ItemTypeId > items_type, ConstArrayView< Int64 > items_connectivity, ArrayView< Int32 > local_ids, bool allow_null=false)=0
Recherche les identifiants locaux des faces à partir de leur connectivité.
Interface d'un service d'écriture d'un maillage.
Definition IMeshWriter.h:36
Exception lorsqu'une erreur d'entrée/sortie est détectée.
Definition IOException.h:32
Interface du gestionnaire de parallélisme pour un sous-domaine.
virtual Int32 commRank() const =0
Rang de cette instance dans le communicateur.
virtual void broadcastMemoryBuffer(ByteArray &bytes, Int32 rank)=0
Effectue un broadcast d'une zone mémoire.
Interface du gestionnaire de traces.
virtual void addAutoDestroyVariable(VariableRef *var)=0
Ajoute la variable à la liste des variables qui sont conservées jusqu'à la fin de l'exécution.
Interface du gestionnaire de variables.
virtual IVariableMngInternal * _internalApi()=0
API interne à Arcane.
static IXmlDocumentHolder * loadFromBuffer(Span< const Byte > buffer, const String &name, ITraceMng *tm)
Charge un document XML.
Definition DomUtils.cc:426
Groupe d'entités de maillage.
Definition ItemGroup.h:49
const String & name() const
Nom du groupe.
Definition ItemGroup.h:76
Type d'une entité (Item).
Definition ItemTypeId.h:32
Int16 dimension() const
Dimension de l'élément (<0 si inconnu)
Gestionnaire des types d'entités d'un maillage.
Definition ItemTypeMng.h:65
ItemTypeInfo * typeFromId(Integer id) const
Type correspondant au numéro id.
Elément de maillage s'appuyant sur des noeuds (Edge,Face,Cell).
Definition Item.h:727
NodeConnectedListViewType nodes() const
Liste des noeuds de l'entité
Definition Item.h:785
Int32 nbNode() const
Nombre de noeuds de l'entité
Definition Item.h:779
constexpr Int32 localId() const
Identifiant local de l'entité dans le sous-domaine du processeur.
Definition Item.h:219
ItemUniqueId uniqueId() const
Identifiant unique sur tous les domaines.
Definition Item.h:225
Paramètres nécessaires à la construction d'un maillage.
Noeud d'un maillage.
Definition Item.h:576
Flot de sortie lié à une String.
Classe gérant un vecteur de réel de dimension 3.
Definition Real3.h:132
Référence à une instance.
Encapsulation d'un pointeur qui se détruit automatiquement.
Definition ScopedPtr.h:44
Structure contenant les informations pour créer un service.
Propriétés de création d'un service.
Chaîne de caractères unicode.
String lower() const
Transforme tous les caractères de la chaîne en minuscules.
Definition String.cc:478
const char * localstr() const
Retourne la conversion de l'instance dans l'encodage UTF-8.
Definition String.cc:227
ByteConstArrayView utf8() const
Retourne la conversion de l'instance dans l'encodage UTF-8.
Definition String.cc:275
bool endsWith(const String &s) const
Indique si la chaîne se termine par les caractères de s.
Definition String.cc:1083
String substring(Int64 pos) const
Sous-chaîne commençant à la position pos.
Definition String.cc:1114
TraceAccessor(ITraceMng *m)
Construit un accesseur via le gestionnaire de trace m.
TraceMessageDbg debug(Trace::eDebugLevel=Trace::Medium) const
Flot pour un message de debug.
TraceMessage fatal() const
Flot pour un message d'erreur fatale.
TraceMessage info() const
Flot pour un message d'information.
TraceMessage error() const
Flot pour un message d'erreur.
ITraceMng * traceMng() const
Gestionnaire de trace.
Vecteur 1D de données avec sémantique par valeur (style STL).
Informations pour allouer les entités d'un maillage non structuré.
void allocateMesh()
Alloue le maillage avec les mailles ajoutées lors de l'appel à addCell().
void setMeshDimension(Int32 v)
Positionne la dimension du maillage.
void preAllocate(Int32 nb_cell, Int64 nb_connectivity_node)
Pre-alloue la mémoire.
void addCell(ItemTypeId type_id, Int64 cell_uid, SmallSpan< const Int64 > nodes_uid)
Ajoute une maille au maillage.
Paramètres nécessaires à la construction d'une variable.
const char * getCurrentLine()
Permet de retourner la ligne présente dans le buffer.
bool isEmptyNextLine()
Permet de voir si la prochaine ligne est vide.
bool m_is_init
Y'a-t-il eu au moins une ligne lue.
std::istream * m_stream
Le stream.
static bool isEqualString(const String &current_value, const String &expected_value)
Permet de vérifier si expected_value == current_value.
double getDouble()
Permet de récupérer le double qui suit.
bool m_is_eof
Est-on à la fin du fichier.
char m_buf[BUFSIZE]
Le buffer contenant la ligne lue.
void getBinary(T &type)
Permet de récupérer le nombre binaire qui suit.
int getInt()
Permet de récupérer le int qui suit.
void checkString(const String &current_value, const String &expected_value)
Permet de vérifier si expected_value == current_value.
bool m_is_binary_file
Est-ce un fichier contenant des données en binaire.
bool m_need_reread_current_line
Doit-on relire la même ligne.
const char * getNextLine()
Permet de récupérer la prochaine ligne du fichier.
float getFloat()
Permet de récupérer le float qui suit.
void allocateMeshItems(IPrimaryMesh *pm) override
Alloue les entités du maillage géré par ce service.
void fillMeshBuildInfo(MeshBuildInfo &build_info) override
Remplit build_info avec les informations nécessaires pour créer le maillage.
Ref< IMeshBuilder > createBuilder(const CaseMeshReaderReadInfo &read_info) const override
Retourne un builder pour créer et lire le maillage dont les informations sont spécifiées dans read_in...
bool writeMeshToFile(IMesh *mesh, const String &file_name) override
Ecriture du maillage au Vtk.
void _writeMeshToFile(IMesh *mesh, const String &file_name, eItemKind cell_kind)
Ecrit le maillage au format Vtk.
void build() override
Construction de niveau build du service.
Lecteur des fichiers de maillage au format Vtk historique (legacy).
void _createFaceGroup(IMesh *mesh, const String &name, Int32ConstArrayView faces_lid)
Permet de créer un groupe de face de nom "name" et composé des faces ayant les ids inclus dans "faces...
bool _readMetadata(IMesh *mesh, VtkFile &vtk_file)
Lecture des metadata.
void _readCellVariable(IMesh *mesh, VtkFile &vtk_file, const String &name_str, Integer nb_cell)
Permet de créer une variable aux mailles à partir des infos du fichier vtk.
bool _readStructuredGrid(IPrimaryMesh *mesh, VtkFile &, bool use_internal_partition)
Permet de lire un fichier vtk contenant une STRUCTURED_GRID.
void _readNodesUnstructuredGrid(IMesh *mesh, VtkFile &vtk_file, Array< Real3 > &node_coords)
Lecture des noeuds et de leur coordonnées.
void _readCellsUnstructuredGrid(IMesh *mesh, VtkFile &vtk_file, Array< Integer > &cells_nb_node, Array< ItemTypeId > &cells_type, Array< Int64 > &cells_connectivity)
Lecture des mailles et de leur connectivité.
void _readItemGroup(IMesh *mesh, VtkFile &vtk_file, const String &name_str, Integer nb_item, eItemKind ik, ConstArrayView< Int32 > local_id)
Permet de créer un groupe d'item.
void _readFacesMesh(IMesh *mesh, const String &file_name, const String &dir_name, bool use_internal_partition)
Permet de lire le fichier truc.vtkfaces.vtk (s'il existe).
bool _readData(IMesh *mesh, VtkFile &vtk_file, bool use_internal_partition, eItemKind cell_kind, Int32ConstArrayView local_id, Integer nb_node)
Permet de lire les données complémentaires (POINT_DATA / CELL_DATA).
bool readMesh(IPrimaryMesh *mesh, const String &file_name, const String &dir_name, bool use_internal_partition)
Permet de débuter la lecture d'un fichier vtk.
void _readNodeGroup(IMesh *mesh, VtkFile &vtk_file, const String &name, Integer nb_item)
Permet de créer un groupe de node.
bool _readUnstructuredGrid(IPrimaryMesh *mesh, VtkFile &vtk_file, bool use_internal_partition)
Permet de lire un fichier vtk contenant une UNSTRUCTURED_GRID.
eReturnType readMeshFromFile(IPrimaryMesh *mesh, const XmlNode &mesh_node, const String &file_name, const String &dir_name, bool use_internal_partition) override
Lit un maillage à partir d'un fichier.
bool allowExtension(const String &str) override
Vérifie si le service supporte les fichiers avec l'extension str.
Liste de noeuds d'un arbre DOM.
Definition XmlNodeList.h:33
Noeud d'un arbre DOM.
Definition XmlNode.h:51
XmlNode documentElement() const
Retourne le noeud élément du document.
Definition XmlNode.cc:556
XmlNodeList children(const String &name) const
Ensemble des noeuds fils de ce noeud ayant pour nom name.
Definition XmlNode.cc:93
#define ARCANE_REGISTER_SERVICE(aclass, a_service_property,...)
Macro pour enregistrer un service.
MeshVariableScalarRefT< Cell, Real > VariableCellReal
Grandeur au centre des mailles de type réel.
MeshVariableScalarRefT< Node, Real3 > VariableNodeReal3
Grandeur au noeud de type coordonnées.
@ ReduceMax
Maximum des valeurs.
-*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
ArrayView< Byte > ByteArrayView
Equivalent C d'un tableau à une dimension de caractères.
Definition UtilsTypes.h:534
std::int64_t Int64
Type entier signé sur 64 bits.
Int32 Integer
Type représentant un entier.
ConstArrayView< Int32 > Int32ConstArrayView
Equivalent C d'un tableau à une dimension d'entiers 32 bits.
Definition UtilsTypes.h:569
@ ST_SubDomain
Le service s'utilise au niveau du sous-domaine.
UniqueArray< Byte > ByteUniqueArray
Tableau dynamique à une dimension de caractères.
Definition UtilsTypes.h:422
UniqueArray< Int32 > Int32UniqueArray
Tableau dynamique à une dimension d'entiers 32 bits.
Definition UtilsTypes.h:428
eItemKind
Genre d'entité de maillage.
@ IK_Node
Entité de maillage de genre noeud.
@ IK_Cell
Entité de maillage de genre maille.
@ IK_Face
Entité de maillage de genre face.
std::int16_t Int16
Type entier signé sur 16 bits.
double Real
Type représentant un réel.
ConstArrayView< Byte > ByteConstArrayView
Equivalent C d'un tableau à une dimension de caractères.
Definition UtilsTypes.h:563
unsigned char Byte
Type d'un octet.
Definition BaseTypes.h:43
auto makeRef(InstanceType *t) -> Ref< InstanceType >
Créé une référence sur un pointeur.
ArrayView< Real > RealArrayView
Equivalent C d'un tableau à une dimension de réels.
Definition UtilsTypes.h:546
std::int32_t Int32
Type entier signé sur 32 bits.
Real y
deuxième composante du triplet
Definition Real3.h:36
Real z
troisième composante du triplet
Definition Real3.h:37
Real x
première composante du triplet
Definition Real3.h:35