20470828925

Committed 23 Dec 2025 08:24PM UTC coverage: 81.898% (-0.04%) from 81.94%

Build # 20470828925

Build Type

Pull #3550

github

Committed by

web-flow

Commit Message

Merge b067e7f24 into 3f06a42ab

Pull Request Pull Request #3550: [Point Detector] Add distribution get_pdf functionality

Run Details

17078 of 23776 branches covered (71.83%)

Branch coverage included in aggregate %.

48 of 212 new or added lines in 13 files covered. (22.64%)

4 existing lines in 2 files now uncovered.

55208 of 64487 relevant lines covered (85.61%)

43370386.44 hits per line

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

75.65

/src/distribution_multi.cpp

#include "openmc/distribution_multi.h"

#include <algorithm> // for move
#include <cmath>     // for sqrt, sin, cos, max

#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/math_functions.h"
#include "openmc/random_dist.h"
#include "openmc/random_lcg.h"
#include "openmc/xml_interface.h"

namespace openmc {

unique_ptr<UnitSphereDistribution> UnitSphereDistribution::create(
  pugi::xml_node node)
{
  // Check for type of angular distribution
  std::string type;
  if (check_for_node(node, "type"))
    type = get_node_value(node, "type", true, true);
  if (type == "isotropic") {
    return UPtrAngle {new Isotropic()};
  } else if (type == "monodirectional") {
    return UPtrAngle {new Monodirectional(node)};
  } else if (type == "mu-phi") {
    return UPtrAngle {new PolarAzimuthal(node)};
  } else {
    fatal_error(fmt::format(
      "Invalid angular distribution for external source: {}", type));
  }
}

//==============================================================================
// UnitSphereDistribution implementation
//==============================================================================

UnitSphereDistribution::UnitSphereDistribution(pugi::xml_node node)
{
  // Read reference directional unit vector
  if (check_for_node(node, "reference_uvw")) {
    auto u_ref = get_node_array<double>(node, "reference_uvw");
    if (u_ref.size() != 3)
      fatal_error("Angular distribution reference direction must have "
                  "three parameters specified.");
    u_ref_ = Direction(u_ref.data());
    u_ref_ /= u_ref_.norm();
  }
}

//==============================================================================
// PolarAzimuthal implementation
//==============================================================================

PolarAzimuthal::PolarAzimuthal(Direction u, UPtrDist mu, UPtrDist phi)
  : UnitSphereDistribution {u}, mu_ {std::move(mu)}, phi_ {std::move(phi)}
{}

PolarAzimuthal::PolarAzimuthal(pugi::xml_node node)
  : UnitSphereDistribution {node}
{
  // Read reference directional unit vector
  if (check_for_node(node, "reference_vwu")) {
    auto v_ref = get_node_array<double>(node, "reference_vwu");
    if (v_ref.size() != 3)
      fatal_error("Angular distribution reference v direction must have "
                  "three parameters specified.");
    v_ref_ = Direction(v_ref.data());
    v_ref_ /= v_ref_.norm();
  }
  w_ref_ = u_ref_.cross(v_ref_);
  if (check_for_node(node, "mu")) {
    pugi::xml_node node_dist = node.child("mu");
    mu_ = distribution_from_xml(node_dist);
  } else {
    mu_ = UPtrDist {new Uniform(-1., 1.)};
  }

  if (check_for_node(node, "phi")) {
    pugi::xml_node node_dist = node.child("phi");
    phi_ = distribution_from_xml(node_dist);
  } else {
    phi_ = UPtrDist {new Uniform(0.0, 2.0 * PI)};
  }
}

Direction PolarAzimuthal::sample(uint64_t* seed) const
{
  // Sample cosine of polar angle
  double mu = mu_->sample(seed);
  if (mu == 1.0)
    return u_ref_;
  if (mu == -1.0)
    return -u_ref_;

  // Sample azimuthal angle
  double phi = phi_->sample(seed);

  double f = std::sqrt(1 - mu * mu);

  return mu * u_ref_ + f * std::cos(phi) * v_ref_ + f * std::sin(phi) * w_ref_;
}

double PolarAzimuthal::evaluate(Direction u) const
{
  double mu = u.dot(u_ref_);
  double phi = std::acos(u.dot(v_ref_) / std::sqrt(1 - mu * mu));
  return mu_->evaluate(mu) * phi_->evaluate(phi);
}

//==============================================================================
// Isotropic implementation
//==============================================================================

Direction isotropic_direction(uint64_t* seed)
{
  double phi = uniform_distribution(0., 2.0 * PI, seed);
  double mu = uniform_distribution(-1., 1., seed);
  return {mu, std::sqrt(1.0 - mu * mu) * std::cos(phi),
    std::sqrt(1.0 - mu * mu) * std::sin(phi)};
}

Direction Isotropic::sample(uint64_t* seed) const
{
  return isotropic_direction(seed);
}

double Isotropic::evaluate(Direction u) const
{
  return 1.0 / (4.0 * PI);
}

//==============================================================================
// Monodirectional implementation
//==============================================================================

Direction Monodirectional::sample(uint64_t* seed) const
{
  return u_ref_;
}

} // namespace openmc

