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.9

/src/distribution_angle.cpp

#include "openmc/distribution_angle.h"

#include <cmath> // for abs, copysign

#include "xtensor/xarray.hpp"
#include "xtensor/xview.hpp"

#include "openmc/endf.h"
#include "openmc/hdf5_interface.h"
#include "openmc/math_functions.h"
#include "openmc/random_lcg.h"
#include "openmc/search.h"
#include "openmc/vector.h" // for vector

namespace openmc {

//==============================================================================
// AngleDistribution implementation
//==============================================================================

AngleDistribution::AngleDistribution(hid_t group)
{
  // Get incoming energies
  read_dataset(group, "energy", energy_);
  int n_energy = energy_.size();

  // Get outgoing energy distribution data
  vector<int> offsets;
  vector<int> interp;
  hid_t dset = open_dataset(group, "mu");
  read_attribute(dset, "offsets", offsets);
  read_attribute(dset, "interpolation", interp);
  xt::xarray<double> temp;
  read_dataset(dset, temp);
  close_dataset(dset);

  for (int i = 0; i < n_energy; ++i) {
    // Determine number of outgoing energies
    int j = offsets[i];
    int n;
    if (i < n_energy - 1) {
      n = offsets[i + 1] - j;
    } else {
      n = temp.shape()[1] - j;
    }

    // Create and initialize tabular distribution
    auto xs = xt::view(temp, 0, xt::range(j, j + n));
    auto ps = xt::view(temp, 1, xt::range(j, j + n));
    auto cs = xt::view(temp, 2, xt::range(j, j + n));
    vector<double> x {xs.begin(), xs.end()};
    vector<double> p {ps.begin(), ps.end()};
    vector<double> c {cs.begin(), cs.end()};

    // To get answers that match ACE data, for now we still use the tabulated
    // CDF values that were passed through to the HDF5 library. At a later
    // time, we can remove the CDF values from the HDF5 library and
    // reconstruct them using the PDF
    Tabular* mudist =
      new Tabular {x.data(), p.data(), n, int2interp(interp[i]), c.data()};

    distribution_.emplace_back(mudist);
  }
}

double AngleDistribution::sample(double E, uint64_t* seed) const
{
  // Find energy bin and calculate interpolation factor
  int i;
  double r;
  get_energy_index(energy_, E, i, r);

  // Sample between the ith and (i+1)th bin
  if (r > prn(seed))
    ++i;

  // Sample i-th distribution
  double mu = distribution_[i]->sample(seed);

  // Make sure mu is in range [-1,1] and return
  if (std::abs(mu) > 1.0)
    mu = std::copysign(1.0, mu);
  return mu;
}

double AngleDistribution::evaluate(double E, double mu) const
{
  // Find energy bin and calculate interpolation factor
  int i;
  double r;
  get_energy_index(energy_, E, i, r);

  double pdf = 0.0;
  if (r > 0.0)
    pdf += r * distribution_[i + 1]->evaluate(mu);
  if (r < 1.0)
    pdf += (1.0 - r) * distribution_[i]->evaluate(mu);
  return pdf;
}

} // namespace openmc

1	#include "openmc/distribution_angle.h"
2
3	#include <cmath> // for abs, copysign
4
5	#include "xtensor/xarray.hpp"
6	#include "xtensor/xview.hpp"
7
8	#include "openmc/endf.h"
9	#include "openmc/hdf5_interface.h"
10	#include "openmc/math_functions.h"
11	#include "openmc/random_lcg.h"
12	#include "openmc/search.h"
13	#include "openmc/vector.h" // for vector
14
15	namespace openmc {
16
17	//==============================================================================
18	// AngleDistribution implementation
19	//==============================================================================
20
21	AngleDistribution::AngleDistribution(hid_t group)	2,963,266✔
22	{
23	// Get incoming energies
24	read_dataset(group, "energy", energy_);	2,963,266✔
25	int n_energy = energy_.size();	2,963,266✔
26
27	// Get outgoing energy distribution data
28	vector<int> offsets;	2,963,266✔
29	vector<int> interp;	2,963,266✔
30	hid_t dset = open_dataset(group, "mu");	2,963,266✔
31	read_attribute(dset, "offsets", offsets);	2,963,266✔
32	read_attribute(dset, "interpolation", interp);	2,963,266✔
33	xt::xarray<double> temp;	2,963,266✔
34	read_dataset(dset, temp);	2,963,266✔
35	close_dataset(dset);	2,963,266✔
36
37	for (int i = 0; i < n_energy; ++i) {	28,475,883✔
38	// Determine number of outgoing energies
39	int j = offsets[i];	25,512,617✔
40	int n;
41	if (i < n_energy - 1) {	25,512,617✔
42	n = offsets[i + 1] - j;	22,549,351✔
43	} else {
44	n = temp.shape()[1] - j;	2,963,266✔
45	}
46
47	// Create and initialize tabular distribution
48	auto xs = xt::view(temp, 0, xt::range(j, j + n));	25,512,617✔
49	auto ps = xt::view(temp, 1, xt::range(j, j + n));	25,512,617✔
50	auto cs = xt::view(temp, 2, xt::range(j, j + n));	25,512,617✔
51	vector<double> x {xs.begin(), xs.end()};	25,512,617✔
52	vector<double> p {ps.begin(), ps.end()};	25,512,617✔
53	vector<double> c {cs.begin(), cs.end()};	25,512,617✔
54
55	// To get answers that match ACE data, for now we still use the tabulated
56	// CDF values that were passed through to the HDF5 library. At a later
57	// time, we can remove the CDF values from the HDF5 library and
58	// reconstruct them using the PDF
59	Tabular* mudist =
60	new Tabular {x.data(), p.data(), n, int2interp(interp[i]), c.data()};	25,512,617✔
61
62	distribution_.emplace_back(mudist);	25,512,617✔
63	}	25,512,617✔
64	}	2,963,266✔
65
66	double AngleDistribution::sample(double E, uint64_t* seed) const	691,691,651✔
67	{
68	// Find energy bin and calculate interpolation factor
69	int i;
70	double r;
71	get_energy_index(energy_, E, i, r);	691,691,651✔
72
73	// Sample between the ith and (i+1)th bin
74	if (r > prn(seed))	691,691,651✔
75	++i;	280,574,831✔
76
77	// Sample i-th distribution
78	double mu = distribution_[i]->sample(seed);	691,691,651✔
79
80	// Make sure mu is in range [-1,1] and return
81	if (std::abs(mu) > 1.0)	691,691,651!
82	mu = std::copysign(1.0, mu);	×
83	return mu;	691,691,651✔
84	}
85
NEW 86	double AngleDistribution::evaluate(double E, double mu) const	×
87	{
88	// Find energy bin and calculate interpolation factor
89	int i;
90	double r;
NEW 91	get_energy_index(energy_, E, i, r);	×
92
NEW 93	double pdf = 0.0;	×
NEW 94	if (r > 0.0)	×
NEW 95	pdf += r * distribution_[i + 1]->evaluate(mu);	×
NEW 96	if (r < 1.0)	×
NEW 97	pdf += (1.0 - r) * distribution_[i]->evaluate(mu);	×
NEW 98	return pdf;	×
99	}
100
101	} // 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