21773869314

Committed 06 Feb 2026 07:55PM UTC coverage: 74.829% (+0.3%) from 74.553%

Build # 21773869314

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #919 from wdhawkins/rz

Numerous fixes to RZ solver

Run Details

190 of 230 new or added lines in 8 files covered. (82.61%)

2 existing lines in 1 file now uncovered.

19579 of 26165 relevant lines covered (74.83%)

69416660.18 hits per line

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

79.2

/framework/math/spatial_discretization/spatial_discretization.cc

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

#include "framework/math/spatial_discretization/spatial_discretization.h"
#include "framework/math/petsc_utils/petsc_utils.h"
#include "framework/logging/log.h"

namespace opensn
{

SpatialDiscretization::SpatialDiscretization(const std::shared_ptr<MeshContinuum> grid,
                                             SpatialDiscretizationType sdm_type)
  : UNITARY_UNKNOWN_MANAGER({std::make_pair(UnknownType::SCALAR, 0)}), grid_(grid), type_(sdm_type)
{
}

SpatialDiscretizationType
SpatialDiscretization::GetType() const
{
  return type_;
}

std::shared_ptr<MeshContinuum>
SpatialDiscretization::GetGrid() const
{
  return grid_;
}

size_t
SpatialDiscretization::GetNumLocalNodes() const
{
  return local_base_block_size_;
}

size_t
SpatialDiscretization::GetNumLocalDOFs(const UnknownManager& unknown_manager) const
{
  unsigned int N = unknown_manager.GetTotalUnknownStructureSize();

  return local_base_block_size_ * N;
}

size_t
SpatialDiscretization::GetNumGlobalNodes() const
{
  return global_base_block_size_;
}

size_t
SpatialDiscretization::GetNumGlobalDOFs(const UnknownManager& unknown_manager) const
{
  unsigned int N = unknown_manager.GetTotalUnknownStructureSize();

  return global_base_block_size_ * N;
}

size_t
SpatialDiscretization::GetNumLocalAndGhostDOFs(const UnknownManager& unknown_manager) const
{
  return GetNumLocalDOFs(unknown_manager) + GetNumGhostDOFs(unknown_manager);
}

size_t
SpatialDiscretization::GetCellNumNodes(const Cell& cell) const
{
  return GetCellMapping(cell).GetNumNodes();
}

const std::vector<Vector3>&
SpatialDiscretization::GetCellNodeLocations(const Cell& cell) const
{
  return GetCellMapping(cell).GetNodeLocations();
}

std::pair<std::set<uint32_t>, std::set<uint32_t>>
SpatialDiscretization::MakeCellInternalAndBndryNodeIDs(const Cell& cell) const
{
  const auto& cell_mapping = GetCellMapping(cell);
  const size_t num_faces = cell.faces.size();
  const size_t num_nodes = cell_mapping.GetNumNodes();

  // Determine which nodes are on the boundary
  std::set<uint32_t> boundary_nodes;
  for (size_t f = 0; f < num_faces; ++f)
  {
    if (not cell.faces[f].has_neighbor)
    {
      const size_t num_face_nodes = cell_mapping.GetNumFaceNodes(f);
      for (size_t fi = 0; fi < num_face_nodes; ++fi)
        boundary_nodes.insert(cell_mapping.MapFaceNode(f, fi));
    }
  } // for f

  // Determine non-boundary nodes
  std::set<uint32_t> internal_nodes;
  for (size_t i = 0; i < num_nodes; ++i)
    if (boundary_nodes.find(i) == boundary_nodes.end())
      internal_nodes.insert(i);

  return {internal_nodes, boundary_nodes};
}

std::vector<std::vector<std::vector<int>>>
SpatialDiscretization::MakeInternalFaceNodeMappings(const double tolerance) const
{
  std::vector<std::vector<std::vector<int>>> cell_adj_mapping;
  for (const auto& cell : grid_->local_cells)
  {
    const auto& cell_mapping = this->GetCellMapping(cell);
    const auto& node_locations = cell_mapping.GetNodeLocations();
    const size_t num_faces = cell.faces.size();

    std::vector<std::vector<int>> per_face_adj_mapping;

    for (size_t f = 0; f < num_faces; ++f)
    {
      const auto& face = cell.faces[f];
      const auto num_face_nodes = cell_mapping.GetNumFaceNodes(f);
      std::vector<int> face_adj_mapping(num_face_nodes, -1);
      if (face.has_neighbor)
      {
        const auto& adj_cell = grid_->cells[face.neighbor_id];
        const auto& adj_cell_mapping = this->GetCellMapping(adj_cell);
        const auto& adj_node_locations = adj_cell_mapping.GetNodeLocations();
        const size_t adj_num_nodes = adj_cell_mapping.GetNumNodes();

        for (size_t fi = 0; fi < num_face_nodes; ++fi)
        {
          const int i = cell_mapping.MapFaceNode(f, fi);
          const auto& ivec3 = node_locations[i];

          for (size_t ai = 0; ai < adj_num_nodes; ++ai)
          {
            const auto& aivec3 = adj_node_locations[ai];
            if ((ivec3 - aivec3).NormSquare() < tolerance)
            {
              face_adj_mapping[fi] = static_cast<int>(ai);
              break;
            }
          } // for ai
          if (face_adj_mapping[fi] < 0)
            throw std::logic_error("Face node mapping failed");
        } // for fi
      } // if internal face

      per_face_adj_mapping.push_back(std::move(face_adj_mapping));
    } // for face

    cell_adj_mapping.push_back(std::move(per_face_adj_mapping));
  } // for cell

  return cell_adj_mapping;
}

void
SpatialDiscretization::CopyVectorWithUnknownScope(const std::vector<double>& from_vector,
                                                  std::vector<double>& to_vector,
                                                  const UnknownManager& from_vec_uk_structure,
                                                  const unsigned int from_vec_uk_id,
                                                  const UnknownManager& to_vec_uk_structure,
                                                  const unsigned int to_vec_uk_id) const
{
  const std::string fname = "spatial_discretization::"
                            "CopyVectorWithUnknownScope";
  const auto& ukmanF = from_vec_uk_structure;
  const auto& ukmanT = to_vec_uk_structure;
  const auto& ukidF = from_vec_uk_id;
  const auto& ukidT = to_vec_uk_id;
  try
  {
    const auto& ukA = from_vec_uk_structure.unknowns.at(from_vec_uk_id);
    const auto& ukB = to_vec_uk_structure.unknowns.at(to_vec_uk_id);

    if (ukA.num_components != ukB.num_components)
      throw std::logic_error(fname + " Unknowns do not have the "
                                     "same number of components");

    const size_t num_comps = ukA.num_components;

    for (const auto& cell : grid_->local_cells)
    {
      const auto& cell_mapping = this->GetCellMapping(cell);
      const size_t num_nodes = cell_mapping.GetNumNodes();

      for (size_t i = 0; i < num_nodes; ++i)
      {
        for (size_t c = 0; c < num_comps; ++c)
        {
          const auto fmap = MapDOFLocal(cell, i, ukmanF, ukidF, c);
          const auto imap = MapDOFLocal(cell, i, ukmanT, ukidT, c);

          to_vector[imap] = from_vector[fmap];
        } // for component c
      } // for node i
    } // for cell
  }
  catch (const std::out_of_range& oor)
  {
    throw std::out_of_range(fname + ": either from_vec_uk_id or to_vec_uk_id is "
                                    "out of range for its respective "
                                    "unknown manager.");
  }
}

void
SpatialDiscretization::LocalizePETScVector(Vec petsc_vector,
                                           std::vector<double>& local_vector,
                                           const UnknownManager& unknown_manager) const
{
  size_t num_local_dofs = GetNumLocalDOFs(unknown_manager);

  CopyVecToSTLvector(petsc_vector, local_vector, num_local_dofs);
}

void
SpatialDiscretization::LocalizePETScVectorWithGhosts(Vec petsc_vector,
                                                     std::vector<double>& local_vector,
                                                     const UnknownManager& unknown_manager) const
{
  size_t num_local_dofs = GetNumLocalAndGhostDOFs(unknown_manager);

  CopyVecToSTLvectorWithGhosts(petsc_vector, local_vector, num_local_dofs);
}

double
SpatialDiscretization::CartesianSpatialWeightFunction(const Vector3& point)
{
  return 1.0;
}

double
SpatialDiscretization::CylindricalRZSpatialWeightFunction(const Vector3& point)
{
  return 2.0 * M_PI * point[0];
}

double
SpatialDiscretization::Spherical1DSpatialWeightFunction(const Vector3& point)
{
  const double r = point[2];
  return 4.0 * M_PI * r * r;
}

SpatialDiscretization::SpatialWeightFunction
SpatialDiscretization::GetSpatialWeightingFunction() const
{
  switch (grid_->GetCoordinateSystem())
  {
    case CoordinateSystemType::CARTESIAN:
      return CartesianSpatialWeightFunction;
    case CoordinateSystemType::CYLINDRICAL:
      return CylindricalRZSpatialWeightFunction;
    case CoordinateSystemType::SPHERICAL:
      return Spherical1DSpatialWeightFunction;
    case CoordinateSystemType::UNDEFINED:
    default:
      OpenSnLogicalError("Coordinate system undefined.");
  }
}

} // namespace opensn

