22654533986

Committed 04 Mar 2026 12:48AM UTC coverage: 74.268% (+0.08%) from 74.192%

Build # 22654533986

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #960 from andrsd/fv-removal

Removing finite-volume spatial discretization

Run Details

0 of 2 new or added lines in 1 file covered. (0.0%)

183 existing lines in 9 files now uncovered.

19993 of 26920 relevant lines covered (74.27%)

67230185.92 hits per line

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

33.33

/framework/field_functions/field_function_grid_based.cc

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

#include "framework/field_functions/field_function_grid_based.h"
#include "framework/math/spatial_discretization/finite_element/piecewise_linear/piecewise_linear_continuous.h"
#include "framework/math/spatial_discretization/finite_element/piecewise_linear/piecewise_linear_discontinuous.h"
#include "framework/math/spatial_discretization/spatial_discretization.h"
#include "framework/mesh/mesh.h"
#include "framework/mesh/mesh_continuum/mesh_continuum.h"
#include "framework/mesh/mesh_continuum/grid_vtk_utils.h"
#include "framework/object_factory.h"
#include "framework/logging/log.h"
#include "framework/runtime.h"
#include <petsc.h>
#include <vtkUnstructuredGrid.h>
#include <vtkCellData.h>
#include <vtkPointData.h>
#include <vtkDoubleArray.h>

namespace opensn
{

InputParameters
FieldFunctionGridBased::GetInputParameters()
{
  InputParameters params = FieldFunction::GetInputParameters();

  params.AddOptionalParameter(
    "discretization", "PWLC", "The spatial discretization type to be used");
  params.AddOptionalParameter(
    "coord_sys", "cartesian", "Coordinate system to apply to element mappings");
  params.AddOptionalParameter("quadrature_order",
                              0,
                              "If supplied, will overwrite the default for the "
                              "specific discretization-coordinate system combination.");

  params.AddOptionalParameter(
    "initial_value", 0.0, "The initial value to assign to the field function");

  params.ConstrainParameterRange("discretization", AllowableRangeList::New({"PWLC", "PWLD"}));
  params.ConstrainParameterRange(
    "coord_sys", AllowableRangeList::New({"cartesian", "cylindrical", "spherical"}));

  return params;
}

FieldFunctionGridBased::FieldFunctionGridBased(const InputParameters& params)
  : FieldFunction(params),
    discretization_(MakeSpatialDiscretization(params)),
    ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),
    local_grid_bounding_box_(params.GetSharedPtrParam<MeshContinuum>("mesh")->GetLocalBoundingBox())
{
  ghosted_field_vector_->Set(params.GetParamValue<double>("initial_value"));
}

FieldFunctionGridBased::FieldFunctionGridBased(
  std::string name, std::shared_ptr<SpatialDiscretization>& discretization_ptr, Unknown unknown)
  : FieldFunction(std::move(name), std::move(unknown)),
    discretization_(discretization_ptr),
    ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),
    local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())
{
}

FieldFunctionGridBased::FieldFunctionGridBased(
  std::string name,
  std::shared_ptr<SpatialDiscretization>& discretization_ptr,
  Unknown unknown,
  const std::vector<double>& field_vector)
  : FieldFunction(std::move(name), std::move(unknown)),
    discretization_(discretization_ptr),
    ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),
    local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())
{
  OpenSnInvalidArgumentIf(field_vector.size() != ghosted_field_vector_->GetLocalSize(),
                          "Constructor called with incompatible size field vector.");

  ghosted_field_vector_->Set(field_vector);
}

FieldFunctionGridBased::FieldFunctionGridBased(
  std::string name,
  std::shared_ptr<SpatialDiscretization>& discretization_ptr,
  Unknown unknown,
  double field_value)
  : FieldFunction(std::move(name), std::move(unknown)),
    discretization_(discretization_ptr),
    ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),
    local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())
{
  ghosted_field_vector_->Set(field_value);
}

const SpatialDiscretization&
FieldFunctionGridBased::GetSpatialDiscretization() const
{
  return *discretization_;
}

std::vector<double>&
FieldFunctionGridBased::GetLocalFieldVector()
{
  return ghosted_field_vector_->GetLocalSTLData();
}

