18928565313

Committed 22 Oct 2025 07:55PM UTC coverage: 74.771%. Remained the same

Build # 18928565313

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #807 from wdhawkins/clang-tidy-init-variables

Fixing clang-tidy init-variables warnings

Coverage Stats

132 of 179 new or added lines in 52 files covered. (73.74%)

184 existing lines in 3 files now uncovered.

18203 of 24345 relevant lines covered (74.77%)

53868061.69 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

61.22

/framework/mesh/mesh_continuum/grid_vtk_utils.cc

// SPDX-FileCopyrightText: 2024 The OpenSn Authors <https://open-sn.github.io/opensn/>
// SPDX-License-Identifier: MIT

#include "framework/mesh/mesh_continuum/grid_vtk_utils.h"
#include "framework/mesh/mesh_continuum/mesh_continuum.h"
#include "framework/runtime.h"
#include "framework/logging/log.h"
#include <vtkPoints.h>
#include <vtkPointData.h>
#include <vtkCellData.h>
#include <vtkUnstructuredGrid.h>
#include <vtkUnsignedIntArray.h>
#include <vtkAppendFilter.h>
#include <vtkXMLUnstructuredGridWriter.h>
#include <vtkXMLPUnstructuredGridWriter.h>
#include <format>

namespace opensn
{

void
UploadCellGeometryDiscontinuous(const std::shared_ptr<MeshContinuum> grid,
                                const Cell& cell,
                                int64_t& node_counter,
                                vtkNew<vtkPoints>& points,
                                vtkNew<vtkUnstructuredGrid>& ugrid)
{
  size_t num_verts = cell.vertex_ids.size();

  std::vector<vtkIdType> cell_vids(num_verts);
  for (size_t v = 0; v < num_verts; ++v)
  {
    uint64_t vgi = cell.vertex_ids[v];
    std::vector<double> d_node(3);
    d_node[0] = grid->vertices[vgi].x;
    d_node[1] = grid->vertices[vgi].y;
    d_node[2] = grid->vertices[vgi].z;

    points->InsertPoint(node_counter, d_node.data());
    cell_vids[v] = node_counter++;
  }

  if (cell.GetType() == CellType::SLAB)
  {
    ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (cell.GetType() == CellType::POLYGON)
  {
    int vtk_subtype = 0;
    switch (cell.GetSubType())
    {
      case CellType::POLYGON:
        vtk_subtype = VTK_POLYGON;
        break;
      case CellType::QUADRILATERAL:
        vtk_subtype = VTK_QUAD;
        break;
      case CellType::TRIANGLE:
        vtk_subtype = VTK_TRIANGLE;
        break;
      default:
        vtk_subtype = VTK_POLYGON;
        break;
    }

    ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (cell.GetType() == CellType::POLYHEDRON)
  {
    // Build polyhedron faces
    std::vector<vtkIdType> faces_vids;

    size_t num_faces = cell.faces.size();
    for (const auto& face : cell.faces)
    {
      size_t num_fverts = face.vertex_ids.size();
      std::vector<vtkIdType> face_info(num_fverts);
      for (size_t fv = 0; fv < num_fverts; ++fv)
      {
        size_t v = 0;
        for (size_t cv = 0; cv < num_verts; ++cv)
          if (cell.vertex_ids[cv] == face.vertex_ids[fv])
          {
            v = cv;
            break;
          }

        face_info[fv] = cell_vids[v];
      }

      faces_vids.push_back(static_cast<vtkIdType>(num_fverts));
      for (auto vid : face_info)
        faces_vids.push_back(vid);
    } // for f

    int vtk_subtype = 0;
    switch (cell.GetSubType())
    {
      case CellType::POLYHEDRON:
        vtk_subtype = VTK_POLYHEDRON;
        break;
      case CellType::PYRAMID:
        vtk_subtype = VTK_PYRAMID;
        break;
      case CellType::WEDGE:
        vtk_subtype = VTK_WEDGE;
        break;
      case CellType::HEXAHEDRON:
        vtk_subtype = VTK_HEXAHEDRON;
        break;
      case CellType::TETRAHEDRON:
        vtk_subtype = VTK_TETRA;
        break;
      default:
        vtk_subtype = VTK_POLYHEDRON;
        break;
    }

    ugrid->InsertNextCell(vtk_subtype,
                          static_cast<vtkIdType>(num_verts),
                          cell_vids.data(),
                          static_cast<vtkIdType>(num_faces),
                          faces_vids.data());
  } // polyhedron
}

void
UploadCellGeometryContinuous(const Cell& cell,
                             const std::vector<uint64_t>& vertex_map,
                             vtkNew<vtkUnstructuredGrid>& ugrid)
{
  size_t num_verts = cell.vertex_ids.size();

  std::vector<vtkIdType> cell_vids(num_verts);
  for (size_t v = 0; v < num_verts; ++v)
    cell_vids[v] = static_cast<vtkIdType>(vertex_map[cell.vertex_ids[v]]);

  if (cell.GetType() == CellType::SLAB)
  {
    ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (cell.GetType() == CellType::POLYGON)
  {
    int vtk_subtype = 0;
    switch (cell.GetSubType())
    {
      case CellType::POLYGON:
        vtk_subtype = VTK_POLYGON;
        break;
      case CellType::QUADRILATERAL:
        vtk_subtype = VTK_QUAD;
        break;
      case CellType::TRIANGLE:
        vtk_subtype = VTK_TRIANGLE;
        break;
      default:
        vtk_subtype = VTK_POLYGON;
        break;
    }

    ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (cell.GetType() == CellType::POLYHEDRON)
  {
    // Build polyhedron faces
    std::vector<vtkIdType> faces_vids;

    int vtk_subtype = 0;
    switch (cell.GetSubType())
    {
      case CellType::POLYHEDRON:
        vtk_subtype = VTK_POLYHEDRON;
        break;
      case CellType::PYRAMID:
        vtk_subtype = VTK_PYRAMID;
        break;
      case CellType::WEDGE:
        vtk_subtype = VTK_WEDGE;
        break;
      case CellType::HEXAHEDRON:
        vtk_subtype = VTK_HEXAHEDRON;
        break;
      case CellType::TETRAHEDRON:
        vtk_subtype = VTK_TETRA;
        break;
      default:
        vtk_subtype = VTK_POLYHEDRON;
        break;
    }

    switch (cell.GetSubType())
    {
      case CellType::POLYHEDRON:
      {
        size_t num_faces = cell.faces.size();
        for (const auto& face : cell.faces)
        {
          size_t num_fverts = face.vertex_ids.size();
          std::vector<vtkIdType> face_info(num_fverts);
          for (size_t fv = 0; fv < num_fverts; ++fv)
          {
            size_t v = 0;
            for (size_t cv = 0; cv < num_verts; ++cv)
              if (cell.vertex_ids[cv] == face.vertex_ids[fv])
              {
                v = cv;
                break;
              }

            face_info[fv] = cell_vids[v];
          }

          faces_vids.push_back(static_cast<vtkIdType>(num_fverts));
          for (auto vid : face_info)
            faces_vids.push_back(vid);
        } // for f

        ugrid->InsertNextCell(vtk_subtype,
                              static_cast<vtkIdType>(num_verts),
                              cell_vids.data(),
                              static_cast<vtkIdType>(num_faces),
                              faces_vids.data());
        break;
      }
      default:
        ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());
    }
  } // polyhedron
}

void
UploadFaceGeometry(const CellFace& cell_face,
                   const std::vector<uint64_t>& vertex_map,
                   vtkNew<vtkUnstructuredGrid>& ugrid)
{
  const size_t num_verts = cell_face.vertex_ids.size();

  std::vector<vtkIdType> cell_vids;
  cell_vids.reserve(cell_face.vertex_ids.size());
  for (uint64_t vid : cell_face.vertex_ids)
    cell_vids.push_back(static_cast<vtkIdType>(vertex_map[vid]));

  if (num_verts == 1)
  {
    ugrid->InsertNextCell(VTK_VERTEX, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (num_verts == 2)
  {
    ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
  if (num_verts >= 3)
  {
    int vtk_subtype = 0;
    switch (num_verts)
    {
      case 3:
        vtk_subtype = VTK_TRIANGLE;
        break;
      case 4:
        vtk_subtype = VTK_QUAD;
        break;
      default:
        vtk_subtype = VTK_POLYGON;
        break;
    }

    ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());
  }
}

int
FindHighestDimension(std::vector<vtkUGridPtrAndName>& ugrid_blocks)
{
  int max_dim = 0;
  for (auto& ugrid : ugrid_blocks)
  {
    const int64_t num_cells = ugrid.first->GetNumberOfCells();
    for (int64_t c = 0; c < num_cells; ++c)
    {
      auto* cell = ugrid.first->GetCell(c);
      OpenSnLogicalErrorIf(not cell, "Failed to obtain VTK-cell pointer");
      max_dim = std::max(max_dim, cell->GetCellDimension());
    }
  } // for ugrid in block

  return max_dim;
}

vtkUGridPtr
ConsolidateGridBlocks(std::vector<vtkUGridPtrAndName>& ugrid_blocks,
                      const std::string& block_id_array_name)
{
  const std::string fname = "ConsolidateGridBlocks";

  // Determine if all blocks have global-ids
  bool has_global_ids = true;
  for (auto& ugrid_name : ugrid_blocks)
  {
    auto& ugrid = ugrid_name.first;
    const bool has_cell_gids = ugrid->GetCellData()->GetGlobalIds();
    const bool has_pnts_gids = ugrid->GetPointData()->GetGlobalIds();
    const bool has_block_ids = ugrid->GetCellData()->GetArray(block_id_array_name.c_str());

    if ((not has_cell_gids) or (not has_pnts_gids))
      has_global_ids = false;

    if (not has_block_ids)
    {
      auto err = std::format("{}: Grid block {} does not have \"{}\" array.",
                             fname,
                             ugrid_name.second,
                             block_id_array_name);
      throw std::logic_error(err);
    }
  } // for grid_name pairs

  if (has_global_ids)
    log.Log() << fname << ": blocks have global-id arrays";

  // Consolidate the blocks
  auto append = vtkSmartPointer<vtkAppendFilter>::New();
  for (auto& ugrid : ugrid_blocks)
    append->AddInputData(ugrid.first);

  append->MergePointsOn();
  append->Update();

  auto consolidated_ugrid =
    vtkSmartPointer<vtkUnstructuredGrid>(vtkUnstructuredGrid::SafeDownCast(append->GetOutput()));

  log.Log0Verbose1() << "Consolidated grid num cells and points: "
                     << consolidated_ugrid->GetNumberOfCells() << " "
                     << consolidated_ugrid->GetNumberOfPoints();

  if (has_global_ids)
  {
    const vtkIdType num_points = consolidated_ugrid->GetNumberOfPoints();
    vtkIdType min_id = num_points;
    vtkIdType max_id = 0;
    for (vtkIdType p = 0; p < num_points; ++p)
    {
      auto* point_gids =
        vtkIdTypeArray::SafeDownCast(consolidated_ugrid->GetPointData()->GetGlobalIds());
      auto point_gid = point_gids->GetValue(p);

      min_id = std::min(min_id, point_gid);
      max_id = std::max(max_id, point_gid);
    }

    log.Log() << "Minimum and Maximum node-ids " << min_id << " " << max_id;
  }

  std::map<std::string, size_t> cell_type_count_map;
  const int64_t num_cells = consolidated_ugrid->GetNumberOfCells();
  for (int64_t c = 0; c < num_cells; ++c)
  {
    auto* cell = consolidated_ugrid->GetCell(c);
    OpenSnLogicalErrorIf(not cell, "Failed to obtain VTK-cell pointer");

    const auto* cell_type_name = cell->GetClassName();
    cell_type_count_map[cell_type_name] += 1;
  }

  if (log.GetVerbosity() >= 1)
  {
    std::stringstream outstr;
    /**Lambda to right pad an entry.*/
    auto RightPad = [](std::string& entry, size_t width)
    {
      const size_t pad_size = width - entry.size();
      entry.append(std::string(pad_size, ' '));
    };

    outstr << "Block statistictics:\n";
    for (const auto& [type_name, count] : cell_type_count_map)
    {
      auto temp_name = type_name;
      RightPad(temp_name, 20);
      outstr << "  " << temp_name << " " << count << "\n";
    }
    log.Log0Verbose1() << outstr.str();
  }

  return consolidated_ugrid;
}

std::vector<vtkUGridPtrAndName>
GetBlocksOfDesiredDimension(std::vector<vtkUGridPtrAndName>& ugrid_blocks, int desired_dimension)
{
  std::vector<vtkUGridPtrAndName> desired_blocks;
  for (auto& ugrid : ugrid_blocks)
  {
    if (ugrid.first->GetNumberOfCells() == 0)
      continue;

    std::vector<vtkUGridPtrAndName> single_grid = {ugrid};
    int block_dimension = FindHighestDimension(single_grid);

    if (block_dimension == desired_dimension)
      desired_blocks.push_back(ugrid);
  }

  return desired_blocks;
}

std::vector<uint64_t>
BuildBlockCellExtents(std::vector<vtkUGridPtrAndName>& ugrid_blocks, const int desired_dimension)
{
  std::vector<uint64_t> block_mat_ids;
  size_t total_cells = 0;

  for (auto& ugrid : ugrid_blocks)
  {
    uint64_t num_cells = ugrid.first->GetNumberOfCells();

    if (num_cells == 0)
      continue;

    if (ugrid.first->GetCell(0)->GetCellDimension() == desired_dimension)
    {
      total_cells += num_cells;
      block_mat_ids.push_back(total_cells);
    }
  }
  return block_mat_ids;
}

void
SetBlockIDArrays(std::vector<vtkUGridPtrAndName>& ugrid_blocks)
{
  int block_id = 0;
  for (auto& ugrid : ugrid_blocks)
  {
    const vtkIdType num_cells = ugrid.first->GetNumberOfCells();

    if (num_cells == 0)
      continue;

    vtkNew<vtkIntArray> block_id_list;
    block_id_list->SetName("BlockID");

    for (vtkIdType c = 0; c < num_cells; ++c)
      block_id_list->InsertNextValue(block_id);

    auto* arr = ugrid.first->GetCellData()->GetArray("BlockID");
    if (not arr)
      ugrid.first->GetCellData()->RemoveArray("BlockID");

    ugrid.first->GetCellData()->AddArray(block_id_list);
    ++block_id;
  }
}

std::vector<int>
BuildCellBlockIDsFromField(vtkUGridPtr& ugrid,
                           const std::string& field_name,
                           const std::string& file_name)
{
  const auto total_cell_count = ugrid->GetNumberOfCells();
  std::vector<int> block_ids(total_cell_count, -1);

  // Determine if reading cell identifiers
  vtkDataArray* cell_id_array_ptr = nullptr;
  if (field_name.empty())
  {
    log.Log0Warning() << "A user-supplied field name from which to recover block identifiers has "
                         "not been found. Block-ids will be left unassigned.";
    return block_ids;
  }
  else
  {
    auto* cell_data = ugrid->GetCellData();
    auto* vtk_abstract_array_ptr = cell_data->GetAbstractArray(field_name.c_str());

    if (not vtk_abstract_array_ptr)
    {
      log.Log0Warning() << "The VTU file : \"" << file_name << "\" "
                        << "does not contain a vtkCellData field of name : \"" << field_name
                        << "\". Block-ids will be left unassigned.";
      return block_ids;
    }

    cell_id_array_ptr = vtkArrayDownCast<vtkDataArray>(vtk_abstract_array_ptr);
    if (not cell_id_array_ptr)
    {
      log.Log0Warning() << "The VTU file : \"" << file_name << "\" "
                        << "with vtkCellData field of name : \"" << field_name << "\" "
                        << "cannot be downcast to vtkDataArray. Block-ids will be left "
                           "unassigned.";
      return block_ids;
    }

    const auto cell_id_n_tup = cell_id_array_ptr->GetNumberOfTuples();
    if (cell_id_n_tup != total_cell_count)
    {
      log.Log0Warning() << "The VTU file : \"" << file_name << "\" "
                        << "with vtkCellData field of name : \"" << field_name
                        << "\" has n. tuples : " << cell_id_n_tup
                        << ", but differs from the value expected : " << total_cell_count
                        << ". Block-ids will be left unassigned.";
      return block_ids;
    }

    const auto cell_id_n_val = cell_id_array_ptr->GetNumberOfValues();
    if (cell_id_n_val != total_cell_count)
    {
      log.Log0Warning() << "The VTU file : \"" << file_name << "\" "
                        << "with vtkCellData field of name : \"" << field_name
                        << "\" has n. values : " << cell_id_n_val
                        << ", but differs from the value expected : " << total_cell_count
                        << ". Block-ids will be left unassigned.";
      return block_ids;
    }
  }

  //  apply cell identifier
  for (vtkIdType c = 0; c < total_cell_count; ++c)
  {
    std::vector<double> cell_id_vec(1);
    cell_id_array_ptr->GetTuple(c, cell_id_vec.data());
    const auto mat_id = (int)cell_id_vec.front();

    block_ids[c] = mat_id;
  }

  return block_ids;
}

vtkNew<vtkUnstructuredGrid>
PrepareVtkUnstructuredGrid(const std::shared_ptr<MeshContinuum> grid, bool discontinuous)
{
  // Instantiate VTK items
  vtkNew<vtkUnstructuredGrid> ugrid;
  vtkNew<vtkPoints> points;
  vtkNew<vtkIntArray> block_array;
  vtkNew<vtkIntArray> partition_id_array;

  points->SetDataType(VTK_DOUBLE);

  // Set names
  block_array->SetName("Block");
  partition_id_array->SetName("Partition");

  std::vector<uint64_t> vertex_map;
  if (not discontinuous)
  {
    vertex_map.assign(grid->GetGlobalVertexCount(), 0);
    const size_t num_verts = grid->GetGlobalVertexCount();
    for (size_t v = 0; v < num_verts; ++v)
      vertex_map[v] = v;
  }

  // Populate cell information
  int64_t node_count = 0;
  for (const auto& cell : grid->local_cells)
  {
    if (discontinuous)
      UploadCellGeometryDiscontinuous(grid, cell, node_count, points, ugrid);
    else
    {
      for (uint64_t vid : cell.vertex_ids)
      {
        const auto& vertex = grid->vertices[vid];
        points->InsertNextPoint(vertex.x, vertex.y, vertex.z);
        vertex_map[vid] = node_count;
        ++node_count;
      }
      UploadCellGeometryContinuous(cell, vertex_map, ugrid);
    }

    block_array->InsertNextValue(cell.block_id);
    partition_id_array->InsertNextValue(cell.partition_id);
  } // for local cells
  ugrid->SetPoints(points);

  ugrid->GetCellData()->AddArray(block_array);
  ugrid->GetCellData()->AddArray(partition_id_array);

  return ugrid;
}

void
WritePVTUFiles(vtkNew<vtkUnstructuredGrid>& ugrid, const std::string& file_base_name)
{
  // Construct file name
  std::string base_filename = std::string(file_base_name);
  std::string location_filename = base_filename + std::string("_") +
                                  std::to_string(opensn::mpi_comm.rank()) + std::string(".vtu");

  // Write master file
  if (opensn::mpi_comm.rank() == 0)
  {
    std::string pvtu_file_name = base_filename + std::string(".pvtu");

    auto pgrid_writer = vtkSmartPointer<vtkXMLPUnstructuredGridWriter>::New();

    pgrid_writer->EncodeAppendedDataOff();
    pgrid_writer->SetFileName(pvtu_file_name.c_str());
    pgrid_writer->SetNumberOfPieces(opensn::mpi_comm.size());
    pgrid_writer->SetStartPiece(opensn::mpi_comm.rank());
    pgrid_writer->SetEndPiece(opensn::mpi_comm.size() - 1);
    pgrid_writer->SetInputData(ugrid);

    pgrid_writer->Write();
  }
  opensn::mpi_comm.barrier();

  // Serial output each piece
  auto grid_writer = vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();

  grid_writer->SetInputData(ugrid);
  grid_writer->SetFileName(location_filename.c_str());

  grid_writer->Write();
}

} // namespace opensn

1	// SPDX-FileCopyrightText: 2024 The OpenSn Authors <https://open-sn.github.io/opensn/>
2	// SPDX-License-Identifier: MIT
3
4	#include "framework/mesh/mesh_continuum/grid_vtk_utils.h"
5	#include "framework/mesh/mesh_continuum/mesh_continuum.h"
6	#include "framework/runtime.h"
7	#include "framework/logging/log.h"
8	#include <vtkPoints.h>
9	#include <vtkPointData.h>
10	#include <vtkCellData.h>
11	#include <vtkUnstructuredGrid.h>
12	#include <vtkUnsignedIntArray.h>
13	#include <vtkAppendFilter.h>
14	#include <vtkXMLUnstructuredGridWriter.h>
15	#include <vtkXMLPUnstructuredGridWriter.h>
16	#include <format>
17
18	namespace opensn
19	{
20
21	void
22	UploadCellGeometryDiscontinuous(const std::shared_ptr<MeshContinuum> grid,	40,951✔
23	const Cell& cell,
24	int64_t& node_counter,
25	vtkNew<vtkPoints>& points,
26	vtkNew<vtkUnstructuredGrid>& ugrid)
27	{
28	size_t num_verts = cell.vertex_ids.size();	40,951✔
29
30	std::vector<vtkIdType> cell_vids(num_verts);	40,951✔
31	for (size_t v = 0; v < num_verts; ++v)	198,915✔
32	{
33	uint64_t vgi = cell.vertex_ids[v];	157,964✔
34	std::vector<double> d_node(3);	157,964✔
35	d_node[0] = grid->vertices[vgi].x;	157,964✔
36	d_node[1] = grid->vertices[vgi].y;	157,964✔
37	d_node[2] = grid->vertices[vgi].z;	157,964✔
38
39	points->InsertPoint(node_counter, d_node.data());	157,964✔
40	cell_vids[v] = node_counter++;	157,964✔
41	}	157,964✔
42
43	if (cell.GetType() == CellType::SLAB)	40,951✔
44	{
45	ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());	2,920✔
46	}
47	if (cell.GetType() == CellType::POLYGON)	40,951✔
48	{
49	int vtk_subtype = 0;	38,031✔
50	switch (cell.GetSubType())	38,031✔
51	{
52	case CellType::POLYGON:
53	vtk_subtype = VTK_POLYGON;
54	break;
55	case CellType::QUADRILATERAL:	38,031✔
56	vtk_subtype = VTK_QUAD;	38,031✔
57	break;	38,031✔
58	case CellType::TRIANGLE:	×
59	vtk_subtype = VTK_TRIANGLE;	×
60	break;	×
61	default:
62	vtk_subtype = VTK_POLYGON;
63	break;
64	}
65
66	ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());	38,031✔
67	}
68	if (cell.GetType() == CellType::POLYHEDRON)	40,951✔
69	{
70	// Build polyhedron faces
71	std::vector<vtkIdType> faces_vids;	×
72
73	size_t num_faces = cell.faces.size();	×
74	for (const auto& face : cell.faces)	×
75	{
76	size_t num_fverts = face.vertex_ids.size();	×
77	std::vector<vtkIdType> face_info(num_fverts);	×
78	for (size_t fv = 0; fv < num_fverts; ++fv)	×
79	{
80	size_t v = 0;	×
81	for (size_t cv = 0; cv < num_verts; ++cv)	×
82	if (cell.vertex_ids[cv] == face.vertex_ids[fv])	×
83	{
84	v = cv;
85	break;
86	}
87
88	face_info[fv] = cell_vids[v];	×
89	}
90
91	faces_vids.push_back(static_cast<vtkIdType>(num_fverts));	×
92	for (auto vid : face_info)	×
93	faces_vids.push_back(vid);	×
94	} // for f	×
95
NEW 96	int vtk_subtype = 0;	×
97	switch (cell.GetSubType())	×
98	{
99	case CellType::POLYHEDRON:
100	vtk_subtype = VTK_POLYHEDRON;
101	break;
102	case CellType::PYRAMID:	×
103	vtk_subtype = VTK_PYRAMID;	×
104	break;	×
105	case CellType::WEDGE:	×
106	vtk_subtype = VTK_WEDGE;	×
107	break;	×
108	case CellType::HEXAHEDRON:	×
109	vtk_subtype = VTK_HEXAHEDRON;	×
110	break;	×
111	case CellType::TETRAHEDRON:	×
112	vtk_subtype = VTK_TETRA;	×
113	break;	×
114	default:
115	vtk_subtype = VTK_POLYHEDRON;
116	break;
117	}
118
119	ugrid->InsertNextCell(vtk_subtype,	×
120	static_cast<vtkIdType>(num_verts),
121	cell_vids.data(),	×
122	static_cast<vtkIdType>(num_faces),
123	faces_vids.data());	×
124	} // polyhedron	×
125	}	40,951✔
126
127	void
128	UploadCellGeometryContinuous(const Cell& cell,	248,446✔
129	const std::vector<uint64_t>& vertex_map,
130	vtkNew<vtkUnstructuredGrid>& ugrid)
131	{
132	size_t num_verts = cell.vertex_ids.size();	248,446✔
133
134	std::vector<vtkIdType> cell_vids(num_verts);	248,446✔
135	for (size_t v = 0; v < num_verts; ++v)	1,769,514✔
136	cell_vids[v] = static_cast<vtkIdType>(vertex_map[cell.vertex_ids[v]]);	1,521,068✔
137
138	if (cell.GetType() == CellType::SLAB)	248,446✔
139	{
140	ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());	×
141	}
142	if (cell.GetType() == CellType::POLYGON)	248,446✔
143	{
144	int vtk_subtype = 0;	83,600✔
145	switch (cell.GetSubType())	83,600✔
146	{
147	case CellType::POLYGON:
148	vtk_subtype = VTK_POLYGON;
149	break;
150	case CellType::QUADRILATERAL:	83,464✔
151	vtk_subtype = VTK_QUAD;	83,464✔
152	break;	83,464✔
153	case CellType::TRIANGLE:	136✔
154	vtk_subtype = VTK_TRIANGLE;	136✔
155	break;	136✔
156	default:
157	vtk_subtype = VTK_POLYGON;
158	break;
159	}
160
161	ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());	83,600✔
162	}
163	if (cell.GetType() == CellType::POLYHEDRON)	248,446✔
164	{
165	// Build polyhedron faces
166	std::vector<vtkIdType> faces_vids;	164,846✔
167
168	int vtk_subtype = 0;	164,846✔
169	switch (cell.GetSubType())	164,846✔
170	{
171	case CellType::POLYHEDRON:
172	vtk_subtype = VTK_POLYHEDRON;
173	break;
174	case CellType::PYRAMID:	×
175	vtk_subtype = VTK_PYRAMID;	×
176	break;	×
177	case CellType::WEDGE:	7,710✔
178	vtk_subtype = VTK_WEDGE;	7,710✔
179	break;	7,710✔
180	case CellType::HEXAHEDRON:	128,000✔
181	vtk_subtype = VTK_HEXAHEDRON;	128,000✔
182	break;	128,000✔
183	case CellType::TETRAHEDRON:	29,136✔
184	vtk_subtype = VTK_TETRA;	29,136✔
185	break;	29,136✔
186	default:
187	vtk_subtype = VTK_POLYHEDRON;
188	break;
189	}
190
191	switch (cell.GetSubType())	164,846✔
192	{
193	case CellType::POLYHEDRON:	×
194	{	×
195	size_t num_faces = cell.faces.size();	×
196	for (const auto& face : cell.faces)	×
197	{
198	size_t num_fverts = face.vertex_ids.size();	×
199	std::vector<vtkIdType> face_info(num_fverts);	×
200	for (size_t fv = 0; fv < num_fverts; ++fv)	×
201	{
202	size_t v = 0;	×
203	for (size_t cv = 0; cv < num_verts; ++cv)	×
204	if (cell.vertex_ids[cv] == face.vertex_ids[fv])	×
205	{
206	v = cv;
207	break;
208	}
209
210	face_info[fv] = cell_vids[v];	×
211	}
212
213	faces_vids.push_back(static_cast<vtkIdType>(num_fverts));	×
214	for (auto vid : face_info)	×
215	faces_vids.push_back(vid);	×
216	} // for f	×
217
218	ugrid->InsertNextCell(vtk_subtype,	×
219	static_cast<vtkIdType>(num_verts),
220	cell_vids.data(),	×
221	static_cast<vtkIdType>(num_faces),
222	faces_vids.data());	×
223	break;
224	}
225	default:	164,846✔
226	ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());	164,846✔
227	}
228	} // polyhedron	164,846✔
229	}	248,446✔
230
231	void
232	UploadFaceGeometry(const CellFace& cell_face,	×
233	const std::vector<uint64_t>& vertex_map,
234	vtkNew<vtkUnstructuredGrid>& ugrid)
235	{
236	const size_t num_verts = cell_face.vertex_ids.size();	×
237
238	std::vector<vtkIdType> cell_vids;	×
239	cell_vids.reserve(cell_face.vertex_ids.size());	×
240	for (uint64_t vid : cell_face.vertex_ids)	×
241	cell_vids.push_back(static_cast<vtkIdType>(vertex_map[vid]));	×
242
243	if (num_verts == 1)	×
244	{
245	ugrid->InsertNextCell(VTK_VERTEX, static_cast<vtkIdType>(num_verts), cell_vids.data());	×
246	}
247	if (num_verts == 2)	×
248	{
249	ugrid->InsertNextCell(VTK_LINE, static_cast<vtkIdType>(num_verts), cell_vids.data());	×
250	}
251	if (num_verts >= 3)	×
252	{
NEW 253	int vtk_subtype = 0;	×
254	switch (num_verts)	×
255	{
256	case 3:
257	vtk_subtype = VTK_TRIANGLE;
258	break;
259	case 4:	×
260	vtk_subtype = VTK_QUAD;	×
261	break;	×
262	default:	×
263	vtk_subtype = VTK_POLYGON;	×
264	break;	×
265	}
266
267	ugrid->InsertNextCell(vtk_subtype, static_cast<vtkIdType>(num_verts), cell_vids.data());	×
268	}
269	}	×
270
271	int
272	FindHighestDimension(std::vector<vtkUGridPtrAndName>& ugrid_blocks)	239✔
273	{
274	int max_dim = 0;	239✔
275	for (auto& ugrid : ugrid_blocks)	542✔
276	{
277	const int64_t num_cells = ugrid.first->GetNumberOfCells();	303✔
278	for (int64_t c = 0; c < num_cells; ++c)	1,492,239✔
279	{
280	auto* cell = ugrid.first->GetCell(c);	1,491,936✔
281	OpenSnLogicalErrorIf(not cell, "Failed to obtain VTK-cell pointer");	1,491,936✔
282	max_dim = std::max(max_dim, cell->GetCellDimension());	1,491,936✔
283	}
284	} // for ugrid in block
285
286	return max_dim;	239✔
287	}
288
289	vtkUGridPtr
290	ConsolidateGridBlocks(std::vector<vtkUGridPtrAndName>& ugrid_blocks,	37✔
291	const std::string& block_id_array_name)
292	{
293	const std::string fname = "ConsolidateGridBlocks";	37✔
294
295	// Determine if all blocks have global-ids
296	bool has_global_ids = true;	37✔
297	for (auto& ugrid_name : ugrid_blocks)	90✔
298	{
299	auto& ugrid = ugrid_name.first;	53✔
300	const bool has_cell_gids = ugrid->GetCellData()->GetGlobalIds();	53✔
301	const bool has_pnts_gids = ugrid->GetPointData()->GetGlobalIds();	53✔
302	const bool has_block_ids = ugrid->GetCellData()->GetArray(block_id_array_name.c_str());	53✔
303
304	if ((not has_cell_gids) or (not has_pnts_gids))	53✔
305	has_global_ids = false;	50✔
306
307	if (not has_block_ids)	53✔
308	{
309	auto err = std::format("{}: Grid block {} does not have \"{}\" array.",	×
310	fname,
311	ugrid_name.second,	×
312	block_id_array_name);	×
313	throw std::logic_error(err);	×
314	}	×
315	} // for grid_name pairs
316
317	if (has_global_ids)	37✔
318	log.Log() << fname << ": blocks have global-id arrays";	9✔
319
320	// Consolidate the blocks
321	auto append = vtkSmartPointer<vtkAppendFilter>::New();	37✔
322	for (auto& ugrid : ugrid_blocks)	90✔
323	append->AddInputData(ugrid.first);	53✔
324
325	append->MergePointsOn();	37✔
326	append->Update();	37✔
327
328	auto consolidated_ugrid =	37✔
329	vtkSmartPointer<vtkUnstructuredGrid>(vtkUnstructuredGrid::SafeDownCast(append->GetOutput()));	37✔
330
331	log.Log0Verbose1() << "Consolidated grid num cells and points: "	37✔
332	<< consolidated_ugrid->GetNumberOfCells() << " "	37✔
333	<< consolidated_ugrid->GetNumberOfPoints();	74✔
334
335	if (has_global_ids)	37✔
336	{
337	const vtkIdType num_points = consolidated_ugrid->GetNumberOfPoints();	3✔
338	vtkIdType min_id = num_points;	3✔
339	vtkIdType max_id = 0;	3✔
340	for (vtkIdType p = 0; p < num_points; ++p)	72✔
341	{
342	auto* point_gids =	69✔
343	vtkIdTypeArray::SafeDownCast(consolidated_ugrid->GetPointData()->GetGlobalIds());	69✔
344	auto point_gid = point_gids->GetValue(p);	69✔
345
346	min_id = std::min(min_id, point_gid);	69✔
347	max_id = std::max(max_id, point_gid);	69✔
348	}
349
350	log.Log() << "Minimum and Maximum node-ids " << min_id << " " << max_id;	9✔
351	}
352
353	std::map<std::string, size_t> cell_type_count_map;	37✔
354	const int64_t num_cells = consolidated_ugrid->GetNumberOfCells();	37✔
355	for (int64_t c = 0; c < num_cells; ++c)	437,349✔
356	{
357	auto* cell = consolidated_ugrid->GetCell(c);	437,312✔
358	OpenSnLogicalErrorIf(not cell, "Failed to obtain VTK-cell pointer");	437,312✔
359
360	const auto* cell_type_name = cell->GetClassName();	437,312✔
361	cell_type_count_map[cell_type_name] += 1;	437,312✔
362	}
363
364	if (log.GetVerbosity() >= 1)	37✔
365	{
366	std::stringstream outstr;	×
367	/*Lambda to right pad an entry./
368	auto RightPad = [](std::string& entry, size_t width)	×
369	{
370	const size_t pad_size = width - entry.size();	×
371	entry.append(std::string(pad_size, ' '));	×
372	};	×
373
374	outstr << "Block statistictics:\n";	×
375	for (const auto& [type_name, count] : cell_type_count_map)	×
376	{
377	auto temp_name = type_name;	×
378	RightPad(temp_name, 20);	×
379	outstr << " " << temp_name << " " << count << "\n";	×
380	}	×
381	log.Log0Verbose1() << outstr.str();	×
382	}	×
383
384	return consolidated_ugrid;	37✔
385	}	37✔
386
387	std::vector<vtkUGridPtrAndName>
388	GetBlocksOfDesiredDimension(std::vector<vtkUGridPtrAndName>& ugrid_blocks, int desired_dimension)	74✔
389	{
390	std::vector<vtkUGridPtrAndName> desired_blocks;	74✔
391	for (auto& ugrid : ugrid_blocks)	276✔
392	{
393	if (ugrid.first->GetNumberOfCells() == 0)	202✔
394	continue;	×
395
396	std::vector<vtkUGridPtrAndName> single_grid = {ugrid};	404✔
397	int block_dimension = FindHighestDimension(single_grid);	202✔
398
399	if (block_dimension == desired_dimension)	202✔
400	desired_blocks.push_back(ugrid);	101✔
401	}	202✔
402
403	return desired_blocks;	74✔
404	}	202✔
405
406	std::vector<uint64_t>
407	BuildBlockCellExtents(std::vector<vtkUGridPtrAndName>& ugrid_blocks, const int desired_dimension)	×
408	{
409	std::vector<uint64_t> block_mat_ids;	×
410	size_t total_cells = 0;	×
411
412	for (auto& ugrid : ugrid_blocks)	×
413	{
414	uint64_t num_cells = ugrid.first->GetNumberOfCells();	×
415
416	if (num_cells == 0)	×
417	continue;	×
418
419	if (ugrid.first->GetCell(0)->GetCellDimension() == desired_dimension)	×
420	{
421	total_cells += num_cells;	×
422	block_mat_ids.push_back(total_cells);	×
423	}
424	}
425	return block_mat_ids;	×
426	}	×
427
428	void
429	SetBlockIDArrays(std::vector<vtkUGridPtrAndName>& ugrid_blocks)	19✔
430	{
431	int block_id = 0;	19✔
432	for (auto& ugrid : ugrid_blocks)	54✔
433	{
434	const vtkIdType num_cells = ugrid.first->GetNumberOfCells();	35✔
435
436	if (num_cells == 0)	35✔
437	continue;	×
438
439	vtkNew<vtkIntArray> block_id_list;	35✔
440	block_id_list->SetName("BlockID");	35✔
441
442	for (vtkIdType c = 0; c < num_cells; ++c)	153,855✔
443	block_id_list->InsertNextValue(block_id);	153,820✔
444
445	auto* arr = ugrid.first->GetCellData()->GetArray("BlockID");	35✔
446	if (not arr)	35✔
447	ugrid.first->GetCellData()->RemoveArray("BlockID");	35✔
448
449	ugrid.first->GetCellData()->AddArray(block_id_list);	35✔
450	++block_id;	35✔
451	}	35✔
452	}	19✔
453
454	std::vector<int>
455	BuildCellBlockIDsFromField(vtkUGridPtr& ugrid,	18✔
456	const std::string& field_name,
457	const std::string& file_name)
458	{
459	const auto total_cell_count = ugrid->GetNumberOfCells();	18✔
460	std::vector<int> block_ids(total_cell_count, -1);	18✔
461
462	// Determine if reading cell identifiers
463	vtkDataArray* cell_id_array_ptr = nullptr;	18✔
464	if (field_name.empty())	18✔
465	{
466	log.Log0Warning() << "A user-supplied field name from which to recover block identifiers has "	×
467	"not been found. Block-ids will be left unassigned.";	×
468	return block_ids;	×
469	}
470	else
471	{
472	auto* cell_data = ugrid->GetCellData();	18✔
473	auto* vtk_abstract_array_ptr = cell_data->GetAbstractArray(field_name.c_str());	18✔
474
475	if (not vtk_abstract_array_ptr)	18✔
476	{
477	log.Log0Warning() << "The VTU file : \"" << file_name << "\" "	×
478	<< "does not contain a vtkCellData field of name : \"" << field_name	×
479	<< "\". Block-ids will be left unassigned.";	×
480	return block_ids;	18✔
481	}
482
483	cell_id_array_ptr = vtkArrayDownCast<vtkDataArray>(vtk_abstract_array_ptr);	18✔
484	if (not cell_id_array_ptr)	18✔
485	{
486	log.Log0Warning() << "The VTU file : \"" << file_name << "\" "	×
487	<< "with vtkCellData field of name : \"" << field_name << "\" "	×
488	<< "cannot be downcast to vtkDataArray. Block-ids will be left "
489	"unassigned.";	×
490	return block_ids;	×
491	}
492
493	const auto cell_id_n_tup = cell_id_array_ptr->GetNumberOfTuples();	18✔
494	if (cell_id_n_tup != total_cell_count)	18✔
495	{
496	log.Log0Warning() << "The VTU file : \"" << file_name << "\" "	×
497	<< "with vtkCellData field of name : \"" << field_name	×
498	<< "\" has n. tuples : " << cell_id_n_tup	×
499	<< ", but differs from the value expected : " << total_cell_count	×
500	<< ". Block-ids will be left unassigned.";	×
501	return block_ids;	×
502	}
503
504	const auto cell_id_n_val = cell_id_array_ptr->GetNumberOfValues();	18✔
505	if (cell_id_n_val != total_cell_count)	18✔
506	{
507	log.Log0Warning() << "The VTU file : \"" << file_name << "\" "	×
508	<< "with vtkCellData field of name : \"" << field_name	×
509	<< "\" has n. values : " << cell_id_n_val	×
510	<< ", but differs from the value expected : " << total_cell_count	×
511	<< ". Block-ids will be left unassigned.";	×
512	return block_ids;	×
513	}
514	}
515
516	// apply cell identifier
517	for (vtkIdType c = 0; c < total_cell_count; ++c)	283,510✔
518	{
519	std::vector<double> cell_id_vec(1);	283,492✔
520	cell_id_array_ptr->GetTuple(c, cell_id_vec.data());	283,492✔
521	const auto mat_id = (int)cell_id_vec.front();	283,492✔
522
523	block_ids[c] = mat_id;	283,492✔
524	}	283,492✔
525
526	return block_ids;
527	}	×
528
529	vtkNew<vtkUnstructuredGrid>
530	PrepareVtkUnstructuredGrid(const std::shared_ptr<MeshContinuum> grid, bool discontinuous)	86✔
531	{
532	// Instantiate VTK items
533	vtkNew<vtkUnstructuredGrid> ugrid;	86✔
534	vtkNew<vtkPoints> points;	86✔
535	vtkNew<vtkIntArray> block_array;	86✔
536	vtkNew<vtkIntArray> partition_id_array;	86✔
537
538	points->SetDataType(VTK_DOUBLE);	86✔
539
540	// Set names
541	block_array->SetName("Block");	86✔
542	partition_id_array->SetName("Partition");	86✔
543
544	std::vector<uint64_t> vertex_map;	86✔
545	if (not discontinuous)	86✔
546	{
547	vertex_map.assign(grid->GetGlobalVertexCount(), 0);	19✔
548	const size_t num_verts = grid->GetGlobalVertexCount();	19✔
549	for (size_t v = 0; v < num_verts; ++v)	232,633✔
550	vertex_map[v] = v;	232,614✔
551	}
552
553	// Populate cell information
554	int64_t node_count = 0;	86✔
555	for (const auto& cell : grid->local_cells)	289,483✔
556	{
557	if (discontinuous)	289,397✔
558	UploadCellGeometryDiscontinuous(grid, cell, node_count, points, ugrid);	122,853✔
559	else
560	{
561	for (uint64_t vid : cell.vertex_ids)	1,769,514✔
562	{
563	const auto& vertex = grid->vertices[vid];	1,521,068✔
564	points->InsertNextPoint(vertex.x, vertex.y, vertex.z);	1,521,068✔
565	vertex_map[vid] = node_count;	1,521,068✔
566	++node_count;	1,521,068✔
567	}
568	UploadCellGeometryContinuous(cell, vertex_map, ugrid);	248,446✔
569	}
570
571	block_array->InsertNextValue(cell.block_id);	289,397✔
572	partition_id_array->InsertNextValue(cell.partition_id);	289,397✔
573	} // for local cells
574	ugrid->SetPoints(points);	86✔
575
576	ugrid->GetCellData()->AddArray(block_array);	86✔
577	ugrid->GetCellData()->AddArray(partition_id_array);	86✔
578
579	return ugrid;	86✔
580	}	86✔
581
582	void
583	WritePVTUFiles(vtkNew<vtkUnstructuredGrid>& ugrid, const std::string& file_base_name)	86✔
584	{
585	// Construct file name
586	std::string base_filename = std::string(file_base_name);	86✔
587	std::string location_filename = base_filename + std::string("_") +	258✔
588	std::to_string(opensn::mpi_comm.rank()) + std::string(".vtu");	344✔
589
590	// Write master file
591	if (opensn::mpi_comm.rank() == 0)	86✔
592	{
593	std::string pvtu_file_name = base_filename + std::string(".pvtu");	30✔
594
595	auto pgrid_writer = vtkSmartPointer<vtkXMLPUnstructuredGridWriter>::New();	30✔
596
597	pgrid_writer->EncodeAppendedDataOff();	30✔
598	pgrid_writer->SetFileName(pvtu_file_name.c_str());	30✔
599	pgrid_writer->SetNumberOfPieces(opensn::mpi_comm.size());	30✔
600	pgrid_writer->SetStartPiece(opensn::mpi_comm.rank());	30✔
601	pgrid_writer->SetEndPiece(opensn::mpi_comm.size() - 1);	30✔
602	pgrid_writer->SetInputData(ugrid);	30✔
603
604	pgrid_writer->Write();	30✔
605	}	30✔
606	opensn::mpi_comm.barrier();	86✔
607
608	// Serial output each piece
609	auto grid_writer = vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();	86✔
610
611	grid_writer->SetInputData(ugrid);	86✔
612	grid_writer->SetFileName(location_filename.c_str());	86✔
613
614	grid_writer->Write();	86✔
615	}	86✔
616
617	} // namespace opensn

Open-Sn / opensn / 18928565313

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous