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"
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"
47#include "arcane/core/internal/IVariableMngInternal.h"
48#include "arcane/core/internal/VtkCellTypes.h"
60using namespace Arcane::VtkUtils;
101class VtkMeshIOService
119 VTK_MT_StructuredGrid,
120 VTK_MT_UnstructuredGrid
159 const String& dir_name,
bool use_internal_partition);
101class VtkMeshIOService {
…};
172 static const int BUFSIZE = 10000;
176 explicit VtkFile(std::istream* stream)
191 const String& expected_value1,
192 const String& expected_value2);
270 throw IOException(
"VtkFile::isEmptyNextLine()",
"Unexpected EndOfFile");
296 bool is_comment =
true;
299 for (
int i = 0; i < BUFSIZE &&
m_buf[i] !=
'\0'; ++i) {
300 if (!isspace(
m_buf[i]) &&
m_buf[i] !=
'#' && is_comment) {
303 if (
m_buf[i] ==
'#') {
312 for (
int i = 0; i < BUFSIZE &&
m_buf[i] !=
'\0'; ++i) {
313 if (
m_buf[i] ==
'\r') {
329 throw IOException(
"VtkFile::isEmptyNextLine()",
"Not Good");
354 throw IOException(
"VtkFile::isEmptyNextLine()",
"Unexpected EndOfFile");
369 bool is_comment =
true;
377 for (
int i = 0; i < BUFSIZE &&
m_buf[i] !=
'\0'; ++i) {
378 if (!isspace(
m_buf[i]) &&
m_buf[i] !=
'#' && is_comment) {
381 if (
m_buf[i] ==
'#') {
390 for (
int i = 0; i < BUFSIZE &&
m_buf[i] !=
'\0'; ++i) {
391 if (
m_buf[i] ==
'\r') {
399 throw IOException(
"VtkFile::getNextLine()",
"Not good");
418 (*m_stream) >> ws >> v;
423 throw IOException(
"VtkFile::getFloat()",
"Bad float");
442 (*m_stream) >> ws >> v;
447 throw IOException(
"VtkFile::getDouble()",
"Bad double");
466 (*m_stream) >> ws >> v;
486 constexpr size_t sizeofT =
sizeof(T);
489 Byte big_endian[sizeofT];
490 Byte little_endian[sizeofT];
493 m_stream->read((
char*)big_endian, sizeofT);
496 for (
size_t i = 0; i < sizeofT; i++) {
497 little_endian[sizeofT - 1 - i] = big_endian[i];
501 T* conv =
new (little_endian) T;
521 String expected_value_low = expected_value.
lower();
523 if (current_value_low != expected_value_low) {
524 String s =
"Expecting chain '" + expected_value +
"', found '" + current_value +
"'";
546 String expected_value1_low = expected_value1.
lower();
547 String expected_value2_low = expected_value2.
lower();
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 +
"'";
571 String expected_value_low = expected_value.
lower();
572 return (current_value_low == expected_value_low);
593 std::ifstream ifile(file_name.
localstr(), std::ifstream::binary);
596 error() <<
"Unable to read file '" << file_name <<
"'";
609 info() <<
"Titre du fichier VTK : " << title.
localstr();
617 vtk_file.setIsBinaryFile(
true);
620 eMeshType mesh_type = VTK_MT_Unknown;
628 std::istringstream mesh_type_line(buf);
629 std::string dataset_str;
630 std::string mesh_type_str;
632 mesh_type_line >> ws >> dataset_str >> ws >> mesh_type_str;
637 mesh_type = VTK_MT_StructuredGrid;
641 mesh_type = VTK_MT_UnstructuredGrid;
644 if (mesh_type == VTK_MT_Unknown) {
645 error() <<
"Support exists only for 'STRUCTURED_GRID' and 'UNSTRUCTURED_GRID' formats (format=" << mesh_type_str <<
"')";
649 debug() <<
"Lecture en-tête OK";
653 case VTK_MT_StructuredGrid:
657 case VTK_MT_UnstructuredGrid:
659 debug() <<
"Lecture _readUnstructuredGrid OK";
663 debug() <<
"Lecture _readFacesMesh OK";
693 const char* buf =
nullptr;
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;
704 error() <<
"Syntax error while reading grid 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;
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;
718 std::string float_str;
721 std::istringstream iline(buf);
722 std::string points_str;
724 iline >> ws >> points_str >> ws >> nb_node_read >> ws >> float_str;
726 error() <<
"Syntax error while reading grid dimensions";
730 if (nb_node_read != nb_node) {
731 error() <<
"Number of invalid nodes: expected=" << nb_node <<
" found=" << nb_node_read;
736 Int32 rank =
mesh->parallelMng()->commRank();
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;
742 if (use_internal_partition && rank != 0) {
751 const Integer nb_node_yz = nb_node_y * nb_node_z;
752 const Integer nb_node_xy = nb_node_x * nb_node_y;
754 Integer nb_cell = nb_cell_x * nb_cell_y * nb_cell_z;
761 info() <<
" NODE YZ = " << nb_node_yz;
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) {
770 Integer node_unique_id = y + (z)*nb_node_y + x * nb_node_y * nb_node_z;
772 nodes_unique_id[node_local_id] = node_unique_id;
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;
806 Int64 cell_unique_id = x + y * nb_cell_x + z * nb_cell_x * nb_cell_y;
808 cells_infos[cells_infos_index] = IT_Hexaedron8;
811 cells_infos[cells_infos_index] = cell_unique_id;
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;
832 mesh->setDimension(3);
833 mesh->allocateCells(nb_cell, cells_infos,
false);
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) {
847 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
848 coords[node_unique_id] =
Real3(nx, ny, nz);
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) {
860 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
861 coords[node_unique_id] =
Real3(nx, ny, nz);
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) {
873 Integer node_unique_id = x + y * nb_node_x + z * nb_node_xy;
874 coords[node_unique_id] =
Real3(nx, ny, nz);
880 throw IOException(
"_readStructuredGrid",
"Invalid type name");
886 nodes_coord_var[inode] = coords[node.
uniqueId().asInt32()];
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;
916 if (node_x != (nb_node_x - 1))
920 if (node_y != (nb_node_y - 1))
924 if (node_z != (nb_node_z - 1))
928 xmin_surface_lid.
add(face_local_id);
930 xmax_surface_lid.
add(face_local_id);
932 ymin_surface_lid.
add(face_local_id);
934 ymax_surface_lid.
add(face_local_id);
936 zmin_surface_lid.
add(face_local_id);
938 zmax_surface_lid.
add(face_local_id);
973 const char* func_name =
"VtkMeshIOService::_readNodesUnstructuredGrid()";
975 std::istringstream iline(buf);
976 std::string points_str;
977 std::string data_type_str;
980 iline >> ws >> points_str >> ws >> nb_node >> ws >> data_type_str;
983 throw IOException(func_name,
"Syntax error while reading number of nodes");
988 throw IOException(A_FUNCINFO, String::format(
"Invalid number of nodes: n={0}", nb_node));
990 info() <<
" Info: " << nb_node;
993 node_coords.
resize(nb_node);
996 for (
Integer i = 0; i < nb_node; ++i) {
1000 node_coords[i] =
Real3(nx, ny, nz);
1004 for (
Integer i = 0; i < nb_node; ++i) {
1008 node_coords[i] =
Real3(nx, ny, nz);
1012 for (
Integer i = 0; i < nb_node; ++i) {
1016 node_coords[i] =
Real3(nx, ny, nz);
1020 throw IOException(func_name,
"Invalid type name");
1046 ARCANE_UNUSED(
mesh);
1048 const char* func_name =
"VtkMeshIOService::_readCellsUnstructuredGrid()";
1054 std::istringstream iline(buf);
1055 std::string cells_str;
1059 iline >> ws >> cells_str >> ws >> nb_cell >> ws >> nb_cell_node;
1062 throw IOException(func_name,
"Syntax error while reading cells");
1066 if (nb_cell < 0 || nb_cell_node < 0) {
1068 String::format(
"Invalid dimensions: nb_cell={0} nb_cell_node={1}",
1069 nb_cell, nb_cell_node));
1072 cells_nb_node.
resize(nb_cell);
1073 cells_type.
resize(nb_cell);
1074 cells_connectivity.
resize(nb_cell_node);
1077 Integer connectivity_index = 0;
1078 for (
Integer i = 0; i < nb_cell; ++i) {
1080 cells_nb_node[i] = n;
1081 for (
Integer j = 0; j < n; ++j) {
1083 cells_connectivity[connectivity_index] = id;
1084 ++connectivity_index;
1094 std::istringstream iline(buf);
1095 std::string cell_types_str;
1097 iline >> ws >> cell_types_str >> ws >> nb_cell_type;
1100 throw IOException(func_name,
"Syntax error while reading cell types");
1103 vtk_file.
checkString(cell_types_str,
"CELL_TYPES");
1104 if (nb_cell_type != nb_cell) {
1106 nb_cell_type, nb_cell);
1110 for (
Integer i = 0; i < nb_cell; ++i) {
1112 Int16 it = vtkToArcaneCellType(vtk_ct, cells_nb_node[i]);
1131 ARCANE_UNUSED(
mesh);
1135 if (vtk_file.isEof())
1142 vtk_file.reReadSameLine();
1242 Int32 sid =
mesh->parallelMng()->commRank();
1247 bool need_read =
true;
1249 if (use_internal_partition)
1250 need_read = (sid == 0);
1252 std::array<Int64, 4> nb_cell_by_dimension = {};
1253 Int32 mesh_dimension = -1;
1260 debug() <<
"Lecture _readNodesUnstructuredGrid OK";
1261 nb_node = node_coords.
size();
1269 debug() <<
"Lecture _readCellsUnstructuredGrid OK";
1271 nb_cell = cells_nb_node.
size();
1272 nb_cell_node = cells_connectivity.
size();
1273 cells_local_id.
resize(nb_cell);
1276 mesh_build_info.
preAllocate(nb_cell, nb_cell_node);
1279 Int32 connectivity_index = 0;
1280 for (
Integer i = 0; i < nb_cell; ++i) {
1281 Int32 current_cell_nb_node = cells_nb_node[i];
1282 Int64 cell_unique_id = i;
1284 cells_local_id[i] = i;
1287 if (cell_dim >= 0 && cell_dim <= 3)
1288 ++nb_cell_by_dimension[cell_dim];
1290 auto cell_nodes = cells_connectivity.
subView(connectivity_index,current_cell_nb_node);
1291 mesh_build_info.
addCell(cells_type[i],cell_unique_id, cell_nodes);
1292 connectivity_index += current_cell_nb_node;
1296 Int32 nb_different_dim = 0;
1297 for (
Int32 i = 0; i < 4; ++i)
1298 if (nb_cell_by_dimension[i] != 0) {
1302 if (nb_different_dim > 1)
1303 ARCANE_FATAL(
"The mesh contains cells of different dimension. nb0={0} nb1={1} nb2={2} nb3={3}",
1304 nb_cell_by_dimension[0], nb_cell_by_dimension[1], nb_cell_by_dimension[2], nb_cell_by_dimension[3]);
1309 Integer wanted_dimension = mesh_dimension;
1311 mesh->setDimension(wanted_dimension);
1322 nodes_coord_var[inode] = node_coords[node.
uniqueId().asInt32()];
1328 debug() <<
"Lecture _readData OK";
1346 bool use_internal_partition)
1348 ARCANE_UNUSED(dir_name);
1350 std::ifstream ifile(file_name.
localstr(), std::ifstream::binary);
1352 info() <<
"No face descriptor file found '" << file_name <<
"'";
1357 const char* buf = 0;
1361 info() <<
"Titre du fichier VTK : " << title;
1365 vtk_file.setIsBinaryFile(
true);
1368 eMeshType mesh_type = VTK_MT_Unknown;
1374 std::istringstream mesh_type_line(buf);
1375 std::string dataset_str;
1376 std::string mesh_type_str;
1377 mesh_type_line >> ws >> dataset_str >> ws >> mesh_type_str;
1381 mesh_type = VTK_MT_UnstructuredGrid;
1384 if (mesh_type == VTK_MT_Unknown) {
1385 error() <<
"Face descriptor file type must be 'UNSTRUCTURED_GRID' (format=" << mesh_type_str <<
"')";
1398 bool need_read =
true;
1399 if (use_internal_partition)
1400 need_read = (sid == 0);
1417 nb_face = faces_nb_node.
size();
1421 faces_local_id.
resize(nb_face);
1464 if ((sid == 0) && vtk_file.isEof())
1468 pm->broadcast(bb, 0);
1475 created_infos_str() <<
"<?xml version='1.0' ?>\n";
1476 created_infos_str() <<
"<infos>";
1479 const char* buf = 0;
1482 bool reading_node =
false;
1483 bool reading_cell =
false;
1484 while (((buf = vtk_file.
getNextLine()) != 0) && !vtk_file.isEof()) {
1485 debug() <<
"Read line";
1486 std::istringstream iline(buf);
1487 std::string data_str;
1493 iline >> ws >> nb_item;
1494 reading_node =
false;
1495 reading_cell =
true;
1496 if (nb_item != nb_cell_kind)
1497 error() <<
"Size expected = " << nb_cell_kind <<
" found = " << nb_item;
1503 iline >> ws >> nb_item;
1504 reading_node =
true;
1505 reading_cell =
false;
1506 if (nb_item != nb_node)
1507 error() <<
"Size expected = " << nb_node <<
" found = " << nb_item;
1512 std::string name_str;
1515 iline >> ws >> name_str >> ws >> nb_fields;
1518 std::string type_str;
1519 std::string s_name_str;
1520 int nb_component = 1;
1521 bool is_group =
false;
1523 for (
Integer i = 0; i < nb_fields; i++) {
1525 std::istringstream iline(buf);
1526 iline >> ws >> s_name_str >> ws >> nb_component >> ws >> nb_item >> ws >> type_str;
1528 if (nb_item != nb_cell_kind && reading_cell && !reading_node)
1529 error() <<
"Size expected = " << nb_cell_kind <<
" found = " << nb_item;
1531 if (nb_item != nb_node && !reading_cell && reading_node)
1532 error() <<
"Size expected = " << nb_node <<
" found = " << nb_item;
1534 String name_str = s_name_str;
1537 if (cstr ==
"GROUP_") {
1540 debug() <<
"** ** ** GROUP ! name=" << new_name;
1541 name_str = new_name;
1546 error() <<
"Group type must be 'int', found=" << type_str;
1551 created_infos_str() <<
"<node-group name='" << name_str <<
"'/>";
1556 created_infos_str() <<
"<cell-group name='" << name_str <<
"'/>";
1564 error() <<
"Expecting 'float' or 'double' data type, found=" << type_str;
1569 fatal() <<
"Unable to read POINT_DATA: feature not implemented";
1572 created_infos_str() <<
"<cell-variable name='" << name_str <<
"'/>";
1575 throw IOException(
"Unable to read face variables: feature not supported");
1585 if (reading_node || reading_cell) {
1586 std::string type_str;
1587 std::string s_name_str;
1589 bool is_group =
false;
1590 int nb_component = 1;
1592 iline >> ws >> s_name_str >> ws >> type_str >> ws >> nb_component;
1593 debug() <<
"** ** ** READNAME: name=" << s_name_str <<
" type=" << type_str;
1595 String name_str = s_name_str;
1598 if (cstr ==
"GROUP_") {
1601 info() <<
"** ** ** GROUP ! name=" << new_name;
1602 name_str = new_name;
1606 error() <<
"Expecting 'SCALARS' data type, found=" << data_str;
1612 error() <<
"Group type must be 'int', found=" << type_str;
1620 created_infos_str() <<
"<node-group name='" << name_str <<
"'/>";
1625 created_infos_str() <<
"<cell-group name='" << name_str <<
"'/>";
1631 error() <<
"Expecting 'float' or 'double' data type, found=" << type_str;
1638 fatal() <<
"Unable to read POINT_DATA: feature not implemented";
1641 created_infos_str() <<
"<cell-variable name='" << name_str <<
"'/>";
1644 throw IOException(
"Unable to read face variables: feature not supported");
1651 error() <<
"Expecting value CELL_DATA or POINT_DATA, found='" << data_str <<
"'";
1657 created_infos_str() <<
"</infos>";
1658 if (use_internal_partition) {
1661 String str = created_infos_str.str();
1671 String str = String::fromUtf8(bytes);
1672 info() <<
"FOUND STR=" << bytes.
size() <<
" " << str;
1674 XmlNode doc_node = doc->documentNode();
1680 String name = xnode.attrValue(
"name");
1681 info() <<
"Building variable: " << name;
1692 String name = xnode.attrValue(
"name");
1693 info() <<
"Building group: " << name;
1703 String name = xnode.attrValue(
"name");
1704 info() <<
"Create node group: " << name;
1726 info() <<
"Building face group '" << name <<
"'"
1727 <<
" size=" << faces_lid.
size();
1729 mesh->faceFamily()->createGroup(name, faces_lid);
1747 info() <<
"Reading values for variable: " << var_name <<
" n=" << nb_cell;
1749 mesh->variableMng()->_internalApi()->addAutoDestroyVariable(var);
1751 for (
Integer i = 0; i < nb_cell; ++i) {
1756 info() <<
"Variable build finished: " << vtk_file.isEof();
1777 info() <<
"Reading group info for group: " << name;
1780 for (
Integer i = 0; i < nb_item; ++i) {
1783 ids.
add(local_id[i]);
1785 info() <<
"Building group: " << name <<
" nb_element=" << ids.
size();
1807 info() <<
"Lecture infos groupes de noeuds pour le groupe: " << name;
1810 for (
Integer i = 0; i < nb_item; ++i) {
1815 info() <<
"Création groupe: " << name <<
" nb_element=" << ids.
size();
1828class VtkLegacyMeshWriter
1829:
public BasicService
1849 void _saveGroups(
IItemFamily* family, std::ostream& ofile);
1828class VtkLegacyMeshWriter {
…};
1876 String fname = file_name;
1893 std::ofstream ofile(file_name.
localstr());
1896 throw IOException(
"VtkMeshIOService::writeMeshToFile(): Unable to open file");
1897 ofile <<
"# vtk DataFile Version 2.0\n";
1898 ofile <<
"Maillage Arcane\n";
1900 ofile <<
"DATASET UNSTRUCTURED_GRID\n";
1902 nodes_local_id_to_current.
fill(NULL_ITEM_ID);
1910 ofile <<
"POINTS " << nb_node <<
" double\n";
1914 const Node& node = *inode;
1915 nodes_local_id_to_current[node.
localId()] = node_index;
1916 Real3 xyz = coords[inode];
1917 ofile << xyz.
x <<
' ' << xyz.
y <<
' ' << xyz.
z <<
'\n';
1924 Integer nb_node_cell_kind = nb_cell_kind;
1925 ENUMERATE_ITEMWITHNODES(iitem, cell_kind_family->
allItems())
1927 nb_node_cell_kind += (*iitem).nbNode();
1929 ofile <<
"CELLS " << nb_cell_kind <<
' ' << nb_node_cell_kind <<
"\n";
1930 ENUMERATE_ITEMWITHNODES(iitem, cell_kind_family->
allItems())
1934 ofile << item_nb_node;
1935 for (NodeLocalId inode_lid : item.
nodes() ) {
1936 ofile <<
' ' << nodes_local_id_to_current[inode_lid];
1941 ofile <<
"CELL_TYPES " << nb_cell_kind <<
"\n";
1943 Int16 arcane_type = (*iitem).type();
1944 int type = arcaneToVtkCellType(arcane_type);
1945 ofile << type <<
'\n';
1951 ofile <<
"POINT_DATA " << nb_node <<
"\n";
1956 ofile <<
"CELL_DATA " << nb_cell_kind <<
"\n";
1957 _saveGroups(
mesh->itemFamily(cell_kind), ofile);
1963void VtkLegacyMeshWriter::
1964_saveGroups(
IItemFamily* family, std::ostream& ofile)
1966 info() <<
"Saving groups for family name=" << family->
name();
1974 if (group.
name() ==
"OuterFaces")
1976 ofile <<
"SCALARS GROUP_" << group.
name() <<
" int 1\n";
1977 ofile <<
"LOOKUP_TABLE default\n";
1978 in_group_list.fill(
'0');
1980 in_group_list[(*iitem).localId()] =
'1';
1983 ofile << in_group_list[(*iitem).localId()] <<
'\n';
2008 const String& dir_name,
bool use_internal_partition)
override
2011 ARCANE_UNUSED(mesh_node);
2013 bool ret = vtk_service.
readMesh(
mesh, file_name, dir_name, use_internal_partition);
1994class VtkMeshReader {
…};
2036class VtkLegacyCaseMeshReader
2049 , m_read_info(read_info)
2056 ARCANE_UNUSED(build_info);
2061 String fname = m_read_info.fileName();
2062 m_trace_mng->info() <<
"VtkLegacy Reader (ICaseMeshReader) file_name=" << fname;
2063 bool ret = vtk_service.
readMesh(pm, fname, m_read_info.directoryName(), m_read_info.isParallelRead());
2085 if (read_info.format() ==
"vtk")
2036class VtkLegacyCaseMeshReader {
…};
#define ARCANE_THROW(exception_class,...)
Macro pour envoyer une exception avec formattage.
#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.
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
Vue constante d'un tableau de type T.
constexpr Integer size() const noexcept
Nombre d'éléments du tableau.
Informations sur le type flottant.
Interface du service de lecture du maillage à partir du jeu de données.
Interface d'une famille d'entités.
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.
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 d'une classe proposant des fonctions utilitaires sur maillage.
virtual void localIdsFromConnectivity(eItemKind item_kind, IntegerConstArrayView items_nb_node, Int64ConstArrayView items_connectivity, Int32ArrayView local_ids, bool allow_null=false)=0
Recherche les identifiants locaux des entités à partir de leur connectivité.
Interface d'un service d'écriture d'un maillage.
Exception lorsqu'une erreur d'entrée/sortie est détectée.
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.
Groupe d'entités de maillage.
const String & name() const
Nom du groupe.
Type d'une entité (Item).
Int16 dimension() const
Dimension de l'élément (<0 si inconnu)
Gestionnaire des types d'entités de maillage.
ItemTypeInfo * typeFromId(Integer id) const
Type correspondant au numéro id.
Elément de maillage s'appuyant sur des noeuds (Edge,Face,Cell).
NodeConnectedListViewType nodes() const
Liste des noeuds de l'entité
Int32 nbNode() const
Nombre de noeuds de l'entité
constexpr Int32 localId() const
Identifiant local de l'entité dans le sous-domaine du processeur.
ItemUniqueId uniqueId() const
Identifiant unique sur tous les domaines.
Paramètres nécessaires à la construction d'un maillage.
Flot de sortie lié à une String.
Classe gérant un vecteur de réel de dimension 3.
Référence à une instance.
Encapsulation d'un pointeur qui se détruit automatiquement.
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.
const char * localstr() const
Retourne la conversion de l'instance dans l'encodage UTF-8.
ByteConstArrayView utf8() const
Retourne la conversion de l'instance dans l'encodage UTF-8.
String substring(Int64 pos) const
Sous-chaîne commençant à la position pos.
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 ¤t_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 ¤t_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.
XmlNode documentElement() const
Retourne le noeud élément du document.
XmlNodeList children(const String &name) const
Ensemble des noeuds fils de ce noeud ayant pour nom name.
#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.
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.
@ ST_SubDomain
Le service s'utilise au niveau du sous-domaine.
UniqueArray< Byte > ByteUniqueArray
Tableau dynamique à une dimension de caractères.
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.
@ 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.
unsigned char Byte
Type d'un octet.
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.
std::int32_t Int32
Type entier signé sur 32 bits.
Real y
deuxième composante du triplet
Real z
troisième composante du triplet
Real x
première composante du triplet