const std::vector<double>&
FieldFunctionGridBased::GetLocalFieldVector() const
{
  return ghosted_field_vector_->GetLocalSTLData();
}

std::vector<double>
FieldFunctionGridBased::GetGhostedFieldVector() const
{
  return ghosted_field_vector_->GetLocalSTLData();
}

void
FieldFunctionGridBased::UpdateFieldVector(const std::vector<double>& field_vector)
{
  OpenSnInvalidArgumentIf(field_vector.size() < ghosted_field_vector_->GetLocalSize(),
                          "Attempted update with a vector of insufficient size.");

  ghosted_field_vector_->Set(field_vector);
  ghosted_field_vector_->CommunicateGhostEntries();
}

void
FieldFunctionGridBased::UpdateFieldVector(const Vec& field_vector)
{
  ghosted_field_vector_->CopyLocalValues(field_vector);

  ghosted_field_vector_->CommunicateGhostEntries();
}

std::vector<double>
FieldFunctionGridBased::GetPointValue(const Vector3& point) const
{
  const auto& uk_man = GetUnknownManager();
  const size_t num_components = uk_man.GetTotalUnknownStructureSize();

  size_t local_num_point_hits = 0;
  std::vector<double> local_point_value(num_components, 0.0);

  const auto& xyz_min = local_grid_bounding_box_.first;
  const auto& xyz_max = local_grid_bounding_box_.second;

  const auto xmin = xyz_min.x;
  const auto ymin = xyz_min.y;
  const auto zmin = xyz_min.z;

  const auto xmax = xyz_max.x;
  const auto ymax = xyz_max.y;
  const auto zmax = xyz_max.z;

  const auto& field_vector = *ghosted_field_vector_;

  if (point.x >= xmin and point.x <= xmax and point.y >= ymin and point.y <= ymax and
      point.z >= zmin and point.z <= zmax)
  {
    const auto& grid = discretization_->GetGrid();
    for (const auto& cell : grid->local_cells)
    {
      if (grid->CheckPointInsideCell(cell, point))
      {
        const auto& cell_mapping = discretization_->GetCellMapping(cell);
        Vector<double> shape_values;
        cell_mapping.ShapeValues(point, shape_values);

        local_num_point_hits += 1;

        const auto num_nodes = cell_mapping.GetNumNodes();
        for (size_t c = 0; c < num_components; ++c)
        {
          for (size_t j = 0; j < num_nodes; ++j)
          {
            const auto dof_map = discretization_->MapDOFLocal(cell, j, uk_man, 0, c);
            const double dof_value = field_vector[dof_map];

            local_point_value[c] += dof_value * shape_values(j);
          } // for node i
        } // for component c
      } // if inside cell
    } // for cell
  } // if in bounding box

  // Communicate number of point hits
  size_t global_num_point_hits = 0;
  mpi_comm.all_reduce(local_num_point_hits, global_num_point_hits, mpi::op::sum<size_t>());

  std::vector<double> global_point_value(num_components, 0.0);
  mpi_comm.all_reduce(
    local_point_value.data(), 1, global_point_value.data(), mpi::op::sum<double>());

  Scale(global_point_value, 1.0 / static_cast<double>(global_num_point_hits));

  return global_point_value;
}

double
FieldFunctionGridBased::Evaluate(const Cell& cell, const Vector3& position, int component) const
{
  const auto& field_vector = *ghosted_field_vector_;

  const auto& cell_mapping = discretization_->GetCellMapping(cell);

  Vector<double> shape_values;
  cell_mapping.ShapeValues(position, shape_values);

  double value = 0.0;
  const size_t num_nodes = cell_mapping.GetNumNodes();
  for (size_t j = 0; j < num_nodes; ++j)
  {
    const auto dof_map = discretization_->MapDOFLocal(cell, j, GetUnknownManager(), 0, component);
    value += field_vector[dof_map] * shape_values(j);
  }

  return value;
}

