diff --git a/Code/Source/solver/VtkData.cpp b/Code/Source/solver/VtkData.cpp index 4e46b2885..a357da38a 100644 --- a/Code/Source/solver/VtkData.cpp +++ b/Code/Source/solver/VtkData.cpp @@ -629,6 +629,30 @@ Array VtkVtpData::get_connectivity() const return conn; } +/// @brief Copy an array of cell data from an polydata mesh into the given Vector. +// +void VtkVtpData::copy_cell_data(const std::string& data_name, Vector& mesh_data) +{ + auto vtk_data = vtkIntArray::SafeDownCast(impl->vtk_polydata->GetCellData()->GetArray(data_name.c_str())); + if (vtk_data == nullptr) { + return; + } + + int num_data = vtk_data->GetNumberOfTuples(); + if (num_data == 0) { + return; + } + + int num_comp = vtk_data->GetNumberOfComponents(); + + for (int i = 0; i < num_data; i++) { + auto tuple = vtk_data->GetTuple(i); + for (int j = 0; j < num_comp; j++) { + mesh_data(i) = tuple[j]; + } + } +} + /// @brief Copy an array of point data from an polydata mesh into the given Array. // void VtkVtpData::copy_point_data(const std::string& data_name, Array& mesh_data) @@ -774,6 +798,19 @@ Array VtkVtpData::get_points() const return points_array; } +bool VtkVtpData::has_cell_data(const std::string& data_name) +{ + int num_arrays = impl->vtk_polydata->GetPointData()->GetNumberOfArrays(); + + for (int i = 0; i < num_arrays; i++) { + if (!strcmp(impl->vtk_polydata->GetCellData()->GetArrayName(i), data_name.c_str())) { + return true; + } + } + + return false; +} + bool VtkVtpData::has_point_data(const std::string& data_name) { int num_arrays = impl->vtk_polydata->GetPointData()->GetNumberOfArrays(); @@ -911,6 +948,30 @@ std::vector VtkVtuData::get_point_data_names() return data_names; } +/// @brief Copy an array of cell data from an unstructured mesh into the given Vector. +// +void VtkVtuData::copy_cell_data(const std::string& data_name, Vector& mesh_data) +{ + auto vtk_data = vtkIntArray::SafeDownCast(impl->vtk_ugrid->GetCellData()->GetArray(data_name.c_str())); + if (vtk_data == nullptr) { + return; + } + + int num_data = vtk_data->GetNumberOfTuples(); + if (num_data == 0) { + return; + } + + int num_comp = vtk_data->GetNumberOfComponents(); + + for (int i = 0; i < num_data; i++) { + auto tuple = vtk_data->GetTuple(i); + for (int j = 0; j < num_comp; j++) { + mesh_data(i) = tuple[j]; + } + } +} + /// @brief Copy an array of point data from an unstructured grid into the given Array. // void VtkVtuData::copy_point_data(const std::string& data_name, Array& mesh_data) @@ -996,6 +1057,19 @@ void VtkVtuData::copy_points(Array& points) return; } +bool VtkVtuData::has_cell_data(const std::string& data_name) +{ + int num_arrays = impl->vtk_ugrid->GetPointData()->GetNumberOfArrays(); + + for (int i = 0; i < num_arrays; i++) { + if (!strcmp(impl->vtk_ugrid->GetCellData()->GetArrayName(i), data_name.c_str())) { + return true; + } + } + + return false; +} + bool VtkVtuData::has_point_data(const std::string& data_name) { int num_arrays = impl->vtk_ugrid->GetPointData()->GetNumberOfArrays(); diff --git a/Code/Source/solver/VtkData.h b/Code/Source/solver/VtkData.h index b8cee9c67..89d788eba 100644 --- a/Code/Source/solver/VtkData.h +++ b/Code/Source/solver/VtkData.h @@ -32,11 +32,20 @@ class VtkData { virtual void set_points(const Array& points) = 0; virtual void set_connectivity(const int nsd, const Array& conn, const int pid = 0) = 0; + virtual bool has_cell_data(const std::string& data_name) = 0; virtual bool has_point_data(const std::string& data_name) = 0; virtual void copy_points(Array& points) = 0; + virtual void copy_point_data(const std::string& data_name, Array& mesh_data) = 0; virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) = 0; + virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) = 0; + + virtual void copy_cell_data(const std::string& data_name, Vector& mesh_data) = 0; + + virtual Array get_point_data(const std::string& data_name) = 0; + virtual std::vector get_point_data_names() = 0; + virtual void write() = 0; static VtkData* create_reader(const std::string& file_name); @@ -64,13 +73,20 @@ class VtkVtpData : public VtkData { virtual int num_points() const override; virtual void read_file(const std::string& file_name) override; - void copy_points(Array& points) override; - void copy_point_data(const std::string& data_name, Array& mesh_data) override; - void copy_point_data(const std::string& data_name, Vector& mesh_data) override; - void copy_point_data(const std::string& data_name, Vector& mesh_data); - Array get_point_data(const std::string& data_name); - std::vector get_point_data_names(); - bool has_point_data(const std::string& data_name) override; + virtual void copy_points(Array& points) override; + + virtual void copy_point_data(const std::string& data_name, Array& mesh_data) override; + virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) override; + virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) override; + + virtual void copy_cell_data(const std::string& data_name, Vector& mesh_data) override; + + virtual Array get_point_data(const std::string& data_name) override; + virtual std::vector get_point_data_names() override; + + virtual bool has_cell_data(const std::string& data_name) override; + virtual bool has_point_data(const std::string& data_name) override; + virtual void set_connectivity(const int nsd, const Array& conn, const int pid = 0) override; virtual void set_element_data(const std::string& data_name, const Array& data) override; @@ -101,15 +117,22 @@ class VtkVtuData : public VtkData { virtual int num_points() const override; virtual void read_file(const std::string& file_name) override; - void copy_points(Array& points) override; - void copy_point_data(const std::string& data_name, Array& mesh_data) override; - void copy_point_data(const std::string& data_name, Vector& mesh_data) override; - void copy_point_data(const std::string& data_name, Vector& mesh_data); + virtual void copy_points(Array& points) override; + + virtual void copy_point_data(const std::string& data_name, Array& mesh_data) override; + virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) override; + virtual void copy_point_data(const std::string& data_name, Vector& mesh_data) override; + + virtual void copy_cell_data(const std::string& data_name, Vector& mesh_data) override; + + virtual Array get_point_data(const std::string& data_name) override; + virtual std::vector get_point_data_names() override; - Array get_point_data(const std::string& data_name); - std::vector get_point_data_names(); virtual Array get_points() const override; - bool has_point_data(const std::string& data_name) override; + + virtual bool has_cell_data(const std::string& data_name) override; + virtual bool has_point_data(const std::string& data_name) override; + virtual void set_connectivity(const int nsd, const Array& conn, const int pid = 0) override; virtual void set_element_data(const std::string& data_name, const Array& data) override; diff --git a/Code/Source/solver/read_msh.cpp b/Code/Source/solver/read_msh.cpp index 0fdb19000..b73b1e09e 100644 --- a/Code/Source/solver/read_msh.cpp +++ b/Code/Source/solver/read_msh.cpp @@ -2086,16 +2086,32 @@ void set_dmn_id_ff(Simulation* simulation, mshType& lM, const std::string& file_ } } -/// @brief Read mesh domains from a vtu/vtp file. +/// @brief Read mesh domains from a vtk vtu/vtp file. /// -/// \todo [NOTE] Not implemented. -// -void set_dmn_id_vtk(Simulation* simulation, mshType& mesh, const std::string& file_name, const std::string& kwrd) +void set_dmn_id_vtk(Simulation* simulation, mshType& lM, const std::string& file_name, const std::string& data_name) { - int btSiz = std::numeric_limits::digits; + #define n_debug_set_dmn_id_vtk + #ifdef debug_set_dmn_id_vtk + DebugMsg dmsg(__func__, simulation->com_mod.cm.idcm()); + dmsg.banner(); + dmsg << "file_name: " << file_name; + dmsg << "lM.gnEl: " << lM.gnEl; + #endif + + if (lM.eId.size() == 0) { + lM.eId.resize(lM.gnEl); + } + + Vector element_data(lM.gnEl); + vtk_xml::read_element_data(lM, file_name, data_name, element_data); + + for (int a = 0; a < lM.gnEl; a++) { + lM.eId(a) = lM.eId(a) | (1UL << element_data(a)); + } + } -/// @brief This routines associates two faces with each other and sets gN. +/// @brief Associate two faces with each other and set gN. /// /// Data set /// \code {.cpp} @@ -2232,4 +2248,4 @@ void set_projector(Simulation* simulation, utils::stackType& avNds) } } -}; \ No newline at end of file +}; diff --git a/Code/Source/solver/read_msh.h b/Code/Source/solver/read_msh.h index 600f174fd..258f826ae 100644 --- a/Code/Source/solver/read_msh.h +++ b/Code/Source/solver/read_msh.h @@ -52,7 +52,7 @@ namespace read_msh_ns { void read_msh(Simulation* simulation); void set_dmn_id_ff(Simulation* simulation, mshType& mesh, const std::string& file_name); - void set_dmn_id_vtk(Simulation* simulation, mshType& mesh, const std::string& file_name, const std::string& kwrd); + void set_dmn_id_vtk(Simulation* simulation, mshType& lM, const std::string& file_name, const std::string& kwrd); void set_projector(Simulation* simulation, utils::stackType& avNds); void set_ris_projector(Simulation* simulation); void set_uris_meshes(Simulation* simulation); diff --git a/Code/Source/solver/vtk_xml.cpp b/Code/Source/solver/vtk_xml.cpp index 8cb79d9b0..cc3d604f8 100644 --- a/Code/Source/solver/vtk_xml.cpp +++ b/Code/Source/solver/vtk_xml.cpp @@ -454,6 +454,39 @@ void read_vtp(const std::string& file_name, faceType& face) } } +/// @brief Read in element data from an integer vector array data +/// named `data_name` stored in VTK VTU file into 'element_data'. +/// +void read_element_data(const mshType& mesh, const std::string& file_name, const std::string& data_name, Vector& element_data) +{ + if (FILE *file = fopen(file_name.c_str(), "r")) { + fclose(file); + } else { + throw std::runtime_error("The VTK element data file '" + file_name + "' can't be read."); + } + + auto vtk_data = VtkData::create_reader(file_name); + if (vtk_data == nullptr) { + throw std::runtime_error("Failed to create VTK reader for file: " + file_name); + } + + int num_elems = vtk_data->num_elems(); + int num_points = vtk_data->num_points(); + + if (num_elems != mesh.gnEl) { + throw std::runtime_error("The number of elements (" + std::to_string(num_elems) + + ") in the file '" + file_name + "' is not equal to the number of elements (" + + std::to_string(mesh.gnEl) + ") for the mesh named '" + mesh.name + "'."); + } + + if (!vtk_data->has_cell_data(data_name)) { + throw std::runtime_error("No CellData DataArray named '" + data_name + + "' found in the VTK file '" + file_name + "' for the mesh named '" + mesh.name + "'."); + } + + vtk_data->copy_cell_data(data_name, element_data); +} + //---------------- // read_vtp_pdata //---------------- diff --git a/Code/Source/solver/vtk_xml.h b/Code/Source/solver/vtk_xml.h index 1dd72d659..1c0c31f88 100644 --- a/Code/Source/solver/vtk_xml.h +++ b/Code/Source/solver/vtk_xml.h @@ -22,6 +22,8 @@ void read_vtu(const std::string& file_name, mshType& mesh); void read_precomputed_solution_vtu(const std::string& file_name, const std::string& field_name, mshType& mesh); +void read_element_data(const mshType& mesh, const std::string& fName, const std::string& kwrd, Vector& tmpR); + void read_vtu_pdata(const std::string& fName, const std::string& kwrd, const int nsd, const int m, const int idx, mshType& mesh); void read_vtus(Simulation* simulation, Array& lA, Array& lY, Array& lD, const std::string& fName);