crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-29 06:35:48 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

wip

Browse source
This commit is contained in:
GuySten 2026-05-10 17:44:50 +03:00
parent c98f0626f7
commit debc5a921a
4 changed files with 14 additions and 20 deletions
Download patch file Download diff file Expand all files Collapse all files

3
include/openmc/angle_energy.h
View file

@ -35,7 +35,8 @@ public:
virtual ~AngleEnergy() = default;
};
double get_jac_and_transform(double E_in, double& mu, double& E_out, uint64_t* seed, double awr);
double get_jac_and_transform(
double E_in, double& mu, double& E_out, uint64_t* seed, double awr);
double get_jac_and_transform(double E_in, double& mu, double& E_out,
uint64_t* seed, double awr, double E_com);

4
src/angle_energy.cpp
View file

@ -2,7 +2,6 @@
#include <cmath> // for sqrt
#include "openmc/error.h"
#include "openmc/random_lcg.h"
namespace openmc {
@ -40,8 +39,7 @@ double get_jac_and_transform(double E_in, double& mu, double& E_out,
}
}
mu = mu_lab * std::sqrt(E_out / E_cm) - std::sqrt(E_com / E_cm);
double jac = mult * E_out / (D * std::sqrt(E_cm * E_com));
return jac;
return mult * E_out / (D * std::sqrt(E_cm * E_com));
}
} // namespace openmc

10
src/physics.cpp
View file

@ -719,8 +719,6 @@ void scatter(Particle& p, int i_nuclide)
if (prob > cutoff && !sampled) {
// =======================================================================
// S(A,B) SCATTERING
fatal_error("moo oo 3");
sab_scatter(i_nuclide, micro.index_sab, p);
@ -798,11 +796,11 @@ void elastic_scatter(int i_nuclide, const Reaction& rx, double kT, Particle& p)
double mu_cm;
auto& d = rx.products_[0].distribution_[0];
auto d_ = dynamic_cast<UncorrelatedAngleEnergy*>(d.get());
//if (!d_->angle().empty()) {
// mu_cm = d_->angle().sample(p.E(), p.current_seed());
//} else {
if (!d_->angle().empty()) {
mu_cm = d_->angle().sample(p.E(), p.current_seed());
} else {
mu_cm = uniform_distribution(-1., 1., p.current_seed());
//}
}
// Determine direction cosines in CM
Direction u_cm = v_n / vel;

17
src/tallies/next_event_scoring.cpp
View file

@ -21,23 +21,23 @@ void score_point_tally_elastic(
double E_com = v_cm.dot(v_cm);
double E_out = (v_n - v_cm).dot(v_n - v_cm);
double E_in = p.E();
auto u_n = p.u();
auto u_cm = v_cm / v_cm.norm();
auto& d = rx.products_[i_product].distribution_[0];
auto d_ = dynamic_cast<UncorrelatedAngleEnergy*>(d.get());
auto pdf = [&](Direction u, double& E) {
double mu = u.dot(u_n);
double mu = u.dot(u_cm);
E = E_out;
double jac =
get_jac_and_transform(E_in, mu, E, p.current_seed(), awr, E_com);
//if (!d_->angle().empty()) {
// return jac * d_->angle().evaluate(p.E(), mu) / (2.0 * PI);
//} else {
if (!d_->angle().empty()) {
return jac * d_->angle().evaluate(p.E(), mu) / (2.0 * PI);
} else {
return jac * 0.5 / (2.0 * PI);
//}
}
};
score_point_tally_impl(p.r(), p.type(), p.time(), pdf);
score_point_tally_impl(p.r(), p.type(), p.time(), pdf);
}
void score_point_tally_inelastic(
@ -56,7 +56,6 @@ void score_point_tally_inelastic(
(2.0 * PI) * yield;
};
score_point_tally_impl(p.r(), p.type(), p.time(), pdf);
fatal_error("moo1");
}
void score_point_tally_fission(
@ -75,7 +74,6 @@ void score_point_tally_fission(
(2.0 * PI);
};
score_point_tally_impl(p.r(), p.type(), p.time(), pdf);
fatal_error("moo2");
}
void score_point_tally_sab(Particle& p, int i_nuclide, const ThermalData& sab,
@ -92,7 +90,6 @@ void score_point_tally_sab(Particle& p, int i_nuclide, const ThermalData& sab,
(2.0 * PI);
};
score_point_tally_impl(p.r(), p.type(), p.time(), pdf);
fatal_error("moo3");
}
void score_point_tally_source(SourceSite& site, int source_index)