void
FieldFunctionGridBased::ExportMultipleToPVTU(
  const std::string& file_base_name,
  const std::vector<std::shared_ptr<const FieldFunctionGridBased>>& ff_list)
{
  const std::string fname = "FieldFunctionGridBased::ExportMultipleToPVTU";
  log.Log() << "Exporting field functions to PVTU with file base \"" << file_base_name << "\"";

  if (ff_list.empty())
    throw std::logic_error(fname + ": Cannot be used with empty field-function"
                                   " list");

  // Setup master and check slaves
  const auto& master_ff_ptr = ff_list.front();
  const auto& master_ff = *master_ff_ptr;

  for (const auto& ff_ptr : ff_list)
    if (ff_ptr != master_ff_ptr)
      if (ff_ptr->discretization_->GetGrid() != master_ff_ptr->discretization_->GetGrid())
        throw std::logic_error(fname +
                               ": Cannot be used with field functions based on different grids.");

  // Get grid
  const auto& grid = master_ff.discretization_->GetGrid();

  auto ugrid = PrepareVtkUnstructuredGrid(grid);

  // Upload cell/point data
  auto* point_data = ugrid->GetPointData();
  for (const auto& ff_ptr : ff_list)
  {
    const auto field_vector = ff_ptr->GetGhostedFieldVector();

    const auto& uk_man = ff_ptr->GetUnknownManager();
    const auto& unknown = ff_ptr->GetUnknown();
    const auto& sdm = ff_ptr->discretization_;
    const size_t num_comps = unknown.GetNumComponents();

    for (uint c = 0; c < num_comps; ++c)
    {
      std::string component_name = ff_ptr->GetName() + unknown.name;
      if (num_comps > 1)
        component_name += unknown.component_names[c];

      vtkNew<vtkDoubleArray> point_array;

      point_array->SetName(component_name.c_str());

      // Populate the array here
      for (const auto& cell : grid->local_cells)
      {
        const size_t num_nodes = sdm->GetCellNumNodes(cell);

        if (num_nodes == cell.vertex_ids.size())
        {
          for (size_t n = 0; n < num_nodes; ++n)
          {
            const auto nmap = sdm->MapDOFLocal(cell, n, uk_man, 0, c);

            const double field_value = field_vector[nmap];

            point_array->InsertNextValue(field_value);
          } // for node
        }
        else
        {
          double node_average = 0.0;
          for (size_t n = 0; n < num_nodes; ++n)
          {
            const auto nmap = sdm->MapDOFLocal(cell, n, uk_man, 0, c);

            const double field_value = field_vector[nmap];
            node_average += field_value;
          } // for node
          node_average /= static_cast<double>(num_nodes);
          for (std::size_t n = 0; n < cell.vertex_ids.size(); ++n)
          {
            point_array->InsertNextValue(node_average);
          } // for vertex
        }

      } // for cell

      point_data->AddArray(point_array);
    } // for component
  } // for ff_ptr

  WritePVTUFiles(ugrid, file_base_name);

  log.Log() << "Done exporting field functions to PVTU.";
  opensn::mpi_comm.barrier();
}

std::shared_ptr<SpatialDiscretization>
FieldFunctionGridBased::MakeSpatialDiscretization(const InputParameters& params)
{
  const auto grid = params.GetSharedPtrParam<MeshContinuum>("mesh");
  const auto sdm_type = params.GetParamValue<std::string>("discretization");

  std::string cs = "cartesian";
  if (params.IsParameterValid("coord_sys"))
    cs = params.GetParamValue<std::string>("coord_sys");

  QuadratureOrder q_order = QuadratureOrder::SECOND;

  if (params.IsParameterValid("quadrature_order"))
  {
    const auto max_order = static_cast<uint32_t>(QuadratureOrder::FORTYTHIRD);
    const auto q_order_int = params.GetParamValue<uint32_t>("quadrature_order");
    OpenSnInvalidArgumentIf(q_order_int > max_order,
                            "Invalid quadrature point order " + std::to_string(q_order_int));
    q_order = static_cast<QuadratureOrder>(q_order_int);
  }
  else // Defaulted
  {
    if (cs == "cartesian")
      q_order = QuadratureOrder::SECOND;
    else if (cs == "cylindrical")
      q_order = QuadratureOrder::THIRD;
    else if (cs == "spherical")
      q_order = QuadratureOrder::FOURTH;
    else
      throw std::logic_error("Unrecognized coordinate system.");
  }

  if (sdm_type == "PWLC")
    return PieceWiseLinearContinuous::New(grid, q_order);
  else if (sdm_type == "PWLD")
    return PieceWiseLinearDiscontinuous::New(grid, q_order);

  // If not returned by now
  OpenSnInvalidArgument("Unsupported discretization \"" + sdm_type + "\"");
}

