14#include "arcane/utils/ArcanePrecomp.h"
15#include "arcane/utils/Iostream.h"
16#include "arcane/utils/StdHeader.h"
17#include "arcane/utils/HashTableMap.h"
18#include "arcane/utils/ValueConvert.h"
19#include "arcane/utils/ScopedPtr.h"
20#include "arcane/utils/ITraceMng.h"
21#include "arcane/utils/String.h"
22#include "arcane/utils/IOException.h"
23#include "arcane/utils/Collection.h"
24#include "arcane/utils/Enumerator.h"
25#include "arcane/utils/NotImplementedException.h"
26#include "arcane/utils/Real3.h"
28#include "arcane/core/AbstractService.h"
29#include "arcane/core/FactoryService.h"
30#include "arcane/core/IMainFactory.h"
31#include "arcane/core/IMeshReader.h"
32#include "arcane/core/ISubDomain.h"
33#include "arcane/core/IMesh.h"
34#include "arcane/core/IMeshSubMeshTransition.h"
35#include "arcane/core/IItemFamily.h"
36#include "arcane/core/Item.h"
38#include "arcane/core/VariableTypes.h"
39#include "arcane/core/IVariableAccessor.h"
40#include "arcane/core/IParallelMng.h"
41#include "arcane/core/IIOMng.h"
42#include "arcane/core/IXmlDocumentHolder.h"
43#include "arcane/core/XmlNodeList.h"
44#include "arcane/core/XmlNode.h"
45#include "arcane/core/IMeshUtilities.h"
46#include "arcane/core/IMeshWriter.h"
47#include "arcane/core/BasicService.h"
48#include "arcane/core/ICaseMeshReader.h"
49#include "arcane/core/IMeshBuilder.h"
50#include "arcane/core/ItemPrinter.h"
51#include "arcane/core/MeshKind.h"
52#include "arcane/core/ServiceBuilder.h"
54#include <vtkXMLUnstructuredGridReader.h>
55#include <vtkUnstructuredGrid.h>
57#include <vtkDataObject.h>
58#include <vtkDataSet.h>
59#include <vtkCellData.h>
60#include <vtkPointData.h>
61#include <vtkLongArray.h>
62#include <vtkIntArray.h>
68#define VTK_FILE_EXT_VTI "vti"
69#define VTK_FILE_EXT_VTP "vtp"
70#define VTK_FILE_EXT_VTR "vtr"
71#define VTK_FILE_EXT_VTS "vts"
72#define VTK_FILE_EXT_VTU "vtu"
73#define VTK_FILE_EXT_PVTI "pvti"
74#define VTK_FILE_EXT_PVTP "pvtp"
75#define VTK_FILE_EXT_PVTR "pvtr"
76#define VTK_FILE_EXT_PVTS "pvts"
77#define VTK_FILE_EXT_PVTU "pvtu"
79#if VTK_MAJOR_VERSION >= 9
80#define CURRENT_VTK_VERSION_LONG_TYPE VTK_LONG_LONG
81#include <vtkLongLongArray.h>
82using vtkLongArrayType = vtkLongLongArray;
84#define CURRENT_VTK_VERSION_LONG_TYPE VTK_LONG
85using vtkLongArrayType = vtkLongArray;
104class VtuMeshReaderBase
108 explicit VtuMeshReaderBase(
ITraceMng* trace_mng);
109 virtual ~VtuMeshReaderBase() =
default;
113 virtual void build() {}
116 const String& file_name,
const String& dir_name,
bool use_internal_partition);
120 bool readGroupsFromFieldData(
IMesh *
mesh, vtkFieldData*,
int);
141 {VTK_FILE_EXT_VTU, &VtuMeshReaderBase::readMeshFromVtuFile},
150VtuMeshReaderBase(
ITraceMng* trace_mng) : m_trace_mng(trace_mng)
161bool VtuMeshReaderBase::
162readGroupsFromFieldData(
IMesh*
mesh,vtkFieldData* allFieldData,
int i)
164 const std::string group_name = allFieldData->GetArrayName(i);
165 vtkDataArray* iFieldData = allFieldData->GetArray(i);
169 if (iFieldData->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
172 auto nb_tuple = iFieldData->GetNumberOfTuples();
173 m_trace_mng->info() <<
"[readGroupsFromFieldData] iFieldData->GetNumberOfTuples="<< nb_tuple;
176 vtkLongArrayType* vtk_array = vtkLongArrayType::SafeDownCast(iFieldData);
180 IItemFamily* family = mesh->itemFamily(kind_type);
181 auto nb_item = nb_tuple - 1;
184 for(
Integer z=0; z<nb_item; ++z )
185 unique_ids[z] = vtk_array->GetValue(z+1);
189 family->itemsUniqueIdToLocalId(local_ids,unique_ids,
false);
194 for(
Integer index = 0; index < nb_item; ++index )
195 if (local_ids[index]!=NULL_ITEM_LOCAL_ID)
196 ids.
add(local_ids[index]);
198 m_trace_mng->info() <<
"Create group family=" << family->name() <<
" name=" << group_name <<
" ids=" << ids.size();
199 family->createGroup(group_name,ids);
212 bool use_internal_partition)
214 ARCANE_UNUSED(dir_name);
215 ARCANE_UNUSED(use_internal_partition);
216 bool itWasAnArcanProduction=
true;
217 IParallelMng* pm = mesh->parallelMng();
218 bool is_parallel = pm->isParallel();
221 m_trace_mng->info() <<
"[readMeshFromVtuFile] Entering";
222 vtkXMLUnstructuredGridReader *reader = vtkXMLUnstructuredGridReader::New();
223 std::string fname = file_name.localstr();
224 reader->SetFileName(fname.c_str());
225 if (!reader->CanReadFile(fname.c_str()))
227 reader->UpdateInformation();
229 vtkUnstructuredGrid *unstructuredGrid = reader->GetOutput();
232 auto nbOfCells = unstructuredGrid->GetNumberOfCells();
233 auto nbOfNodes = unstructuredGrid->GetNumberOfPoints();
239 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Now Fetching Nodes Unique IDs ##";
240 vtkPointData* allPointData=unstructuredGrid->GetPointData();
241 allPointData->Update();
242 vtkDataArray* dataNodeArray = allPointData->GetArray(
"NodesUniqueIDs");
243 vtkLongArrayType *nodesUidArray =
nullptr;
245 m_trace_mng->info() <<
"[readMeshFromVtuFile] Could not be found, creating new one";
246 nodesUidArray = vtkLongArrayType::New();
247 for(
Integer uid=0; uid<nbOfNodes; ++uid)
248 nodesUidArray->InsertNextValue(uid);
249 itWasAnArcanProduction=
false;
252 m_trace_mng->info() <<
"dataNodeArray->GetDataType()" << dataNodeArray->GetDataType();
253 if (dataNodeArray->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
255 nodesUidArray = vtkLongArrayType::SafeDownCast(dataNodeArray);
257 m_trace_mng->info() <<
"[readMeshFromVtuFile] Fetched";
259 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Now Fetching Cells Unique IDs ##";
260 vtkCellData* allCellData=unstructuredGrid->GetCellData();
261 allCellData->Update();
262 vtkDataArray* dataCellArray = allCellData->GetArray(
"CellsUniqueIDs");
263 vtkLongArrayType *cellsUidArray =
nullptr;
265 m_trace_mng->info() <<
"[readMeshFromVtuFile] Could not be found, creating new one";
266 cellsUidArray = vtkLongArrayType::New();
267 for(
Integer uid=0; uid<nbOfCells; ++uid)
268 cellsUidArray->InsertNextValue(uid);
269 itWasAnArcanProduction=
false;
272 if (dataCellArray->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
274 cellsUidArray = vtkLongArrayType::SafeDownCast(dataCellArray);
279 vtkDataArray* data_owner_array = allPointData->GetArray(
"NodesOwner");
280 vtkIntArray* nodes_owner_array =
nullptr;
281 if (data_owner_array){
282 if (data_owner_array->GetDataType() == VTK_INT) {
283 nodes_owner_array = vtkIntArray::SafeDownCast(data_owner_array);
286 m_trace_mng->warning() <<
"NodesOwner array is present but has bad type (not Int32)";
289 m_trace_mng->info() <<
"[readMeshFromVtuFile] Fetched";
295 m_trace_mng->info() <<
"[readMeshFromVtuFile] nbOfCells=" << nbOfCells <<
", nbOfNodes=" << nbOfNodes;
297 HashTableMapT<Int64,Real3> nodes_coords(
Integer(nbOfNodes),
true);
298 vtkPoints* vtkAllPoints=unstructuredGrid->GetPoints();
300 mesh->toPrimaryMesh()->setDimension(vtkAllPoints->GetData()->GetNumberOfComponents());
301 for(vtkIdType i=0; i<nbOfNodes; ++i ){
303 vtkAllPoints->GetPoint(i,xyz);
305 coords[i] = Real3(xyz[0],xyz[1],xyz[2]);
306 nodes_coords.nocheckAdd(nodesUidArray->GetValue(i),coords[i]);
310 HashTableMapT<Int64,Int32> nodes_owner_map(
Integer(nbOfNodes),
true);
311 if (nodes_owner_array){
312 for(
Integer i=0; i<nbOfNodes; ++i ){
313 nodes_owner_map.nocheckAdd(nodesUidArray->GetValue(i),nodes_owner_array->GetValue(i));
318 cells_filter.
resize(nbOfCells);
319 for(
Integer i=0; i<nbOfCells; ++i )
323 vtkIdType mesh_nb_cell_node = 0;
324 for(
Integer j=0, js=cells_filter.size(); j<js; ++j )
325 mesh_nb_cell_node += unstructuredGrid->GetCell(j)->GetNumberOfPoints();
326 m_trace_mng->info() <<
"Number of mesh_nb_cell_node = "<<mesh_nb_cell_node;
336 for(
Integer i_cell=0, s_cell=cells_filter.size(); i_cell<s_cell; ++i_cell ){
337 vtkIdType iVtkCell=i_cell;
338 vtkCell *cell = unstructuredGrid->GetCell(iVtkCell);
339 if (!cell->IsLinear())
throw NotImplementedException(A_FUNCINFO);
340 Integer arcItemType = IT_NullType;
341 switch(cell->GetCellType()){
343 case(VTK_TETRA): arcItemType = IT_Tetraedron4;
break;
344 case(VTK_HEXAHEDRON): arcItemType = IT_Hexaedron8;
break;
345 case(VTK_PYRAMID): arcItemType = IT_Pyramid5;
break;
347 case(VTK_WEDGE): arcItemType = IT_Pentaedron6;
break;
348 case(VTK_PENTAGONAL_PRISM): arcItemType = IT_Heptaedron10;
break;
349 case(VTK_HEXAGONAL_PRISM): arcItemType = IT_Octaedron12;
break;
352 if (mesh->dimension() == 2){ arcItemType = IT_Quad4;}
353 else if (mesh->dimension() == 3 && !mesh->meshKind().isMonoDimension()) {arcItemType = IT_Cell3D_Quad4;}
354 else {
ARCANE_FATAL(
"VTK_QUAD is not supported in mono-dimension meshes of dimension {0}",mesh->dimension());}
357 if (mesh->dimension() == 2) {arcItemType = IT_Triangle3;}
358 else if (mesh->dimension() == 3 && !mesh->meshKind().isMonoDimension()) {arcItemType = IT_Cell3D_Triangle3;}
359 else {
ARCANE_FATAL(
"VTK_TRIANGLE is not supported in mono-dimension meshes of dimension {0}",mesh->dimension());}
362 case(VTK_POLY_VERTEX):
363 switch (cell->GetNumberOfPoints()){
364 case (4): arcItemType = IT_Tetraedron4;
break;
365 case (7): arcItemType = IT_HemiHexa7;
break;
366 case (6): arcItemType = IT_HemiHexa6;
break;
367 case (5): arcItemType = IT_HemiHexa5;
break;
368 default:
throw NotImplementedException(A_FUNCINFO);;
372 switch (cell->GetNumberOfPoints()){
373 case (4): arcItemType = IT_Tetraedron4;
break;
374 default:
throw NotImplementedException(A_FUNCINFO);;
377#ifndef NO_USER_WARNING
378#warning IT_Vertex returns 0 nodes in ItemTypeMng.cc@42:type->setInfos(this,IT_Vertex,0,0,0);
379#warning IT_Vertex vs IT_DualNode HACK
381 case(VTK_VERTEX): arcItemType=IT_DualNode;
break;
383 switch (cell->GetNumberOfPoints()){
385 if (mesh->dimension() == 2) {arcItemType = IT_Triangle3;}
386 else if (mesh->dimension() == 3 && !mesh->meshKind().isMonoDimension()) {arcItemType = IT_Cell3D_Triangle3;}
387 else {
ARCANE_FATAL(
"VTK_TRIANGLE is not supported in mono-dimension meshes of dimension {0}",mesh->dimension());}
390 if (mesh->dimension() == 2){ arcItemType = IT_Quad4;}
391 else if (mesh->dimension() == 3 && !mesh->meshKind().isMonoDimension()) {arcItemType = IT_Cell3D_Quad4;}
392 else {
ARCANE_FATAL(
"VTK_QUAD is not supported in mono-dimension meshes of dimension {0}",mesh->dimension());}
394 case (5): arcItemType = IT_Pentagon5;
break;
395 case (6): arcItemType = IT_Hexagon6;
break;
396 default:
throw NotImplementedException(A_FUNCINFO);;
401 switch (cell->GetNumberOfPoints()){
402 case (2): arcItemType = IT_CellLine2;
break;
403 default:
throw NotImplementedException(A_FUNCINFO);;
406 case(VTK_EMPTY_CELL):
407 m_trace_mng->info() <<
"VTK_EMPTY_CELL n="<<cell->GetNumberOfPoints();
409 case(VTK_TRIANGLE_STRIP):
410 m_trace_mng->info() <<
"VTK_TRIANGLE_STRIP n="<<cell->GetNumberOfPoints();
413 m_trace_mng->info() <<
"VTK_VOXEL n="<<cell->GetNumberOfPoints();
414 switch (cell->GetNumberOfPoints()){
415 case (3): arcItemType = IT_Triangle3;
break;
416 case (4): arcItemType = IT_Quad4;
break;
417 case (5): arcItemType = IT_Pentagon5;
break;
418 case (6): arcItemType = IT_Hexagon6;
break;
419 case (7): arcItemType = IT_HemiHexa7;
break;
420 case (8): arcItemType = IT_Hexaedron8;
break;
421 default:
throw NotImplementedException(A_FUNCINFO);;
424 default:
throw NotImplementedException(A_FUNCINFO);;
426 auto nNodes = cell->GetNumberOfPoints();
429 cells_infos[cells_infos_index++]=arcItemType;
432 Integer cell_indirect_id = cells_filter[i_cell];
433 cells_infos[cells_infos_index++] = cellsUidArray->GetValue(cell_indirect_id);
436 vtkIdList* nodeIds=cell->GetPointIds();
437 for(
int iNode=0; iNode<nNodes; ++iNode){
438 auto localId = nodeIds->GetId(iNode);
439 long uniqueUid=nodesUidArray->GetValue(localId);
441 cells_infos[cells_infos_index++] = uniqueUid;
449 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Mesh 3D ##";
450 PRIMARYMESH_CAST(mesh)->setDimension(3);
451 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Allocating " << cells_filter.size() <<
" cells ##";
452 PRIMARYMESH_CAST(mesh)->allocateCells(cells_filter.size(), cells_infos,
false);
457 m_trace_mng->info() <<
"[readMeshFromVtuFile] internalNodes.size()="<<internalNodes.size();
458 if (nodes_owner_array){
459 m_trace_mng->info() <<
"Set nodes owners from vtu file";
460 for(
Integer i=0, is=internalNodes.size(); i<is; ++i){
461 ItemInternal* internal_node = internalNodes[i];
462 Int32 true_owner = nodes_owner_map[internal_node->uniqueId()];
463 internal_node->setOwner(true_owner,sid);
467 for(
Integer i=0, is=internalNodes.size(); i<is; ++i)
468 internalNodes[i]->setOwner(sid,sid);
471 m_trace_mng->info() <<
"[readMeshFromVtuFile] internalCells.size()="<<internalCells.size();
472 for(
Integer i=0, is=internalCells.size(); i<is; ++i)
473 internalCells[i]->setOwner(sid,sid);
479 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Ending with endAllocate ##";
480 PRIMARYMESH_CAST(mesh)->endAllocate();
483 mesh->utilities()->changeOwnersFromCells();
486 m_trace_mng->info() <<
"\n[readMeshFromVtuFile] ## Now dealing with ghost's layer ##";
487 m_trace_mng->info() <<
"[readMeshFromVtuFile] mesh.nbNode=" <<mesh->nbNode() <<
" mesh.nbCell="<< mesh->nbCell();
493 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Now Fetching Other Fields ##";
494 vtkFieldData* allFieldData=unstructuredGrid->GetFieldData();
495 int nbrOfArrays = allFieldData->GetNumberOfArrays();
496 m_trace_mng->info() <<
"[readMeshFromVtuFile] nbrOfArrays = " << nbrOfArrays;
497 for(
int i=0;i<nbrOfArrays;++i){
498 if (itWasAnArcanProduction==
false)
continue;
499 m_trace_mng->info() <<
"[readMeshFromVtuFile] Focussing on \"" << allFieldData->GetArrayName(i) <<
"\" (i="<<i<<
")";
500 if (readGroupsFromFieldData(mesh, allFieldData, i) !=
true)
508 m_trace_mng->info() <<
"[readMeshFromVtuFile] ## Now insert coords ##";
512 nodes_coord_var[inode] = nodes_coords[inode->uniqueId()];
519 mesh->synchronizeGroupsAndVariables();
525 m_trace_mng->info() <<
"[readMeshFromVtuFile] RTOk";
529 nodesUidArray->Delete();
531 cellsUidArray->Delete();
556 const String& file_name,
const String& dir_name,
bool use_internal_partition)
override;
572 for(vtkFileExtIdx=0;vtkFileExtReader[vtkFileExtIdx].ext!=
nullptr;++vtkFileExtIdx){
574 if (str == vtkFileExtReader[vtkFileExtIdx].ext){
593 bool use_internal_partition)
595 ARCANE_UNUSED(mesh_node);
596 info() <<
"[readMeshFromFile] Forwarding to vtkFileExtReader[" << vtkFileExtIdx <<
"].reader";
599 return vtu_mesh_reader.readMeshFromVtuFile(
mesh, file_name, dir_name, use_internal_partition);
619class VtuCaseMeshReader
632 , m_read_info(read_info)
639 ARCANE_UNUSED(build_info);
644 String fname = m_read_info.fileName();
645 m_trace_mng->info() <<
"Vtu Reader (ICaseMeshReader) file_name=" << fname;
646 bool ret = vtu_service.readMeshFromVtuFile(pm, fname, m_read_info.directoryName(), m_read_info.isParallelRead());
668 if (read_info.format() ==
"vtu")
#define ARCANE_FATAL(...)
Macro envoyant une exception FatalErrorException.
#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.
Classe de base d'un service.
AbstractService(const ServiceBuildInfo &)
Constructeur à partir d'un ServiceBuildInfo.
void resize(Int64 s)
Change le nombre d'éléments du tableau à s.
void add(ConstReferenceType val)
Ajoute l'élément val à la fin du tableau.
Informations nécessaires pour la lecture d'un fichier de maillage.
Interface du service de lecture du maillage à partir du jeu de données.
virtual ITraceMng * traceMng()=0
Gestionnaire de message associé
Interface d'un service de création/lecture du maillage.
Interface du service gérant la lecture d'un maillage.
eReturnType
Types des codes de retour d'une lecture ou écriture.
@ RTError
Erreur lors de l'opération.
@ RTOk
Opération effectuée avec succès.
Interface du gestionnaire d'un sous-domaine.
Interface du gestionnaire de traces.
Paramètres nécessaires à la construction d'un maillage.
Référence à une instance.
ISubDomain * subDomain() const
Accès au ISubDomain associé.
Structure contenant les informations pour créer un service.
Propriétés de création d'un service.
Chaîne de caractères unicode.
TraceMessage info() const
Flot pour un message d'information.
ITraceMng * traceMng() const
Gestionnaire de trace.
void fillMeshBuildInfo(MeshBuildInfo &build_info) override
Remplit build_info avec les informations nécessaires pour créer le maillage.
void allocateMeshItems(IPrimaryMesh *pm) override
Alloue les entités du maillage géré par ce service.
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...
Lecteur des fichiers de maillage aux format Vtk.
bool allowExtension(const String &str) override
Vérifie si le service supporte les fichiers avec l'extension str.
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.
#define ARCANE_REGISTER_SERVICE(aclass, a_service_property,...)
Macro pour enregistrer un service.
MeshVariableScalarRefT< Node, Real3 > VariableNodeReal3
Grandeur au noeud de type coordonnées.
-*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
UniqueArray< Int64 > Int64UniqueArray
Tableau dynamique à une dimension d'entiers 64 bits.
Int32 Integer
Type représentant un entier.
UniqueArray< Real3 > Real3UniqueArray
Tableau dynamique à une dimension de vecteurs de rang 3.
ConstArrayView< ItemInternal * > ItemInternalList
Type de la liste interne des entités.
@ ST_SubDomain
Le service s'utilise au niveau du sous-domaine.
UniqueArray< Int32 > Int32UniqueArray
Tableau dynamique à une dimension d'entiers 32 bits.
eItemKind
Genre d'entité de maillage.
@ IK_Node
Entité de maillage de genre noeud.
@ IK_Cell
Entité de maillage de genre maille.
auto makeRef(InstanceType *t) -> Ref< InstanceType >
Créé une référence sur un pointeur.
UniqueArray< Integer > IntegerUniqueArray
Tableau dynamique à une dimension d'entiers.
std::int32_t Int32
Type entier signé sur 32 bits.