16766938693

Committed 05 Aug 2025 03:26PM UTC coverage: 73.386% (+0.1%) from 73.282%

Build # 16766938693

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #693 from andrsd/move-bnds

Moving sweep boundaries into `DiscreteOrdinatesProblem`

Coverage Stats

247 of 311 new or added lines in 19 files covered. (79.42%)

625 existing lines in 49 files now uncovered.

18320 of 24964 relevant lines covered (73.39%)

43106214.05 hits per line

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

48.08

/framework/math/math.cc

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

#include "framework/math/math.h"
#include "framework/runtime.h"
#include <cassert>

namespace opensn
{

double
Factorial(const int x)
{
  double factorial_value = 1.0;
  for (int i = 2; i <= x; ++i)
    factorial_value *= i;

  return factorial_value;
}

std::pair<double, double>
OmegaToPhiThetaSafe(const Vector3& omega)
{
  // Notes: asin maps [-1,+1] to [-pi/2,+pi/2]
  //        acos maps [-1,+1] to [0,pi]
  // This mapping requires some logic for determining the azimuthal angle.
  //
  const auto omega_hat = omega.Normalized();

  double mu = omega_hat.z;
  mu = std::min(mu, 1.0);
  mu = std::max(mu, -1.0);

  double theta = acos(mu);

  // Handling omega aligned to k_hat
  if (std::fabs(omega_hat.z) < 1.0e-16)
    return {0.0, theta};

  double cos_phi = omega_hat.x / sin(theta);
  cos_phi = std::min(cos_phi, 1.0);
  cos_phi = std::max(cos_phi, -1.0);

  // Computing varphi for NE and NW quadrant
  if (omega_hat.y >= 0.0)
    return {acos(cos_phi), theta};
  // Computing varphi for SE and SW quadrant
  else
    return {2.0 * M_PI - acos(cos_phi), theta};
}

void
PrintVector(const std::vector<double>& x)
{
  for (auto& xi : x)
    std::cout << xi << ' ';
  std::cout << std::endl;
}

void
Scale(std::vector<double>& x, const double& val)
{
  for (double& xi : x)
    xi *= val;
}

void
Set(std::vector<double>& x, const double& val)
{
  for (double& xi : x)
    xi = val;
}

double
Dot(const std::vector<double>& x, const std::vector<double>& y)
{
  // Error Checks
  assert(x.size() > 0);
  assert(y.size() > 0);
  assert(x.size() == y.size());
  // Local Variables
  size_t n = x.size();
  double val = 0.0;

  for (size_t i = 0; i != n; ++i)
    val += x[i] * y[i];

  return val;
}

std::vector<double>
Mult(const std::vector<double>& x, const double& val)
{
  size_t n = x.size();
  std::vector<double> y(n);

  for (size_t i = 0; i != n; ++i)
    y[i] = val * x[i];

  return y;
}

double
L1Norm(const std::vector<double>& x)
{
  // Local Variables
  size_t n = x.size();
  double val = 0.0;

  for (size_t i = 0; i != n; ++i)
    val += std::fabs(x[i]);

  return val;
}

double
L2Norm(const std::vector<double>& x)
{
  // Local Variables
  size_t n = x.size();
  double val = 0.0;

  for (size_t i = 0; i != n; ++i)
    val += x[i] * x[i];

  return std::sqrt(val);
}

double
LInfNorm(const std::vector<double>& x)
{
  // Local Variables
  size_t n = x.size();
  double val = 0.0;

  for (size_t i = 0; i != n; ++i)
    val += std::max(std::fabs(x[i]), val);

  return val;
}

double
LpNorm(const std::vector<double>& x, const double& p)
{
  // Local Variables
  size_t n = x.size();
  double val = 0.0;

  for (size_t i = 0; i != n; ++i)
    val += std::pow(std::fabs(x[i]), p);

  return std::pow(val, 1.0 / p);
}

std::vector<double>
operator+(const std::vector<double>& a, const std::vector<double>& b)
{
  assert(a.size() == b.size());
  std::vector<double> result(a.size(), 0.0);

  for (size_t i = 0; i < a.size(); ++i)
    result[i] = a[i] + b[i];

  return result;
}

std::vector<double>
operator-(const std::vector<double>& a, const std::vector<double>& b)
{
  assert(a.size() == b.size());
  std::vector<double> result(a.size(), 0.0);

  for (size_t i = 0; i < a.size(); ++i)
    result[i] = a[i] - b[i];

  return result;
}

double
ComputeL2Change(std::vector<double>& x, std::vector<double>& y)
{
  assert(x.size() == y.size());

  double norm = 0.0;
  for (auto i = 0; i < x.size(); ++i)
  {
    double val = x[i] - y[i];
    norm += val * val;
  }

  double global_norm = 0.0;
  mpi_comm.all_reduce<double>(norm, global_norm, mpi::op::sum<double>());

  return std::sqrt(global_norm);
}

double
ComputePointwiseChange(std::vector<double>& x, std::vector<double>& y)
{
  assert(x.size() == y.size());

  double pw_change = 0.0;
  for (auto i = 0; i < x.size(); ++i)
  {
    double max = std::max(x[i], y[i]);
    double delta = std::fabs(x[i] - y[i]);
    if (max >= std::numeric_limits<double>::min())
      pw_change = std::max(delta / max, pw_change);
    else
      pw_change = std::max(delta, pw_change);
  }

  double global_pw_change = 0.0;
  mpi_comm.all_reduce<double>(pw_change, global_pw_change, mpi::op::max<double>());

  return global_pw_change;
}

} // 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/math.h"
5	#include "framework/runtime.h"
6	#include <cassert>
7
8	namespace opensn
9	{
10
11	double
12	Factorial(const int x)	29,501,968✔
13	{
14	double factorial_value = 1.0;	29,501,968✔
15	for (int i = 2; i <= x; ++i)	44,266,392✔
16	factorial_value *= i;	14,764,424✔
17
18	return factorial_value;	29,501,968✔
19	}
20
21	std::pair<double, double>
22	OmegaToPhiThetaSafe(const Vector3& omega)	839,857✔
23	{
24	// Notes: asin maps [-1,+1] to [-pi/2,+pi/2]
25	// acos maps [-1,+1] to [0,pi]
26	// This mapping requires some logic for determining the azimuthal angle.
27	//
28	const auto omega_hat = omega.Normalized();	839,857✔
29
30	double mu = omega_hat.z;	839,857✔
31	mu = std::min(mu, 1.0);	839,857✔
32	mu = std::max(mu, -1.0);	839,857✔
33
34	double theta = acos(mu);	839,857✔
35
36	// Handling omega aligned to k_hat
37	if (std::fabs(omega_hat.z) < 1.0e-16)	839,857✔
UNCOV 38	return {0.0, theta};	×
39
40	double cos_phi = omega_hat.x / sin(theta);	839,857✔
41	cos_phi = std::min(cos_phi, 1.0);	839,857✔
42	cos_phi = std::max(cos_phi, -1.0);	839,857✔
43
44	// Computing varphi for NE and NW quadrant
45	if (omega_hat.y >= 0.0)	839,857✔
46	return {acos(cos_phi), theta};	421,089✔
47	// Computing varphi for SE and SW quadrant
48	else
49	return {2.0 * M_PI - acos(cos_phi), theta};	418,768✔
50	}
51
52	void
UNCOV 53	PrintVector(const std::vector<double>& x)	×
54	{
55	for (auto& xi : x)	×
56	std::cout << xi << ' ';	×
57	std::cout << std::endl;	×
UNCOV 58	}	×
59
60	void
61	Scale(std::vector<double>& x, const double& val)	7,248✔
62	{
63	for (double& xi : x)	500,444,404✔
64	xi *= val;	500,437,156✔
65	}	7,248✔
66
67	void
68	Set(std::vector<double>& x, const double& val)	496,348✔
69	{
70	for (double& xi : x)	1,016,326,363✔
71	xi = val;	1,015,830,015✔
72	}	496,348✔
73
74	double
UNCOV 75	Dot(const std::vector<double>& x, const std::vector<double>& y)	×
76	{
77	// Error Checks
78	assert(x.size() > 0);	×
79	assert(y.size() > 0);	×
UNCOV 80	assert(x.size() == y.size());	×
81	// Local Variables
82	size_t n = x.size();	×
UNCOV 83	double val = 0.0;	×
84
85	for (size_t i = 0; i != n; ++i)	×
UNCOV 86	val += x[i] * y[i];	×
87
UNCOV 88	return val;	×
89	}
90
91	std::vector<double>
UNCOV 92	Mult(const std::vector<double>& x, const double& val)	×
93	{
94	size_t n = x.size();	×
UNCOV 95	std::vector<double> y(n);	×
96
97	for (size_t i = 0; i != n; ++i)	×
UNCOV 98	y[i] = val * x[i];	×
99
UNCOV 100	return y;	×
101	}
102
103	double
UNCOV 104	L1Norm(const std::vector<double>& x)	×
105	{
106	// Local Variables
107	size_t n = x.size();	×
UNCOV 108	double val = 0.0;	×
109
110	for (size_t i = 0; i != n; ++i)	×
UNCOV 111	val += std::fabs(x[i]);	×
112
UNCOV 113	return val;	×
114	}
115
116	double
UNCOV 117	L2Norm(const std::vector<double>& x)	×
118	{
119	// Local Variables
120	size_t n = x.size();	×
UNCOV 121	double val = 0.0;	×
122
123	for (size_t i = 0; i != n; ++i)	×
UNCOV 124	val += x[i] * x[i];	×
125
UNCOV 126	return std::sqrt(val);	×
127	}
128
129	double
UNCOV 130	LInfNorm(const std::vector<double>& x)	×
131	{
132	// Local Variables
133	size_t n = x.size();	×
UNCOV 134	double val = 0.0;	×
135
136	for (size_t i = 0; i != n; ++i)	×
UNCOV 137	val += std::max(std::fabs(x[i]), val);	×
138
UNCOV 139	return val;	×
140	}
141
142	double
UNCOV 143	LpNorm(const std::vector<double>& x, const double& p)	×
144	{
145	// Local Variables
146	size_t n = x.size();	×
UNCOV 147	double val = 0.0;	×
148
149	for (size_t i = 0; i != n; ++i)	×
UNCOV 150	val += std::pow(std::fabs(x[i]), p);	×
151
UNCOV 152	return std::pow(val, 1.0 / p);	×
153	}
154
155	std::vector<double>
156	operator+(const std::vector<double>& a, const std::vector<double>& b)	9,852✔
157	{
158	assert(a.size() == b.size());	9,852✔
159	std::vector<double> result(a.size(), 0.0);	9,852✔
160
161	for (size_t i = 0; i < a.size(); ++i)	44,157,052✔
162	result[i] = a[i] + b[i];	44,147,200✔
163
164	return result;	9,852✔
165	}
166
167	std::vector<double>
168	operator-(const std::vector<double>& a, const std::vector<double>& b)	2,528✔
169	{
170	assert(a.size() == b.size());	2,528✔
171	std::vector<double> result(a.size(), 0.0);	2,528✔
172
173	for (size_t i = 0; i < a.size(); ++i)	8,092,128✔
174	result[i] = a[i] - b[i];	8,089,600✔
175
176	return result;	2,528✔
177	}
178
179	double
UNCOV 180	ComputeL2Change(std::vector<double>& x, std::vector<double>& y)	×
181	{
UNCOV 182	assert(x.size() == y.size());	×
183
184	double norm = 0.0;	×
UNCOV 185	for (auto i = 0; i < x.size(); ++i)	×
186	{
187	double val = x[i] - y[i];	×
UNCOV 188	norm += val * val;	×
189	}
190
191	double global_norm = 0.0;	×
UNCOV 192	mpi_comm.all_reduce<double>(norm, global_norm, mpi::op::sum<double>());	×
193
UNCOV 194	return std::sqrt(global_norm);	×
195	}
196
197	double
198	ComputePointwiseChange(std::vector<double>& x, std::vector<double>& y)	5,268✔
199	{
200	assert(x.size() == y.size());	5,268✔
201
202	double pw_change = 0.0;	5,268✔
203	for (auto i = 0; i < x.size(); ++i)	109,284,000✔
204	{
205	double max = std::max(x[i], y[i]);	109,278,732✔
206	double delta = std::fabs(x[i] - y[i]);	109,278,732✔
207	if (max >= std::numeric_limits<double>::min())	109,278,732✔
208	pw_change = std::max(delta / max, pw_change);	82,527,561✔
209	else
210	pw_change = std::max(delta, pw_change);	26,751,171✔
211	}
212
213	double global_pw_change = 0.0;	5,268✔
214	mpi_comm.all_reduce<double>(pw_change, global_pw_change, mpi::op::max<double>());	5,268✔
215
216	return global_pw_change;	5,268✔
217	}
218
219	} // namespace opensn

Open-Sn / opensn / 16766938693

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