std::unique_ptr<GhostedParallelSTLVector>
FieldFunctionGridBased::MakeFieldVector(const SpatialDiscretization& discretization,
                                        const UnknownManager& uk_man)
{
  auto field = std::make_unique<GhostedParallelSTLVector>(discretization.GetNumLocalDOFs(uk_man),
                                                          discretization.GetNumGlobalDOFs(uk_man),
                                                          discretization.GetGhostDOFIndices(uk_man),
                                                          mpi_comm);

  return field;
}

} // namespace opensn

1	// SPDX-FileCopyrightText: 2024 The OpenSn Authors <https://open-sn.github.io/opensn/>
2	// SPDX-License-Identifier: MIT
3
4	#include "framework/field_functions/field_function_grid_based.h"
5	#include "framework/math/spatial_discretization/finite_element/piecewise_linear/piecewise_linear_continuous.h"
6	#include "framework/math/spatial_discretization/finite_element/piecewise_linear/piecewise_linear_discontinuous.h"
7	#include "framework/math/spatial_discretization/spatial_discretization.h"
8	#include "framework/mesh/mesh.h"
9	#include "framework/mesh/mesh_continuum/mesh_continuum.h"
10	#include "framework/mesh/mesh_continuum/grid_vtk_utils.h"
11	#include "framework/object_factory.h"
12	#include "framework/logging/log.h"
13	#include "framework/runtime.h"
14	#include <petsc.h>
15	#include <vtkUnstructuredGrid.h>
16	#include <vtkCellData.h>
17	#include <vtkPointData.h>
18	#include <vtkDoubleArray.h>
19
20	namespace opensn
21	{
22
23	InputParameters
24	FieldFunctionGridBased::GetInputParameters()	×
25	{
26	InputParameters params = FieldFunction::GetInputParameters();	×
27
NEW 28	params.AddOptionalParameter(	×
29	"discretization", "PWLC", "The spatial discretization type to be used");
UNCOV 30	params.AddOptionalParameter(	×
31	"coord_sys", "cartesian", "Coordinate system to apply to element mappings");
32	params.AddOptionalParameter("quadrature_order",	×
33	0,
34	"If supplied, will overwrite the default for the "
35	"specific discretization-coordinate system combination.");
36
37	params.AddOptionalParameter(	×
38	"initial_value", 0.0, "The initial value to assign to the field function");
39
NEW 40	params.ConstrainParameterRange("discretization", AllowableRangeList::New({"PWLC", "PWLD"}));	×
41	params.ConstrainParameterRange(	×
42	"coord_sys", AllowableRangeList::New({"cartesian", "cylindrical", "spherical"}));	×
43
44	return params;	×
45	}	×
46
47	FieldFunctionGridBased::FieldFunctionGridBased(const InputParameters& params)	×
48	: FieldFunction(params),
49	discretization_(MakeSpatialDiscretization(params)),	×
50	ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),	×
51	local_grid_bounding_box_(params.GetSharedPtrParam<MeshContinuum>("mesh")->GetLocalBoundingBox())	×
52	{
53	ghosted_field_vector_->Set(params.GetParamValue<double>("initial_value"));	×
54	}	×
55
56	FieldFunctionGridBased::FieldFunctionGridBased(	78,004✔
57	std::string name, std::shared_ptr<SpatialDiscretization>& discretization_ptr, Unknown unknown)	78,004✔
58	: FieldFunction(std::move(name), std::move(unknown)),
59	discretization_(discretization_ptr),	156,008✔
60	ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),	78,004✔
61	local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())	234,012✔
62	{
63	}	78,004✔
64
65	FieldFunctionGridBased::FieldFunctionGridBased(	×
66	std::string name,
67	std::shared_ptr<SpatialDiscretization>& discretization_ptr,
68	Unknown unknown,
69	const std::vector<double>& field_vector)	×
70	: FieldFunction(std::move(name), std::move(unknown)),
71	discretization_(discretization_ptr),	×
72	ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),	×
73	local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())	×
74	{
75	OpenSnInvalidArgumentIf(field_vector.size() != ghosted_field_vector_->GetLocalSize(),	×
76	"Constructor called with incompatible size field vector.");
77
78	ghosted_field_vector_->Set(field_vector);	×
79	}	×
80
81	FieldFunctionGridBased::FieldFunctionGridBased(	×
82	std::string name,
83	std::shared_ptr<SpatialDiscretization>& discretization_ptr,
84	Unknown unknown,
85	double field_value)	×
86	: FieldFunction(std::move(name), std::move(unknown)),
87	discretization_(discretization_ptr),	×
88	ghosted_field_vector_(MakeFieldVector(*discretization_, GetUnknownManager())),	×
89	local_grid_bounding_box_(discretization_->GetGrid()->GetLocalBoundingBox())	×
90	{
91	ghosted_field_vector_->Set(field_value);	×
92	}	×
93
94	const SpatialDiscretization&
95	FieldFunctionGridBased::GetSpatialDiscretization() const	14,394✔
96	{
97	return *discretization_;	14,394✔
98	}
99
100	std::vector<double>&
101	FieldFunctionGridBased::GetLocalFieldVector()	29✔
102	{
103	return ghosted_field_vector_->GetLocalSTLData();	29✔
104	}
105
106	const std::vector<double>&
107	FieldFunctionGridBased::GetLocalFieldVector() const	×
108	{
109	return ghosted_field_vector_->GetLocalSTLData();	×
110	}
111
112	std::vector<double>
113	FieldFunctionGridBased::GetGhostedFieldVector() const	13,156✔
114	{
115	return ghosted_field_vector_->GetLocalSTLData();	13,156✔
116	}
117
118	void
119	FieldFunctionGridBased::UpdateFieldVector(const std::vector<double>& field_vector)	85,037✔
120	{
121	OpenSnInvalidArgumentIf(field_vector.size() < ghosted_field_vector_->GetLocalSize(),	85,037✔
122	"Attempted update with a vector of insufficient size.");
123
124	ghosted_field_vector_->Set(field_vector);	85,037✔
125	ghosted_field_vector_->CommunicateGhostEntries();	85,037✔
126	}	85,037✔
127
128	void
129	FieldFunctionGridBased::UpdateFieldVector(const Vec& field_vector)	×
130	{
131	ghosted_field_vector_->CopyLocalValues(field_vector);	×
132
133	ghosted_field_vector_->CommunicateGhostEntries();	×
134	}	×
135
136	std::vector<double>
137	FieldFunctionGridBased::GetPointValue(const Vector3& point) const	×
138	{
139	const auto& uk_man = GetUnknownManager();	×
140	const size_t num_components = uk_man.GetTotalUnknownStructureSize();	×
141
142	size_t local_num_point_hits = 0;	×
143	std::vector<double> local_point_value(num_components, 0.0);	×
144
145	const auto& xyz_min = local_grid_bounding_box_.first;	×
146	const auto& xyz_max = local_grid_bounding_box_.second;	×
147
148	const auto xmin = xyz_min.x;	×
149	const auto ymin = xyz_min.y;	×
150	const auto zmin = xyz_min.z;	×
151
152	const auto xmax = xyz_max.x;	×
153	const auto ymax = xyz_max.y;	×
154	const auto zmax = xyz_max.z;	×
155
156	const auto& field_vector = *ghosted_field_vector_;	×
157
158	if (point.x >= xmin and point.x <= xmax and point.y >= ymin and point.y <= ymax and	×
159	point.z >= zmin and point.z <= zmax)	×
160	{
161	const auto& grid = discretization_->GetGrid();	×
162	for (const auto& cell : grid->local_cells)	×
163	{
164	if (grid->CheckPointInsideCell(cell, point))	×
165	{
166	const auto& cell_mapping = discretization_->GetCellMapping(cell);	×
167	Vector<double> shape_values;	×
168	cell_mapping.ShapeValues(point, shape_values);	×
169
170	local_num_point_hits += 1;	×
171
172	const auto num_nodes = cell_mapping.GetNumNodes();	×
173	for (size_t c = 0; c < num_components; ++c)	×
174	{
175	for (size_t j = 0; j < num_nodes; ++j)	×
176	{
177	const auto dof_map = discretization_->MapDOFLocal(cell, j, uk_man, 0, c);	×
178	const double dof_value = field_vector[dof_map];	×
179
180	local_point_value[c] += dof_value * shape_values(j);	×
181	} // for node i
182	} // for component c
183	} // if inside cell	×
184	} // for cell
185	} // if in bounding box	×
186
187	// Communicate number of point hits
188	size_t global_num_point_hits = 0;	×
189	mpi_comm.all_reduce(local_num_point_hits, global_num_point_hits, mpi::op::sum<size_t>());	×
190
191	std::vector<double> global_point_value(num_components, 0.0);	×
192	mpi_comm.all_reduce(	×
193	local_point_value.data(), 1, global_point_value.data(), mpi::op::sum<double>());	×
194
195	Scale(global_point_value, 1.0 / static_cast<double>(global_num_point_hits));	×
196
197	return global_point_value;	×
198	}	×
199
200	double
201	FieldFunctionGridBased::Evaluate(const Cell& cell, const Vector3& position, int component) const	×
202	{
203	const auto& field_vector = *ghosted_field_vector_;	×
204
205	const auto& cell_mapping = discretization_->GetCellMapping(cell);	×
206
207	Vector<double> shape_values;	×
208	cell_mapping.ShapeValues(position, shape_values);	×
209
210	double value = 0.0;	×
211	const size_t num_nodes = cell_mapping.GetNumNodes();	×
212	for (size_t j = 0; j < num_nodes; ++j)	×
213	{
214	const auto dof_map = discretization_->MapDOFLocal(cell, j, GetUnknownManager(), 0, component);	×
215	value += field_vector[dof_map] * shape_values(j);	×
216	}
217
218	return value;	×
219	}	×
220
221	void
222	FieldFunctionGridBased::ExportMultipleToPVTU(	70✔
223	const std::string& file_base_name,
224	const std::vector<std::shared_ptr<const FieldFunctionGridBased>>& ff_list)
225	{
226	const std::string fname = "FieldFunctionGridBased::ExportMultipleToPVTU";	70✔
227	log.Log() << "Exporting field functions to PVTU with file base \"" << file_base_name << "\"";	210✔
228
229	if (ff_list.empty())	70✔
230	throw std::logic_error(fname + ": Cannot be used with empty field-function"	×
231	" list");	×
232
233	// Setup master and check slaves
234	const auto& master_ff_ptr = ff_list.front();	70✔
235	const auto& master_ff = *master_ff_ptr;	70✔
236
237	for (const auto& ff_ptr : ff_list)	6,024✔
238	if (ff_ptr != master_ff_ptr)	5,954✔
239	if (ff_ptr->discretization_->GetGrid() != master_ff_ptr->discretization_->GetGrid())	17,652✔
240	throw std::logic_error(fname +	×
241	": Cannot be used with field functions based on different grids.");	×
242
243	// Get grid
244	const auto& grid = master_ff.discretization_->GetGrid();	70✔
245
246	auto ugrid = PrepareVtkUnstructuredGrid(grid);	140✔
247
248	// Upload cell/point data
249	auto* point_data = ugrid->GetPointData();	70✔
250	for (const auto& ff_ptr : ff_list)	6,024✔
251	{
252	const auto field_vector = ff_ptr->GetGhostedFieldVector();	5,954✔
253
254	const auto& uk_man = ff_ptr->GetUnknownManager();	5,954✔
255	const auto& unknown = ff_ptr->GetUnknown();	5,954✔
256	const auto& sdm = ff_ptr->discretization_;	5,954✔
257	const size_t num_comps = unknown.GetNumComponents();	5,954✔
258
259	for (uint c = 0; c < num_comps; ++c)	11,984✔
260	{
261	std::string component_name = ff_ptr->GetName() + unknown.name;	6,030✔
262	if (num_comps > 1)	6,030✔
263	component_name += unknown.component_names[c];	80✔
264
265	vtkNew<vtkDoubleArray> point_array;	6,030✔
266
267	point_array->SetName(component_name.c_str());	6,030✔
268
269	// Populate the array here
270	for (const auto& cell : grid->local_cells)	1,336,680✔
271	{
272	const size_t num_nodes = sdm->GetCellNumNodes(cell);	1,330,650✔
273
274	if (num_nodes == cell.vertex_ids.size())	1,330,650✔
275	{
276	for (size_t n = 0; n < num_nodes; ++n)	6,540,490✔
277	{
278	const auto nmap = sdm->MapDOFLocal(cell, n, uk_man, 0, c);	5,209,840✔
279
280	const double field_value = field_vector[nmap];	5,209,840✔
281
282	point_array->InsertNextValue(field_value);	5,209,840✔
283	} // for node
284	}
285	else
286	{
287	double node_average = 0.0;
UNCOV 288	for (size_t n = 0; n < num_nodes; ++n)	×
289	{
UNCOV 290	const auto nmap = sdm->MapDOFLocal(cell, n, uk_man, 0, c);	×
291
UNCOV 292	const double field_value = field_vector[nmap];	×
UNCOV 293	node_average += field_value;	×
294	} // for node
UNCOV 295	node_average /= static_cast<double>(num_nodes);	×
UNCOV 296	for (std::size_t n = 0; n < cell.vertex_ids.size(); ++n)	×
297	{
UNCOV 298	point_array->InsertNextValue(node_average);	×
299	} // for vertex
300	}
301
302	} // for cell
303
304	point_data->AddArray(point_array);	6,030✔
305	} // for component	6,030✔
306	} // for ff_ptr	5,954✔
307
308	WritePVTUFiles(ugrid, file_base_name);	70✔
309
310	log.Log() << "Done exporting field functions to PVTU.";	140✔
311	opensn::mpi_comm.barrier();	70✔
312	}	140✔
313
314	std::shared_ptr<SpatialDiscretization>
315	FieldFunctionGridBased::MakeSpatialDiscretization(const InputParameters& params)	×
316	{
317	const auto grid = params.GetSharedPtrParam<MeshContinuum>("mesh");	×
318	const auto sdm_type = params.GetParamValue<std::string>("discretization");	×
319
320	std::string cs = "cartesian";	×
321	if (params.IsParameterValid("coord_sys"))	×
322	cs = params.GetParamValue<std::string>("coord_sys");	×
323
324	QuadratureOrder q_order = QuadratureOrder::SECOND;	×
325
326	if (params.IsParameterValid("quadrature_order"))	×
327	{
328	const auto max_order = static_cast<uint32_t>(QuadratureOrder::FORTYTHIRD);	×
329	const auto q_order_int = params.GetParamValue<uint32_t>("quadrature_order");	×
330	OpenSnInvalidArgumentIf(q_order_int > max_order,	×
331	"Invalid quadrature point order " + std::to_string(q_order_int));
332	q_order = static_cast<QuadratureOrder>(q_order_int);	×
333	}
334	else // Defaulted
335	{
336	if (cs == "cartesian")	×
337	q_order = QuadratureOrder::SECOND;
338	else if (cs == "cylindrical")	×
339	q_order = QuadratureOrder::THIRD;
340	else if (cs == "spherical")	×
341	q_order = QuadratureOrder::FOURTH;
342	else
343	throw std::logic_error("Unrecognized coordinate system.");	×
344	}
345
346	if (sdm_type == "PWLC")	×
347	return PieceWiseLinearContinuous::New(grid, q_order);	×
348	else if (sdm_type == "PWLD")	×
349	return PieceWiseLinearDiscontinuous::New(grid, q_order);	×
350
351	// If not returned by now
352	OpenSnInvalidArgument("Unsupported discretization \"" + sdm_type + "\"");	×
353	}	×
354
355	std::unique_ptr<GhostedParallelSTLVector>
356	FieldFunctionGridBased::MakeFieldVector(const SpatialDiscretization& discretization,	78,004✔
357	const UnknownManager& uk_man)
358	{
359	auto field = std::make_unique<GhostedParallelSTLVector>(discretization.GetNumLocalDOFs(uk_man),	312,016✔
360	discretization.GetNumGlobalDOFs(uk_man),	156,008✔
361	discretization.GetGhostDOFIndices(uk_man),	78,004✔
362	mpi_comm);	78,004✔
363
364	return field;	78,004✔
365	}
366
367	} // namespace opensn

Open-Sn / opensn / 22654533986

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