1	#include "openmc/distribution_multi.h"
2
3	#include <algorithm> // for move
4	#include <cmath> // for sqrt, sin, cos, max
5
6	#include "openmc/constants.h"
7	#include "openmc/error.h"
8	#include "openmc/math_functions.h"
9	#include "openmc/random_dist.h"
10	#include "openmc/random_lcg.h"
11	#include "openmc/xml_interface.h"
12
13	namespace openmc {
14
15	unique_ptr<UnitSphereDistribution> UnitSphereDistribution::create(	3,164✔
16	pugi::xml_node node)
17	{
18	// Check for type of angular distribution
19	std::string type;	3,164✔
20	if (check_for_node(node, "type"))	3,164!
21	type = get_node_value(node, "type", true, true);	3,164✔
22	if (type == "isotropic") {	3,164✔
23	return UPtrAngle {new Isotropic()};	256✔
24	} else if (type == "monodirectional") {	2,908✔
25	return UPtrAngle {new Monodirectional(node)};	2,870✔
26	} else if (type == "mu-phi") {	38!
27	return UPtrAngle {new PolarAzimuthal(node)};	38✔
28	} else {
29	fatal_error(fmt::format(	×
30	"Invalid angular distribution for external source: {}", type));
31	}
32	}	3,164✔
33
34	//==============================================================================
35	// UnitSphereDistribution implementation
36	//==============================================================================
37
38	UnitSphereDistribution::UnitSphereDistribution(pugi::xml_node node)	2,908✔
39	{
40	// Read reference directional unit vector
41	if (check_for_node(node, "reference_uvw")) {	2,908!
42	auto u_ref = get_node_array<double>(node, "reference_uvw");	2,908✔
43	if (u_ref.size() != 3)	2,908!
44	fatal_error("Angular distribution reference direction must have "	×
45	"three parameters specified.");
46	u_ref_ = Direction(u_ref.data());	2,908✔
47	u_ref_ /= u_ref_.norm();	2,908✔
48	}	2,908✔
49	}	2,908✔
50
51	//==============================================================================
52	// PolarAzimuthal implementation
53	//==============================================================================
54
55	PolarAzimuthal::PolarAzimuthal(Direction u, UPtrDist mu, UPtrDist phi)	×
56	: UnitSphereDistribution {u}, mu_ {std::move(mu)}, phi_ {std::move(phi)}	×
57	{}	×
58
59	PolarAzimuthal::PolarAzimuthal(pugi::xml_node node)	38✔
60	: UnitSphereDistribution {node}	38✔
61	{
62	// Read reference directional unit vector
63	if (check_for_node(node, "reference_vwu")) {	38!
64	auto v_ref = get_node_array<double>(node, "reference_vwu");	38✔
65	if (v_ref.size() != 3)	38!
66	fatal_error("Angular distribution reference v direction must have "	×
67	"three parameters specified.");
68	v_ref_ = Direction(v_ref.data());	38✔
69	v_ref_ /= v_ref_.norm();	38✔
70	}	38✔
71	w_ref_ = u_ref_.cross(v_ref_);	38✔
72	if (check_for_node(node, "mu")) {	38!
73	pugi::xml_node node_dist = node.child("mu");	38✔
74	mu_ = distribution_from_xml(node_dist);	38✔
75	} else {
76	mu_ = UPtrDist {new Uniform(-1., 1.)};	×
77	}
78
79	if (check_for_node(node, "phi")) {	38!
80	pugi::xml_node node_dist = node.child("phi");	38✔
81	phi_ = distribution_from_xml(node_dist);	38✔
82	} else {
83	phi_ = UPtrDist {new Uniform(0.0, 2.0 * PI)};	×
84	}
85	}	38✔
86
87	Direction PolarAzimuthal::sample(uint64_t* seed) const	25,300✔
88	{
89	// Sample cosine of polar angle
90	double mu = mu_->sample(seed);	25,300✔
91	if (mu == 1.0)	25,300✔
92	return u_ref_;	759✔
93	if (mu == -1.0)	24,541✔
94	return -u_ref_;	1,628✔
95
96	// Sample azimuthal angle
97	double phi = phi_->sample(seed);	22,913✔
98
99	double f = std::sqrt(1 - mu * mu);	22,913✔
100
101	return mu * u_ref_ + f * std::cos(phi) * v_ref_ + f * std::sin(phi) * w_ref_;	22,913✔
102	}
103
NEW 104	double PolarAzimuthal::evaluate(Direction u) const	×
105	{
NEW 106	double mu = u.dot(u_ref_);	×
NEW 107	double phi = std::acos(u.dot(v_ref_) / std::sqrt(1 - mu * mu));	×
NEW 108	return mu_->evaluate(mu) * phi_->evaluate(phi);	×
109	}
110
111	//==============================================================================
112	// Isotropic implementation
113	//==============================================================================
114
115	Direction isotropic_direction(uint64_t* seed)	112,895,415✔
116	{
117	double phi = uniform_distribution(0., 2.0 * PI, seed);	112,895,415✔
118	double mu = uniform_distribution(-1., 1., seed);	112,895,415✔
119	return {mu, std::sqrt(1.0 - mu * mu) * std::cos(phi),	112,895,415✔
120	std::sqrt(1.0 - mu * mu) * std::sin(phi)};	112,895,415✔
121	}
122
123	Direction Isotropic::sample(uint64_t* seed) const	19,463,958✔
124	{
125	return isotropic_direction(seed);	19,463,958✔
126	}
127
NEW 128	double Isotropic::evaluate(Direction u) const	×
129	{
NEW 130	return 1.0 / (4.0 * PI);	×
131	}
132
133	//==============================================================================
134	// Monodirectional implementation
135	//==============================================================================
136
137	Direction Monodirectional::sample(uint64_t* seed) const	12,465,854✔
138	{
139	return u_ref_;	12,465,854✔
140	}
141
142	} // namespace openmc

openmc-dev / openmc / 20470828925

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