1	// SPDX-FileCopyrightText: 2024 The OpenSn Authors <https://open-sn.github.io/opensn/>
2	// SPDX-License-Identifier: MIT
3
4	#include "framework/math/spatial_discretization/spatial_discretization.h"
5	#include "framework/math/petsc_utils/petsc_utils.h"
6	#include "framework/logging/log.h"
7
8	namespace opensn
9	{
10
11	SpatialDiscretization::SpatialDiscretization(const std::shared_ptr<MeshContinuum> grid,	715✔
12	SpatialDiscretizationType sdm_type)	715✔
13	: UNITARY_UNKNOWN_MANAGER({std::make_pair(UnknownType::SCALAR, 0)}), grid_(grid), type_(sdm_type)	715✔
14	{
15	}	715✔
16
17	SpatialDiscretizationType
18	SpatialDiscretization::GetType() const	237✔
19	{
20	return type_;	237✔
21	}
22
23	std::shared_ptr<MeshContinuum>
24	SpatialDiscretization::GetGrid() const	2,147,483,647✔
25	{
26	return grid_;	2,147,483,647✔
27	}
28
29	size_t
30	SpatialDiscretization::GetNumLocalNodes() const	644✔
31	{
32	return local_base_block_size_;	644✔
33	}
34
35	size_t
36	SpatialDiscretization::GetNumLocalDOFs(const UnknownManager& unknown_manager) const	114,323✔
37	{
38	unsigned int N = unknown_manager.GetTotalUnknownStructureSize();	114,323✔
39
40	return local_base_block_size_ * N;	114,323✔
41	}
42
43	size_t
44	SpatialDiscretization::GetNumGlobalNodes() const	564✔
45	{
46	return global_base_block_size_;	564✔
47	}
48
49	size_t
50	SpatialDiscretization::GetNumGlobalDOFs(const UnknownManager& unknown_manager) const	77,799✔
51	{
52	unsigned int N = unknown_manager.GetTotalUnknownStructureSize();	77,799✔
53
54	return global_base_block_size_ * N;	77,799✔
55	}
56
57	size_t
58	SpatialDiscretization::GetNumLocalAndGhostDOFs(const UnknownManager& unknown_manager) const	133✔
59	{
60	return GetNumLocalDOFs(unknown_manager) + GetNumGhostDOFs(unknown_manager);	133✔
61	}
62
63	size_t
64	SpatialDiscretization::GetCellNumNodes(const Cell& cell) const	7,286,752✔
65	{
66	return GetCellMapping(cell).GetNumNodes();	7,286,752✔
67	}
68
69	const std::vector<Vector3>&
70	SpatialDiscretization::GetCellNodeLocations(const Cell& cell) const	5,016,962✔
71	{
72	return GetCellMapping(cell).GetNodeLocations();	5,016,962✔
73	}
74
75	std::pair<std::set<uint32_t>, std::set<uint32_t>>
76	SpatialDiscretization::MakeCellInternalAndBndryNodeIDs(const Cell& cell) const	141,912✔
77	{
78	const auto& cell_mapping = GetCellMapping(cell);	141,912✔
79	const size_t num_faces = cell.faces.size();	141,912✔
80	const size_t num_nodes = cell_mapping.GetNumNodes();	141,912✔
81
82	// Determine which nodes are on the boundary
83	std::set<uint32_t> boundary_nodes;	141,912✔
84	for (size_t f = 0; f < num_faces; ++f)	892,896✔
85	{
86	if (not cell.faces[f].has_neighbor)	750,984✔
87	{
88	const size_t num_face_nodes = cell_mapping.GetNumFaceNodes(f);	21,540✔
89	for (size_t fi = 0; fi < num_face_nodes; ++fi)	98,044✔
90	boundary_nodes.insert(cell_mapping.MapFaceNode(f, fi));	76,504✔
91	}
92	} // for f
93
94	// Determine non-boundary nodes
95	std::set<uint32_t> internal_nodes;	141,912✔
96	for (size_t i = 0; i < num_nodes; ++i)	1,077,872✔
97	if (boundary_nodes.find(i) == boundary_nodes.end())	935,960✔
98	internal_nodes.insert(i);	862,022✔
99
100	return {internal_nodes, boundary_nodes};	141,912✔
101	}	141,912✔
102
103	std::vector<std::vector<std::vector<int>>>
104	SpatialDiscretization::MakeInternalFaceNodeMappings(const double tolerance) const	1✔
105	{
106	std::vector<std::vector<std::vector<int>>> cell_adj_mapping;	1✔
107	for (const auto& cell : grid_->local_cells)	626✔
108	{
109	const auto& cell_mapping = this->GetCellMapping(cell);	625✔
110	const auto& node_locations = cell_mapping.GetNodeLocations();	625✔
111	const size_t num_faces = cell.faces.size();	625✔
112
113	std::vector<std::vector<int>> per_face_adj_mapping;	625✔
114
115	for (size_t f = 0; f < num_faces; ++f)	3,125✔
116	{
117	const auto& face = cell.faces[f];	2,500✔
118	const auto num_face_nodes = cell_mapping.GetNumFaceNodes(f);	2,500✔
119	std::vector<int> face_adj_mapping(num_face_nodes, -1);	2,500✔
120	if (face.has_neighbor)	2,500✔
121	{
122	const auto& adj_cell = grid_->cells[face.neighbor_id];	2,400✔
123	const auto& adj_cell_mapping = this->GetCellMapping(adj_cell);	2,400✔
124	const auto& adj_node_locations = adj_cell_mapping.GetNodeLocations();	2,400✔
125	const size_t adj_num_nodes = adj_cell_mapping.GetNumNodes();	2,400✔
126
127	for (size_t fi = 0; fi < num_face_nodes; ++fi)	7,200✔
128	{
129	const int i = cell_mapping.MapFaceNode(f, fi);	4,800✔
130	const auto& ivec3 = node_locations[i];	4,800✔
131
132	for (size_t ai = 0; ai < adj_num_nodes; ++ai)	12,000✔
133	{
134	const auto& aivec3 = adj_node_locations[ai];	12,000✔
135	if ((ivec3 - aivec3).NormSquare() < tolerance)	12,000✔
136	{
137	face_adj_mapping[fi] = static_cast<int>(ai);	4,800✔
138	break;	4,800✔
139	}
140	} // for ai
141	if (face_adj_mapping[fi] < 0)	4,800✔
142	throw std::logic_error("Face node mapping failed");	×
143	} // for fi
144	} // if internal face
145
146	per_face_adj_mapping.push_back(std::move(face_adj_mapping));	2,500✔
147	} // for face	2,500✔
148
149	cell_adj_mapping.push_back(std::move(per_face_adj_mapping));	625✔
150	} // for cell	625✔
151
152	return cell_adj_mapping;	1✔
153	}	×
154
155	void
156	SpatialDiscretization::CopyVectorWithUnknownScope(const std::vector<double>& from_vector,	6✔
157	std::vector<double>& to_vector,
158	const UnknownManager& from_vec_uk_structure,
159	const unsigned int from_vec_uk_id,
160	const UnknownManager& to_vec_uk_structure,
161	const unsigned int to_vec_uk_id) const
162	{
163	const std::string fname = "spatial_discretization::"	6✔
164	"CopyVectorWithUnknownScope";	6✔
165	const auto& ukmanF = from_vec_uk_structure;	6✔
166	const auto& ukmanT = to_vec_uk_structure;	6✔
167	const auto& ukidF = from_vec_uk_id;	6✔
168	const auto& ukidT = to_vec_uk_id;	6✔
169	try	6✔
170	{
171	const auto& ukA = from_vec_uk_structure.unknowns.at(from_vec_uk_id);	6✔
172	const auto& ukB = to_vec_uk_structure.unknowns.at(to_vec_uk_id);	6✔
173
174	if (ukA.num_components != ukB.num_components)	6✔
175	throw std::logic_error(fname + " Unknowns do not have the "	×
176	"same number of components");	×
177
178	const size_t num_comps = ukA.num_components;	6✔
179
180	for (const auto& cell : grid_->local_cells)	2,727✔
181	{
182	const auto& cell_mapping = this->GetCellMapping(cell);	2,721✔
183	const size_t num_nodes = cell_mapping.GetNumNodes();	2,721✔
184
185	for (size_t i = 0; i < num_nodes; ++i)	13,305✔
186	{
187	for (size_t c = 0; c < num_comps; ++c)	213,068✔
188	{
189	const auto fmap = MapDOFLocal(cell, i, ukmanF, ukidF, c);	202,484✔
190	const auto imap = MapDOFLocal(cell, i, ukmanT, ukidT, c);	202,484✔
191
192	to_vector[imap] = from_vector[fmap];	202,484✔
193	} // for component c
194	} // for node i
195	} // for cell
196	}
197	catch (const std::out_of_range& oor)	×
198	{
199	throw std::out_of_range(fname + ": either from_vec_uk_id or to_vec_uk_id is "	×
200	"out of range for its respective "
201	"unknown manager.");	×
202	}	×
203	}	6✔
204
205	void
206	SpatialDiscretization::LocalizePETScVector(Vec petsc_vector,	9,037✔
207	std::vector<double>& local_vector,
208	const UnknownManager& unknown_manager) const
209	{
210	size_t num_local_dofs = GetNumLocalDOFs(unknown_manager);	9,037✔
211
212	CopyVecToSTLvector(petsc_vector, local_vector, num_local_dofs);	9,037✔
213	}	9,037✔
214
215	void
216	SpatialDiscretization::LocalizePETScVectorWithGhosts(Vec petsc_vector,	130✔
217	std::vector<double>& local_vector,
218	const UnknownManager& unknown_manager) const
219	{
220	size_t num_local_dofs = GetNumLocalAndGhostDOFs(unknown_manager);	130✔
221
222	CopyVecToSTLvectorWithGhosts(petsc_vector, local_vector, num_local_dofs);	130✔
223	}	130✔
224
225	double
226	SpatialDiscretization::CartesianSpatialWeightFunction(const Vector3& point)	×
227	{
228	return 1.0;	×
229	}
230
231	double
UNCOV 232	SpatialDiscretization::CylindricalRZSpatialWeightFunction(const Vector3& point)	×
233	{
UNCOV 234	return 2.0 * M_PI * point[0];	×
235	}
236
237	double
238	SpatialDiscretization::Spherical1DSpatialWeightFunction(const Vector3& point)	×
239	{
240	const double r = point[2];	×
241	return 4.0 * M_PI * r * r;	×
242	}
243
244	SpatialDiscretization::SpatialWeightFunction
245	SpatialDiscretization::GetSpatialWeightingFunction() const	×
246	{
247	switch (grid_->GetCoordinateSystem())	×
248	{
249	case CoordinateSystemType::CARTESIAN:	×
250	return CartesianSpatialWeightFunction;	×
251	case CoordinateSystemType::CYLINDRICAL:	×
252	return CylindricalRZSpatialWeightFunction;	×
253	case CoordinateSystemType::SPHERICAL:	×
254	return Spherical1DSpatialWeightFunction;	×
255	case CoordinateSystemType::UNDEFINED:	×
256	default:	×
257	OpenSnLogicalError("Coordinate system undefined.");	×
258	}
259	}
260
261	} // namespace opensn

Open-Sn / opensn / 21773869314

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