Fix multipole and two other bugs

This commit is contained in:
Paul Romano 2018-12-21 14:30:28 -06:00
parent 620b71372b
commit 79ef612cb5
3 changed files with 21 additions and 5 deletions

View file

@ -193,6 +193,8 @@ extern "C" MaterialMacroXS material_xs;
//==============================================================================
extern "C" void set_micro_xs();
extern "C" void nuclide_multipole_eval(int i_nuclide, double E, double sqrtkT,
double* sig_s, double* sig_a, double* sig_f);
extern "C" bool nuclide_wmp_present(int i_nuclide);
extern "C" double nuclide_wmp_emin(int i_nuclide);
extern "C" double nuclide_wmp_emax(int i_nuclide);

View file

@ -520,7 +520,7 @@ void Nuclide::calculate_xs(int i_sab, double E, int i_log_union,
if (use_mp) {
// Call multipole kernel
double sig_s, sig_a, sig_f;
//multipole_eval(this % multipole, E, sqrtkT, sig_s, sig_a, sig_f)
nuclide_multipole_eval(i_nuclide_, E, sqrtkT, &sig_s, &sig_a, &sig_f);
micro_xs.total = sig_s + sig_a;
micro_xs.elastic = sig_s;
@ -651,11 +651,10 @@ void Nuclide::calculate_xs(int i_sab, double E, int i_log_union,
// If reaction is present and energy is greater than threshold, set the
// reaction xs appropriately
int i_rx = reaction_index_[DEPLETION_RX[j]];
const auto& rx = reactions_[i_rx];
const auto& rx_xs = rx->xs_[i_temp].value;
if (i_rx >= 0) {
const auto& rx = reactions_[i_rx];
const auto& rx_xs = rx->xs_[i_temp].value;
// Physics says that (n,gamma) is not a threshold reaction, so we don't
// need to specifically check its threshold index
if (j == 0) {

View file

@ -1048,6 +1048,9 @@ contains
! Add entry to nuclide dictionary
call nuclide_dict % set(to_lower(name_), n)
n_nuclides = n
! Initialize nuclide grid
call nuclides(n) % init_grid()
else
err = E_DATA
call set_errmsg("Nuclide '" // trim(name_) // "' is not present &
@ -1100,4 +1103,16 @@ contains
E = nuclides(i_nuclide + 1) % multipole % E_max
end function
subroutine nuclide_multipole_eval(i_nuclide, E, sqrtkT, sig_s, sig_a, sig_f) bind(C)
integer(C_INT), value :: i_nuclide
real(C_DOUBLE), value :: E
real(C_DOUBLE), value :: sqrtkT
real(C_DOUBLE), intent(out) :: sig_s
real(C_DOUBLE), intent(out) :: sig_a
real(C_DOUBLE), intent(out) :: sig_f
call multipole_eval(nuclides(i_nuclide + 1) % multipole, E, sqrtkT, &
sig_s, sig_a, sig_f)
end subroutine
end module nuclide_header