crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 13:45:36 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Move ParticleType into Particle and use for type_ member

Browse source
This commit is contained in:
Paul Romano 2019-02-27 09:00:59 -06:00
parent b4ed267d4b
commit b28ee8087c
21 changed files with 104 additions and 97 deletions
Download patch file Download diff file Expand all files Collapse all files

14
include/openmc/particle.h
View file

@ -33,11 +33,6 @@ constexpr int MAX_LOST_PARTICLES {10};
// Maximum number of lost particles, relative to the total number of particles
constexpr double REL_MAX_LOST_PARTICLES {1.0e-6};
//! Particle types
enum class ParticleType {
neutron, photon, electron, positron
};
struct LocalCoord {
int cell {-1};
int universe {-1};
@ -59,8 +54,13 @@ struct LocalCoord {
class Particle {
public:
//! Particle types
enum class Type {
neutron, photon, electron, positron
};
int64_t id_; //!< Unique ID
int type_; //!< Particle type (n, p, e, etc.)
Type type_; //!< Particle type (n, p, e, etc.)
int n_coord_; //!< number of current coordinate levels
int cell_instance_; //!< offset for distributed properties
@ -134,7 +134,7 @@ public:
//! \param E Energy of the secondary particle in [eV]
//! \param type Particle type
//! \param run_CE Whether continuous-energy data is being used
void create_secondary(const double* uvw, double E, int type, bool run_CE);
void create_secondary(const double* uvw, double E, Type type, bool run_CE);
//! sets default attributes for a particle
void initialize();

2
include/openmc/reaction_product.h
View file

@ -44,7 +44,7 @@ public:
//! \param[out] mu Outgoing cosine with respect to current direction
void sample(double E_in, double& E_out, double& mu) const;
ParticleType particle_; //!< Particle type
Particle::Type particle_; //!< Particle type
EmissionMode emission_mode_; //!< Emission mode
double decay_rate_; //!< Decay rate (for delayed neutron precursors) in [1/s]
std::unique_ptr<Function1D> yield_; //!< Yield as a function of energy

2
include/openmc/source.h
View file

@ -45,7 +45,7 @@ public:
// Properties
double strength() const { return strength_; }
private:
ParticleType particle_ {ParticleType::neutron}; //!< Type of particle emitted
Particle::Type particle_ {Particle::Type::neutron}; //!< Type of particle emitted
double strength_ {1.0}; //!< Source strength
UPtrSpace space_; //!< Spatial distribution
UPtrAngle angle_; //!< Angular distribution

3
include/openmc/tallies/filter_particle.h
View file

@ -3,6 +3,7 @@
#include <vector>
#include "openmc/particle.h"
#include "openmc/tallies/filter.h"
namespace openmc {
@ -27,7 +28,7 @@ public:
std::string text_label(int bin) const override;
std::vector<int> particles_;
std::vector<Particle::Type> particles_;
};
} // namespace openmc

7
src/bremsstrahlung.cpp
View file

@ -30,12 +30,12 @@ void thick_target_bremsstrahlung(Particle& p, double* E_lost)
{
if (p.material_ == MATERIAL_VOID) return;
int photon = static_cast<int>(ParticleType::photon);
int photon = static_cast<int>(Particle::Type::photon);
if (p.E_ < settings::energy_cutoff[photon]) return;
// Get bremsstrahlung data for this material and particle type
BremsstrahlungData* mat;
if (p.type_ == static_cast<int>(ParticleType::positron)) {
if (p.type_ == Particle::Type::positron) {
mat = &model::materials[p.material_]->ttb_->positron;
} else {
mat = &model::materials[p.material_]->ttb_->electron;
@ -108,8 +108,7 @@ void thick_target_bremsstrahlung(Particle& p, double* E_lost)
if (w > settings::energy_cutoff[photon]) {
// Create secondary photon
int photon_ = static_cast<int>(ParticleType::photon);
p.create_secondary(p.coord_[0].uvw, w, photon_, true);
p.create_secondary(p.coord_[0].uvw, w, Particle::Type::photon, true);
*E_lost += w;
}
}

10
src/cross_sections.cpp
View file

@ -228,7 +228,7 @@ read_ce_cross_sections(const std::vector<std::vector<double>>& nuc_temps,
// Determine if minimum/maximum energy for this nuclide is greater/less
// than the previous
if (data::nuclides[i_nuclide]->grid_.size() >= 1) {
int neutron = static_cast<int>(ParticleType::neutron);
int neutron = static_cast<int>(Particle::Type::neutron);
data::energy_min[neutron] = std::max(data::energy_min[neutron],
data::nuclides[i_nuclide]->grid_[0].energy.front());
data::energy_max[neutron] = std::min(data::energy_max[neutron],
@ -262,7 +262,7 @@ read_ce_cross_sections(const std::vector<std::vector<double>>& nuc_temps,
// the previous
const auto& elem {data::elements.back()};
if (elem.energy_.size() >= 1) {
int photon = static_cast<int>(ParticleType::photon);
int photon = static_cast<int>(Particle::Type::photon);
int n = elem.energy_.size();
data::energy_min[photon] = std::max(data::energy_min[photon],
std::exp(elem.energy_(1)));
@ -321,7 +321,7 @@ read_ce_cross_sections(const std::vector<std::vector<double>>& nuc_temps,
for (auto& nuc : data::nuclides) {
nuc->init_grid();
}
int neutron = static_cast<int>(ParticleType::neutron);
int neutron = static_cast<int>(Particle::Type::neutron);
simulation::log_spacing = std::log(data::energy_max[neutron] /
data::energy_min[neutron]) / settings::n_log_bins;
@ -329,7 +329,7 @@ read_ce_cross_sections(const std::vector<std::vector<double>>& nuc_temps,
// Determine if minimum/maximum energy for bremsstrahlung is greater/less
// than the current minimum/maximum
if (data::ttb_e_grid.size() >= 1) {
int photon = static_cast<int>(ParticleType::photon);
int photon = static_cast<int>(Particle::Type::photon);
int n_e = data::ttb_e_grid.size();
data::energy_min[photon] = std::max(data::energy_min[photon], data::ttb_e_grid(1));
data::energy_max[photon] = std::min(data::energy_max[photon], data::ttb_e_grid(n_e - 1));
@ -345,7 +345,7 @@ read_ce_cross_sections(const std::vector<std::vector<double>>& nuc_temps,
// grid has not been allocated
if (nuc->grid_.size() > 0) {
double max_E = nuc->grid_[0].energy.back();
int neutron = static_cast<int>(ParticleType::neutron);
int neutron = static_cast<int>(Particle::Type::neutron);
if (max_E == data::energy_max[neutron]) {
write_message("Maximum neutron transport energy: " +
std::to_string(data::energy_max[neutron]) + " eV for " +

7
src/material.cpp
View file

@ -652,18 +652,17 @@ void Material::calculate_xs(const Particle& p) const
simulation::material_xs.fission = 0.0;
simulation::material_xs.nu_fission = 0.0;
if (p.type_ == static_cast<int>(ParticleType::neutron)) {
if (p.type_ == Particle::Type::neutron) {
this->calculate_neutron_xs(p);
} else if (p.type_ == static_cast<int>(ParticleType::photon)) {
} else if (p.type_ == Particle::Type::photon) {
this->calculate_photon_xs(p);
}
}
void Material::calculate_neutron_xs(const Particle& p) const
{
int neutron = static_cast<int>(ParticleType::neutron);
// Find energy index on energy grid
int neutron = static_cast<int>(Particle::Type::neutron);
int i_grid = std::log(p.E_/data::energy_min[neutron])/simulation::log_spacing;
// Determine if this material has S(a,b) tables

2
src/mgxs_interface.cpp
View file

@ -235,7 +235,7 @@ void read_mg_cross_sections_header()
}
// Get the minimum and maximum energies
int neutron = static_cast<int>(ParticleType::neutron);
int neutron = static_cast<int>(Particle::Type::neutron);
data::energy_min[neutron] = data::energy_bins.back();
data::energy_max[neutron] = data::energy_bins.front();

6
src/nuclide.cpp
View file

@ -287,7 +287,7 @@ void Nuclide::create_derived()
auto xs = xt::adapt(rx->xs_[t].value);
for (const auto& p : rx->products_) {
if (p.particle_ == ParticleType::photon) {
if (p.particle_ == Particle::Type::photon) {
auto pprod = xt::view(xs_[t], xt::range(j, j+n), XS_PHOTON_PROD);
for (int k = 0; k < n; ++k) {
double E = grid_[t].energy[k+j];
@ -391,7 +391,7 @@ void Nuclide::create_derived()
void Nuclide::init_grid()
{
int neutron = static_cast<int>(ParticleType::neutron);
int neutron = static_cast<int>(Particle::Type::neutron);
double E_min = data::energy_min[neutron];
double E_max = data::energy_max[neutron];
int M = settings::n_log_bins;
@ -440,7 +440,7 @@ double Nuclide::nu(double E, EmissionMode mode, int group) const
for (int i = 1; i < rx->products_.size(); ++i) {
// Skip any non-neutron products
const auto& product = rx->products_[i];
if (product.particle_ != ParticleType::neutron) continue;
if (product.particle_ != Particle::Type::neutron) continue;
// Evaluate yield
if (product.emission_mode_ == EmissionMode::delayed) {

8
src/output.cpp
View file

@ -143,16 +143,16 @@ extern "C" void print_particle(Particle* p)
{
// Display particle type and ID.
switch (p->type_) {
case static_cast<int>(ParticleType::neutron):
case Particle::Type::neutron:
std::cout << "Neutron ";
break;
case static_cast<int>(ParticleType::photon):
case Particle::Type::photon:
std::cout << "Photon ";
break;
case static_cast<int>(ParticleType::electron):
case Particle::Type::electron:
std::cout << "Electron ";
break;
case static_cast<int>(ParticleType::positron):
case Particle::Type::positron:
std::cout << "Positron ";
break;
default:

26
src/particle.cpp
View file

@ -55,18 +55,18 @@ Particle::clear()
}
void
Particle::create_secondary(const double* uvw, double E, int type, bool run_CE)
Particle::create_secondary(const double* uvw, double E, Type type, bool run_CE)
{
if (n_secondary_ == MAX_SECONDARY) {
fatal_error("Too many secondary particles created.");
}
int64_t n = n_secondary_;
secondary_bank_[n].particle = type_;
secondary_bank_[n].particle = static_cast<int>(type);
secondary_bank_[n].wgt = wgt_;
std::copy(coord_[0].xyz, coord_[0].xyz + 3, secondary_bank_[n].xyz);
std::copy(uvw, uvw + 3, secondary_bank_[n].uvw);
secondary_bank_[n].E = E_;
secondary_bank_[n].E = E;
if (!run_CE) secondary_bank_[n].E = g_;
++n_secondary_;
}
@ -78,7 +78,7 @@ Particle::initialize()
clear();
// Set particle to neutron that's alive
type_ = static_cast<int>(ParticleType::neutron);
type_ = Particle::Type::neutron;
alive_ = true;
// clear attributes
@ -114,9 +114,9 @@ Particle::from_source(const Bank* src)
initialize();
// copy attributes from source bank site
type_ = src->particle;
wgt_ = src->wgt;
last_wgt_ = src->wgt;
type_ = static_cast<Particle::Type>(src->particle);
wgt_ = src->wgt;
last_wgt_ = src->wgt;
std::copy(src->xyz, src->xyz + 3, coord_[0].xyz);
std::copy(src->uvw, src->uvw + 3, coord_[0].uvw);
std::copy(src->xyz, src->xyz + 3, last_xyz_current_);
@ -163,7 +163,7 @@ Particle::transport()
while (true) {
// Set the random number stream
if (type_ == static_cast<int>(ParticleType::neutron)) {
if (type_ == Particle::Type::neutron) {
prn_set_stream(STREAM_TRACKING);
} else {
prn_set_stream(STREAM_PHOTON);
@ -230,8 +230,8 @@ Particle::transport()
// Sample a distance to collision
double d_collision;
if (type_ == static_cast<int>(ParticleType::electron) ||
type_ == static_cast<int>(ParticleType::positron)) {
if (type_ == Particle::Type::electron ||
type_ == Particle::Type::positron) {
d_collision = 0.0;
} else if (simulation::material_xs.total == 0.0) {
d_collision = INFINITY;
@ -257,7 +257,7 @@ Particle::transport()
// Score track-length estimate of k-eff
if (settings::run_mode == RUN_MODE_EIGENVALUE &&
type_ == static_cast<int>(ParticleType::neutron)) {
type_ == Particle::Type::neutron) {
global_tally_tracklength += wgt_ * distance * simulation::material_xs.nu_fission;
}
@ -300,7 +300,7 @@ Particle::transport()
// Score collision estimate of keff
if (settings::run_mode == RUN_MODE_EIGENVALUE &&
type_ == static_cast<int>(ParticleType::neutron)) {
type_ == Particle::Type::neutron) {
global_tally_collision += wgt_ * simulation::material_xs.nu_fission
/ simulation::material_xs.total;
}
@ -687,7 +687,7 @@ Particle::write_restart() const
break;
}
write_dataset(file_id, "id", id_);
write_dataset(file_id, "type", type_);
write_dataset(file_id, "type", static_cast<int>(type_));
int64_t i = simulation::current_work;
write_dataset(file_id, "weight", simulation::source_bank[i-1].wgt);

4
src/particle_restart.cpp
View file

@ -40,7 +40,9 @@ void read_particle_restart(Particle& p, int& previous_run_mode)
previous_run_mode = RUN_MODE_FIXEDSOURCE;
}
read_dataset(file_id, "id", p.id_);
read_dataset(file_id, "type", p.type_);
int type;
read_dataset(file_id, "type", type);
p.type_ = static_cast<Particle::Type>(type);
read_dataset(file_id, "weight", p.wgt_);
read_dataset(file_id, "energy", p.E_);
std::array<double, 3> x;

11
src/photon.cpp
View file

@ -224,7 +224,7 @@ PhotonInteraction::PhotonInteraction(hid_t group, int i_element)
}
// Truncate the bremsstrahlung data at the cutoff energy
int photon = static_cast<int>(ParticleType::photon);
int photon = static_cast<int>(Particle::Type::photon);
const auto& E {electron_energy};
double cutoff = settings::energy_cutoff[photon];
if (cutoff > E(0)) {
@ -658,8 +658,7 @@ void PhotonInteraction::atomic_relaxation(const ElectronSubshell& shell, Particl
uvw[1] = std::sqrt(1.0 - mu*mu)*std::cos(phi);
uvw[2] = std::sqrt(1.0 - mu*mu)*std::sin(phi);
double E = shell.binding_energy;
int photon = static_cast<int>(ParticleType::photon);
p.create_secondary(uvw.data(), E, photon, true);
p.create_secondary(uvw.data(), E, Particle::Type::photon, true);
return;
}
@ -691,8 +690,7 @@ void PhotonInteraction::atomic_relaxation(const ElectronSubshell& shell, Particl
// Non-radiative transition -- Auger/Coster-Kronig effect
// Create auger electron
int electron = static_cast<int>(ParticleType::electron);
p.create_secondary(uvw.data(), E, electron, true);
p.create_secondary(uvw.data(), E, Particle::Type::electron, true);
// Fill hole left by emitted auger electron
int i_hole = shell_map_.at(secondary);
@ -702,8 +700,7 @@ void PhotonInteraction::atomic_relaxation(const ElectronSubshell& shell, Particl
// Radiative transition -- get X-ray energy
// Create fluorescent photon
int photon = static_cast<int>(ParticleType::photon);
p.create_secondary(uvw.data(), E, photon, true);
p.create_secondary(uvw.data(), E, Particle::Type::photon, true);
}
// Fill hole created by electron transitioning to the photoelectron hole

48
src/physics.cpp
View file

@ -36,23 +36,24 @@ void collision(Particle* p)
++(p->n_collision_);
// Sample reaction for the material the particle is in
switch (static_cast<ParticleType>(p->type_)) {
case ParticleType::neutron:
switch (static_cast<Particle::Type>(p->type_)) {
case Particle::Type::neutron:
sample_neutron_reaction(p);
break;
case ParticleType::photon:
case Particle::Type::photon:
sample_photon_reaction(p);
break;
case ParticleType::electron:
case Particle::Type::electron:
sample_electron_reaction(p);
break;
case ParticleType::positron:
case Particle::Type::positron:
sample_positron_reaction(p);
break;
}
// Kill particle if energy falls below cutoff
if (p->E_ < settings::energy_cutoff[p->type_]) {
int type = static_cast<int>(p->type_);
if (p->E_ < settings::energy_cutoff[type]) {
p->alive_ = false;
p->wgt_ = 0.0;
p->last_wgt_ = 0.0;
@ -61,7 +62,7 @@ void collision(Particle* p)
// Display information about collision
if (settings::verbosity >= 10 || simulation::trace) {
std::stringstream msg;
if (static_cast<ParticleType>(p->type_) == ParticleType::neutron) {
if (p->type_ == Particle::Type::neutron) {
msg << " " << reaction_name(p->event_mt_) << " with " <<
data::nuclides[p->event_nuclide_]->name_ << ". Energy = " << p->E_ << " eV.";
} else {
@ -185,7 +186,7 @@ create_fission_sites(Particle* p, int i_nuclide, const Reaction* rx, Bank* bank_
bank_array[i].xyz[2] = p->coord_[0].xyz[2];
// Set that the bank particle is a neutron
bank_array[i].particle = static_cast<int>(ParticleType::neutron);
bank_array[i].particle = static_cast<int>(Particle::Type::neutron);
// Set the weight of the fission bank site
bank_array[i].wgt = 1. / weight;
@ -218,7 +219,7 @@ void sample_photon_reaction(Particle* p)
// Kill photon if below energy cutoff -- an extra check is made here because
// photons with energy below the cutoff may have been produced by neutrons
// reactions or atomic relaxation
int photon = static_cast<int>(ParticleType::photon);
int photon = static_cast<int>(Particle::Type::photon);
if (p->E_ < settings::energy_cutoff[photon]) {
p->E_ = 0.0;
p->alive_ = false;
@ -272,8 +273,7 @@ void sample_photon_reaction(Particle* p)
double uvw[3];
std::copy(p->coord_[0].uvw, p->coord_[0].uvw + 3, uvw);
rotate_angle_c(uvw, mu_electron, &phi);
int electron = static_cast<int>(ParticleType::electron);
p->create_secondary(uvw, E_electron, electron, true);
p->create_secondary(uvw, E_electron, Particle::Type::electron, true);
// TODO: Compton subshell data does not match atomic relaxation data
// Allow electrons to fill orbital and produce auger electrons
@ -332,8 +332,7 @@ void sample_photon_reaction(Particle* p)
uvw[2] = std::sqrt(1.0 - mu*mu)*std::sin(phi);
// Create secondary electron
int electron = static_cast<int>(ParticleType::electron);
p->create_secondary(uvw.data(), E_electron, electron, true);
p->create_secondary(uvw.data(), E_electron, Particle::Type::electron, true);
// Allow electrons to fill orbital and produce auger electrons
// and fluorescent photons
@ -359,14 +358,12 @@ void sample_photon_reaction(Particle* p)
double uvw[3];
std::copy(p->coord_[0].uvw, p->coord_[0].uvw + 3, uvw);
rotate_angle_c(uvw, mu_electron, nullptr);
int electron = static_cast<int>(ParticleType::electron);
p->create_secondary(uvw, E_electron, electron, true);
p->create_secondary(uvw, E_electron, Particle::Type::electron, true);
// Create secondary positron
std::copy(p->coord_[0].uvw, p->coord_[0].uvw + 3, uvw);
rotate_angle_c(uvw, mu_positron, nullptr);
int positron = static_cast<int>(ParticleType::positron);
p->create_secondary(uvw, E_positron, positron, true);
p->create_secondary(uvw, E_positron, Particle::Type::positron, true);
p->event_mt_ = PAIR_PROD;
p->alive_ = false;
@ -405,13 +402,12 @@ void sample_positron_reaction(Particle* p)
uvw[2] = std::sqrt(1.0 - mu*mu)*std::sin(phi);
// Create annihilation photon pair traveling in opposite directions
int photon = static_cast<int>(ParticleType::photon);
p->create_secondary(uvw.data(), MASS_ELECTRON_EV, photon, true);
p->create_secondary(uvw.data(), MASS_ELECTRON_EV, Particle::Type::photon, true);
uvw[0] = -uvw[0];
uvw[1] = -uvw[1];
uvw[2] = -uvw[2];
p->create_secondary(uvw.data(), MASS_ELECTRON_EV, photon, true);
p->create_secondary(uvw.data(), MASS_ELECTRON_EV, Particle::Type::photon, true);
p->E_ = 0.0;
p->alive_ = false;
@ -541,7 +537,7 @@ void sample_photon_product(int i_nuclide, double E, int* i_rx, int* i_product)
+ f*(rx->xs_[i_temp].value[i_grid - threshold + 1]));
for (int j = 0; j < rx->products_.size(); ++j) {
if (rx->products_[j].particle_ == ParticleType::photon) {
if (rx->products_[j].particle_ == Particle::Type::photon) {
// add to cumulative probability
prob += (*rx->products_[j].yield_)(E) * xs;
@ -1029,7 +1025,7 @@ void sample_fission_neutron(int i_nuclide, const Reaction* rx, double E_in, Bank
rx->products_[group].sample(E_in, site->E, mu);
// resample if energy is greater than maximum neutron energy
constexpr int neutron = static_cast<int>(ParticleType::neutron);
constexpr int neutron = static_cast<int>(Particle::Type::neutron);
if (site->E < data::energy_max[neutron]) break;
// check for large number of resamples
@ -1054,7 +1050,7 @@ void sample_fission_neutron(int i_nuclide, const Reaction* rx, double E_in, Bank
rx->products_[0].sample(E_in, site->E, mu);
// resample if energy is greater than maximum neutron energy
constexpr int neutron = static_cast<int>(ParticleType::neutron);
constexpr int neutron = static_cast<int>(Particle::Type::neutron);
if (site->E < data::energy_max[neutron]) break;
// check for large number of resamples
@ -1109,8 +1105,7 @@ void inelastic_scatter(const Nuclide* nuc, const Reaction* rx, Particle* p)
if (std::floor(yield) == yield) {
// If yield is integral, create exactly that many secondary particles
for (int i = 0; i < static_cast<int>(std::round(yield)) - 1; ++i) {
int neutron = static_cast<int>(ParticleType::neutron);
p->create_secondary(p->coord_[0].uvw, p->E_, neutron, true);
p->create_secondary(p->coord_[0].uvw, p->E_, Particle::Type::neutron, true);
}
} else {
// Otherwise, change weight of particle based on yield
@ -1145,8 +1140,7 @@ void sample_secondary_photons(Particle* p, int i_nuclide)
rotate_angle_c(uvw, mu, nullptr);
// Create the secondary photon
int photon = static_cast<int>(ParticleType::photon);
p->create_secondary(uvw, E, photon, true);
p->create_secondary(uvw, E, Particle::Type::photon, true);
}
}

2
src/physics_mg.cpp
View file

@ -157,7 +157,7 @@ create_fission_sites(Particle* p, Bank* bank_array, int64_t* size_bank,
bank_array[i].xyz[2] = p->coord_[0].xyz[2];
// Set that the bank particle is a neutron
bank_array[i].particle = static_cast<int>(ParticleType::neutron);
bank_array[i].particle = static_cast<int>(Particle::Type::neutron);
// Set the weight of the fission bank site
bank_array[i].wgt = 1. / weight;

2
src/reaction.cpp
View file

@ -72,7 +72,7 @@ Reaction::Reaction(hid_t group, const std::vector<int>& temperatures)
// mark fission reactions so that we avoid the angle sampling.
if (is_fission(mt_)) {
for (auto& p : products_) {
if (p.particle_ == ParticleType::neutron) {
if (p.particle_ == Particle::Type::neutron) {
for (auto& d : p.distribution_) {
auto d_ = dynamic_cast<UncorrelatedAngleEnergy*>(d.get());
if (d_) d_->fission() = true;

4
src/reaction_product.cpp
View file

@ -22,9 +22,9 @@ ReactionProduct::ReactionProduct(hid_t group)
std::string temp;
read_attribute(group, "particle", temp);
if (temp == "neutron") {
particle_ = ParticleType::neutron;
particle_ = Particle::Type::neutron;
} else if (temp == "photon") {
particle_ = ParticleType::photon;
particle_ = Particle::Type::photon;
}
// Read emission mode and decay rate

4
src/source.cpp
View file

@ -47,9 +47,9 @@ SourceDistribution::SourceDistribution(pugi::xml_node node)
if (check_for_node(node, "particle")) {
auto temp_str = get_node_value(node, "particle", true, true);
if (temp_str == "neutron") {
particle_ = ParticleType::neutron;
particle_ = Particle::Type::neutron;
} else if (temp_str == "photon") {
particle_ = ParticleType::photon;
particle_ = Particle::Type::photon;
settings::photon_transport = true;
} else {
fatal_error(std::string("Unknown source particle type: ") + temp_str);

23
src/tallies/filter_particle.cpp
View file

@ -7,8 +7,10 @@ namespace openmc {
void
ParticleFilter::from_xml(pugi::xml_node node)
{
particles_ = get_node_array<int>(node, "bins");
for (auto& p : particles_) --p;
auto particles = get_node_array<int>(node, "bins");
for (auto& p : particles) {
particles_.push_back(static_cast<Particle::Type>(p - 1));
}
n_bins_ = particles_.size();
}
@ -28,13 +30,26 @@ void
ParticleFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
write_dataset(filter_group, "bins", particles_);
std::vector<int> particles;
for (auto p : particles_) {
particles.push_back(static_cast<int>(p) + 1);
}
write_dataset(filter_group, "bins", particles);
}
std::string
ParticleFilter::text_label(int bin) const
{
return "Particle " + std::to_string(particles_[bin]);
switch (particles_[bin]) {
case Particle::Type::neutron:
return "Particle: neutron";
case Particle::Type::photon:
return "Particle: photon";
case Particle::Type::electron:
return "Particle: electron";
case Particle::Type::positron:
return "Particle: positron";
}
}
} // namespace openmc

12
src/tallies/tally.cpp
View file

@ -765,9 +765,9 @@ void read_tallies_xml()
} else {
const auto& f = model::tally_filters[particle_filter_index].get();
auto pf = dynamic_cast<ParticleFilter*>(f);
for (int p : pf->particles_) {
if (p == static_cast<int>(ParticleType::electron) ||
p == static_cast<int>(ParticleType::positron)) {
for (auto p : pf->particles_) {
if (p == Particle::Type::electron ||
p == Particle::Type::positron) {
t->estimator_ = ESTIMATOR_ANALOG;
}
}
@ -776,11 +776,11 @@ void read_tallies_xml()
if (particle_filter_index >= 0) {
const auto& f = model::tally_filters[particle_filter_index].get();
auto pf = dynamic_cast<ParticleFilter*>(f);
for (int p : pf->particles_) {
if (p != static_cast<int>(ParticleType::neutron)) {
for (auto p : pf->particles_) {
if (p != Particle::Type::neutron) {
warning("Particle filter other than NEUTRON used with photon "
"transport turned off. All tallies for particle type " +
std::to_string(p) + " will have no scores");
std::to_string(static_cast<int>(p)) + " will have no scores");
}
}
}

4
src/tallies/tally_scoring.cpp
View file

@ -346,8 +346,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
score = p->last_wgt_;
}
if (p->type_ == static_cast<int>(ParticleType::neutron) ||
p->type_ == static_cast<int>(ParticleType::photon)) {
if (p->type_ == Particle::Type::neutron ||
p->type_ == Particle::Type::photon) {
score *= flux / simulation::material_xs.total;
} else {
score = 0.;