Arcane  v3.16.8.0
Documentation développeur
Chargement...
Recherche...
Aucune correspondance
VtuMeshReader.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/* VtuMeshReader.cc (C) 2000-2025 */
9/* */
10/* Lecture/Ecriture d'un fichier au format VtuMeshReader. */
11/*---------------------------------------------------------------------------*/
12/*---------------------------------------------------------------------------*/
13
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"
27
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"
53
54#include <vtkXMLUnstructuredGridReader.h>
55#include <vtkUnstructuredGrid.h>
56#include <vtkCell.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>
63
64
65/**********************************************************************
66 * Convention data types/file type and particular file extension parity *
67 **********************************************************************/
68#define VTK_FILE_EXT_VTI "vti" // Serial structured vtkImageData
69#define VTK_FILE_EXT_VTP "vtp" // Serial UNstructured vtkPolyData
70#define VTK_FILE_EXT_VTR "vtr" // Serial structured vtkRectilinearGrid
71#define VTK_FILE_EXT_VTS "vts" // Serial structured vtkStructuredGrid
72#define VTK_FILE_EXT_VTU "vtu" // Serial UNstructured vtkUnstructuredGrid
73#define VTK_FILE_EXT_PVTI "pvti" // Parallel structured vtkImageData
74#define VTK_FILE_EXT_PVTP "pvtp" // Parallel UNstructured vtkPolyData
75#define VTK_FILE_EXT_PVTR "pvtr" // Parallel structured vtkRectilinearGrid
76#define VTK_FILE_EXT_PVTS "pvts" // Parallel structured vtkStructuredGrid
77#define VTK_FILE_EXT_PVTU "pvtu" // Parallel UNstructured vtkUnstructuredGrid
78
79#if VTK_MAJOR_VERSION >= 9
80#define CURRENT_VTK_VERSION_LONG_TYPE VTK_LONG_LONG
81#include <vtkLongLongArray.h>
82using vtkLongArrayType = vtkLongLongArray;
83#else
84#define CURRENT_VTK_VERSION_LONG_TYPE VTK_LONG
85using vtkLongArrayType = vtkLongArray;
86#endif
87
88
89/*---------------------------------------------------------------------------*/
90/*---------------------------------------------------------------------------*/
91
92namespace Arcane
93{
94
95/*---------------------------------------------------------------------------*/
96/*---------------------------------------------------------------------------*/
97
98
99/*---------------------------------------------------------------------------*/
100/*---------------------------------------------------------------------------*/
104class VtuMeshReaderBase
105{
106 public:
107
108 explicit VtuMeshReaderBase(ITraceMng* trace_mng);
109 virtual ~VtuMeshReaderBase() = default;
110
111 public:
112
113 virtual void build() {}
114
115 IMeshReader::eReturnType readMeshFromVtuFile(IMesh* mesh,
116 const String& file_name, const String& dir_name, bool use_internal_partition);
117
118 // ISubDomain* subDomain() { return m_trace_mng; }
119
120 bool readGroupsFromFieldData(IMesh *mesh, vtkFieldData*, int);
121
122 private:
123
124 int vtkFileExtIdx;
125 ITraceMng* m_trace_mng;
126 // IMeshReader::eReturnType writeMeshToLemFile(IMesh* mesh, const String& filename,const String& dir_name);
127};
128
129/*---------------------------------------------------------------------------*/
130/*---------------------------------------------------------------------------*/
131
132/****************************************************************************
133 * Assignation des data types/file-type and particular file extension readers *
134\****************************************************************************/
135typedef struct{
136 char *ext;
137 IMeshReader::eReturnType (VtuMeshReaderBase::*reader)(IMesh*,const String&,const String&,bool);
139
140vtkExtReader vtkFileExtReader[]={
141 {VTK_FILE_EXT_VTU, &VtuMeshReaderBase::readMeshFromVtuFile},
142 {nullptr,nullptr}
143};
144
145
146/*---------------------------------------------------------------------------*/
147/*---------------------------------------------------------------------------*/
148
149VtuMeshReaderBase::
150VtuMeshReaderBase(ITraceMng* trace_mng) : m_trace_mng(trace_mng)
151{
152 vtkFileExtIdx=0;
153}
154
155/*---------------------------------------------------------------------------*/
156/*---------------------------------------------------------------------------*/
157
158/*****************************************************************
159 * readGroupsFromFieldData
160 *****************************************************************/
161bool VtuMeshReaderBase::
162readGroupsFromFieldData(IMesh* mesh,vtkFieldData* allFieldData,int i)
163{
164 const std::string group_name = allFieldData->GetArrayName(i);
165 vtkDataArray* iFieldData = allFieldData->GetArray(i);
166 if (!iFieldData)
167 return false;
168
169 if (iFieldData->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
170 return false;
171
172 auto nb_tuple = iFieldData->GetNumberOfTuples();
173 m_trace_mng->info() << "[readGroupsFromFieldData] iFieldData->GetNumberOfTuples="<< nb_tuple;
174 if (nb_tuple==0)
175 return false;
176 vtkLongArrayType* vtk_array = vtkLongArrayType::SafeDownCast(iFieldData);
177
178 // Le premier élément du tableau contient son type
179 eItemKind kind_type = (eItemKind)(vtk_array->GetValue(0));
180 IItemFamily* family = mesh->itemFamily(kind_type);
181 auto nb_item = nb_tuple - 1;
182 Int64UniqueArray unique_ids(nb_item);
183 // Les éléments suivant contiennent les uniqueId() des entités du groupe.
184 for( Integer z=0; z<nb_item; ++z )
185 unique_ids[z] = vtk_array->GetValue(z+1);
186
187 // Récupère le localId() correspondant.
188 Int32UniqueArray local_ids(unique_ids.size());
189 family->itemsUniqueIdToLocalId(local_ids,unique_ids,false);
190
191 // Tous les entités ne sont pas forcément dans le maillage actuel et
192 // il faut donc les filtrer.
194 for( Integer index = 0; index < nb_item; ++index )
195 if (local_ids[index]!=NULL_ITEM_LOCAL_ID)
196 ids.add(local_ids[index]);
197
198 m_trace_mng->info() << "Create group family=" << family->name() << " name=" << group_name << " ids=" << ids.size();
199 family->createGroup(group_name,ids);
200 return true;
201}
202
203/*---------------------------------------------------------------------------*/
204/*---------------------------------------------------------------------------*/
205
206/*****************************************************************
207 * readMeshFromVtuFile
208 *****************************************************************/
209IMeshReader::eReturnType VtuMeshReaderBase::
210readMeshFromVtuFile(IMesh* mesh,
211 const String& file_name, const String& dir_name,
212 bool use_internal_partition)
213{
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();
219 Integer sid = pm->commRank();
220
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();
228 reader->Update();// Force reading
229 vtkUnstructuredGrid *unstructuredGrid = reader->GetOutput();
230 // Avec VTK 7, plus besoin du Update.
231 //unstructuredGrid->Update();// La lecture effective du fichier n'a lieu qu'après l'appel à Update().
232 auto nbOfCells = unstructuredGrid->GetNumberOfCells();
233 auto nbOfNodes = unstructuredGrid->GetNumberOfPoints();
234
235
236 /*******************************
237 *Fetching Nodes UID & Cells UID *
238 *******************************/
239 m_trace_mng->info() << "[readMeshFromVtuFile] ## Now Fetching Nodes Unique IDs ##";
240 vtkPointData* allPointData=unstructuredGrid->GetPointData(); // Data associated to Points
241 allPointData->Update();
242 vtkDataArray* dataNodeArray = allPointData->GetArray("NodesUniqueIDs");
243 vtkLongArrayType *nodesUidArray = nullptr;
244 if (!dataNodeArray){
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;
250 }
251 else {
252 m_trace_mng->info() << "dataNodeArray->GetDataType()" << dataNodeArray->GetDataType();
253 if (dataNodeArray->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
255 nodesUidArray = vtkLongArrayType::SafeDownCast(dataNodeArray);
256 }
257 m_trace_mng->info() << "[readMeshFromVtuFile] Fetched";
258
259 m_trace_mng->info() << "[readMeshFromVtuFile] ## Now Fetching Cells Unique IDs ##";
260 vtkCellData* allCellData=unstructuredGrid->GetCellData(); // Data associated to Points
261 allCellData->Update();
262 vtkDataArray* dataCellArray = allCellData->GetArray("CellsUniqueIDs");
263 vtkLongArrayType *cellsUidArray = nullptr;
264 if (!dataCellArray){
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;
270 }
271 else {
272 if (dataCellArray->GetDataType() != CURRENT_VTK_VERSION_LONG_TYPE)
274 cellsUidArray = vtkLongArrayType::SafeDownCast(dataCellArray);
275 }
276
277 // Tableau contenant les numéros des propriétaires des noeuds.
278 // Ce tableau est optionnel et est utilisé par le partitionneur
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);
284 }
285 else {
286 m_trace_mng->warning() << "NodesOwner array is present but has bad type (not Int32)";
287 }
288 }
289 m_trace_mng->info() << "[readMeshFromVtuFile] Fetched";
290
291
292 /************************
293 * Fetch own nodes coords *
294 ************************/
295 m_trace_mng->info() << "[readMeshFromVtuFile] nbOfCells=" << nbOfCells << ", nbOfNodes=" << nbOfNodes;
296 Real3UniqueArray coords(nbOfNodes);
297 HashTableMapT<Int64,Real3> nodes_coords(Integer(nbOfNodes),true);
298 vtkPoints* vtkAllPoints=unstructuredGrid->GetPoints();
299 // Get mesh dimension
300 mesh->toPrimaryMesh()->setDimension(vtkAllPoints->GetData()->GetNumberOfComponents());
301 for(vtkIdType i=0; i<nbOfNodes; ++i ){
302 double xyz[3];
303 vtkAllPoints->GetPoint(i,xyz);
304 //m_trace_mng->info() << "x=" << xyz[0] << " y=" << xyz[1]<< " z=" << xyz[2] << ", nodUid=" << nodesUidArray->GetValue(i);
305 coords[i] = Real3(xyz[0],xyz[1],xyz[2]);
306 nodes_coords.nocheckAdd(nodesUidArray->GetValue(i),coords[i]);
307 }
308
309 // Create hash table for nodes owner.
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));
314 }
315 }
316
317 IntegerUniqueArray cells_filter;
318 cells_filter.resize(nbOfCells);
319 for( Integer i=0; i<nbOfCells; ++i )
320 cells_filter[i] = i;
321
322 // Calcul le nombre de mailles/noeuds
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;
327
328 // Tableau contenant les infos aux mailles (voir IMesh::allocateMesh())
329 Int64UniqueArray cells_infos(mesh_nb_cell_node+cells_filter.size()*2);
330 Integer cells_infos_index = 0;
331
332 /******************
333 * Now insert CELLS *
334 ******************/
335 // For all cells that are discovered
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()){
342 // Linear cells
343 case(VTK_TETRA): arcItemType = IT_Tetraedron4; break;
344 case(VTK_HEXAHEDRON): arcItemType = IT_Hexaedron8; break;
345 case(VTK_PYRAMID): arcItemType = IT_Pyramid5; break;
346 /* Others not yet implemented */
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;
350 /* 2D */
351 case(VTK_QUAD):
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());}
355 break;
356 case(VTK_TRIANGLE):
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());}
360 break;
361 // Cas du poly vertex à parser
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);;
369 }
370 break;
371 case(VTK_PIXEL):
372 switch (cell->GetNumberOfPoints()){
373 case (4): arcItemType = IT_Tetraedron4; break;
374 default:throw NotImplementedException(A_FUNCINFO);;
375 }
376 break;
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
380#endif /* NO_USER_WARNING */
381 case(VTK_VERTEX): arcItemType=IT_DualNode; break;
382 case(VTK_POLYGON):
383 switch (cell->GetNumberOfPoints()){
384 case (3):
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());}
388 break;
389 case (4):
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());}
393 break;
394 case (5): arcItemType = IT_Pentagon5; break;
395 case (6): arcItemType = IT_Hexagon6; break;
396 default:throw NotImplementedException(A_FUNCINFO);;
397 }
398 break;
399 case(VTK_LINE):
400 case(VTK_POLY_LINE):
401 switch (cell->GetNumberOfPoints()){
402 case (2): arcItemType = IT_CellLine2; break;
403 default:throw NotImplementedException(A_FUNCINFO);;
404 }
405 break;
406 case(VTK_EMPTY_CELL):
407 m_trace_mng->info() << "VTK_EMPTY_CELL n="<<cell->GetNumberOfPoints();
408 break;
409 case(VTK_TRIANGLE_STRIP):
410 m_trace_mng->info() << "VTK_TRIANGLE_STRIP n="<<cell->GetNumberOfPoints();
411 break;
412 case(VTK_VOXEL):
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);;
422 }
423 break;
424 default:throw NotImplementedException(A_FUNCINFO);;
425 }
426 auto nNodes = cell->GetNumberOfPoints();// Return the number of points in the cell
427
428 // First is cell's TYPE Stocke le type de la maille
429 cells_infos[cells_infos_index++]=arcItemType;
430
431 // Then comes its UniqueID Stocke le numéro unique de la maille
432 Integer cell_indirect_id = cells_filter[i_cell];
433 cells_infos[cells_infos_index++] = cellsUidArray->GetValue(cell_indirect_id);
434
435 // And finally the Nodes' unique IDs
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);
440 //info() << "working on localId=" << localId << ", uniqueUid=" << uniqueUid;
441 cells_infos[cells_infos_index++] = uniqueUid;
442 }
443 }
444
445
446 /********************************
447 * Setting Dimension & Allocating *
448 ********************************/
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);
453
454
455 // Positionne les propriétaires des noeuds à partir des groupes de noeuds
456 ItemInternalList internalNodes(mesh->itemsInternal(IK_Node));
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);
464 }
465 }
466 else
467 for(Integer i=0, is=internalNodes.size(); i<is; ++i)
468 internalNodes[i]->setOwner(sid,sid);
469
470 ItemInternalList internalCells(mesh->itemsInternal(IK_Cell));
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);
474
475
476 /********************************************
477 * Now finishing & preparing for ghost layout *
478 ********************************************/
479 m_trace_mng->info() << "[readMeshFromVtuFile] ## Ending with endAllocate ##";
480 PRIMARYMESH_CAST(mesh)->endAllocate();
481 if (is_parallel) {
482 // mesh->setOwnersFromCells();
483 mesh->utilities()->changeOwnersFromCells();
484 }
485
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();
488
489
490 /***********************
491 * Fetching Other Groups *
492 ***********************/
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)
502 }
503
504
505 /*******************
506 * Now insert coords *
507 *******************/
508 m_trace_mng->info() << "[readMeshFromVtuFile] ## Now insert coords ##";
509 // Remplit la variable contenant les coordonnées des noeuds
510 VariableNodeReal3& nodes_coord_var(PRIMARYMESH_CAST(mesh)->nodesCoordinates());
511 ENUMERATE_NODE(inode,mesh->ownNodes()){
512 nodes_coord_var[inode] = nodes_coords[inode->uniqueId()];
513 }
514
515
516 /****************************************
517 * Synchronizing groups/variables & nodes *
518 ****************************************/
519 mesh->synchronizeGroupsAndVariables();
520
521 /**************
522 * Finishing up *
523 **************/
524 reader->Delete();
525 m_trace_mng->info() << "[readMeshFromVtuFile] RTOk";
526
527 // TODO: regarder comment detruire automatiquement
528 if (!dataNodeArray)
529 nodesUidArray->Delete();
530 if (!dataCellArray)
531 cellsUidArray->Delete();
532
533 return IMeshReader::RTOk;
534}
535
536/*---------------------------------------------------------------------------*/
537/*---------------------------------------------------------------------------*/
538
539class VtuMeshReader
540: public AbstractService
541, public IMeshReader
542{
543public:
544
545 explicit VtuMeshReader(const ServiceBuildInfo& sbi)
546 : AbstractService(sbi)
547 , vtkFileExtIdx(-1)
548 , m_sub_domain(sbi.subDomain())
549 {}
550
551public:
552
553 bool allowExtension(const String& str) override;
554
556 const String& file_name, const String& dir_name,bool use_internal_partition) override;
557
558private:
559
560 int vtkFileExtIdx;
561 ISubDomain* m_sub_domain;
562 // eReturnType writeMeshToLemFile(IMesh* mesh, const String& filename,const String& dir_name);
563};
564
565/*****************************************************************
566 * allowExtension
567 *****************************************************************/
569allowExtension(const String& str)
570{
571 //m_trace_mng->info() << "[allowExtension] Checking for file extension...";
572 for(vtkFileExtIdx=0;vtkFileExtReader[vtkFileExtIdx].ext!=nullptr;++vtkFileExtIdx){
573 //m_trace_mng->info() << "Testing for '" << vtkFileExtReader[vtkFileExtIdx].ext << "'...";
574 if (str == vtkFileExtReader[vtkFileExtIdx].ext){
575 return true;
576 }
577 }
578 //m_trace_mng->info()<<"Miss for all!";
579 // Sets our index in place for further service management
580 vtkFileExtIdx=0;
581 return false;
582}
583
584/*---------------------------------------------------------------------------*/
585/*---------------------------------------------------------------------------*/
586
587/*****************************************************************
588 * readMeshFromFile switch
589 *****************************************************************/
591readMeshFromFile(IPrimaryMesh* mesh,const XmlNode& mesh_node,
592 const String& file_name,const String& dir_name,
593 bool use_internal_partition)
594{
595 ARCANE_UNUSED(mesh_node);
596 info() << "[readMeshFromFile] Forwarding to vtkFileExtReader[" << vtkFileExtIdx << "].reader";
597 VtuMeshReaderBase vtu_mesh_reader{traceMng()};
598 // return (vtu_mesh_reader.*vtkFileExtReader[vtkFileExtIdx].reader)(mesh, file_name, dir_name, use_internal_partition);
599 return vtu_mesh_reader.readMeshFromVtuFile(mesh, file_name, dir_name, use_internal_partition);
600}
601
602/*---------------------------------------------------------------------------*/
603/*---------------------------------------------------------------------------*/
604
605
607
609 ServiceProperty("VtuNewMeshReader", ST_SubDomain),
611
612
613/*---------------------------------------------------------------------------*/
614/*---------------------------------------------------------------------------*/
615
616/*---------------------------------------------------------------------------*/
617/*---------------------------------------------------------------------------*/
618
619class VtuCaseMeshReader
620: public AbstractService
621, public ICaseMeshReader
622{
623public:
624
625 class Builder
626 : public IMeshBuilder
627 {
628 public:
629
630 explicit Builder(ITraceMng* tm, const CaseMeshReaderReadInfo& read_info)
631 : m_trace_mng(tm)
632 , m_read_info(read_info)
633 {}
634
635 public:
636
637 void fillMeshBuildInfo(MeshBuildInfo& build_info) override
638 {
639 ARCANE_UNUSED(build_info);
640 }
642 {
643 VtuMeshReaderBase vtu_service(pm->traceMng());
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());
647 if (ret)
648 ARCANE_FATAL("Can not read VTK File");
649 }
650
651 private:
652
653 ITraceMng* m_trace_mng;
654 CaseMeshReaderReadInfo m_read_info;
655 };
656
657public:
658
659 explicit VtuCaseMeshReader(const ServiceBuildInfo& sbi)
660 : AbstractService(sbi)
661 {}
662
663public:
664
666 {
667 IMeshBuilder* builder = nullptr;
668 if (read_info.format() == "vtu")
669 builder = new Builder(traceMng(), read_info);
670 return makeRef(builder);
671 }
672};
673
674/*---------------------------------------------------------------------------*/
675/*---------------------------------------------------------------------------*/
676
677ARCANE_REGISTER_SERVICE(VtuCaseMeshReader,
678 ServiceProperty("VtuCaseMeshReader", ST_SubDomain),
679 ARCANE_SERVICE_INTERFACE(ICaseMeshReader));
680
681/*---------------------------------------------------------------------------*/
682/*---------------------------------------------------------------------------*/
683
684
685} // End namespace Arcane
686
687/*---------------------------------------------------------------------------*/
688/*---------------------------------------------------------------------------*/
#define ARCANE_FATAL(...)
Macro envoyant une exception FatalErrorException.
Types et macros pour itérer sur les entités du maillage.
#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.
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.
Definition IMeshReader.h:32
eReturnType
Types des codes de retour d'une lecture ou écriture.
Definition IMeshReader.h:37
@ RTError
Erreur lors de l'opération.
Definition IMeshReader.h:39
@ RTOk
Opération effectuée avec succès.
Definition IMeshReader.h:38
Interface du gestionnaire d'un sous-domaine.
Definition ISubDomain.h:74
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.
Noeud d'un arbre DOM.
Definition XmlNode.h:51
#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.
Definition UtilsTypes.h:426
Int32 Integer
Type représentant un entier.
UniqueArray< Real3 > Real3UniqueArray
Tableau dynamique à une dimension de vecteurs de rang 3.
Definition UtilsTypes.h:450
ConstArrayView< ItemInternal * > ItemInternalList
Type de la liste interne des entités.
Definition ItemTypes.h:466
@ ST_SubDomain
Le service s'utilise au niveau du sous-domaine.
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.
auto makeRef(InstanceType *t) -> Ref< InstanceType >
Créé une référence sur un pointeur.
UniqueArray< Integer > IntegerUniqueArray
Tableau dynamique à une dimension d'entiers.
Definition UtilsTypes.h:434
std::int32_t Int32
Type entier signé sur 32 bits.