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

Cleaned up the C and Fortran interfaces; includes removing interfaces which are no longer needed and removing the _c from C functions which no longer have Fortran links

Browse source
This commit is contained in:
Adam G Nelson 2018-10-12 11:51:23 -04:00
parent e3e388bb6d
commit a9641b8e2a
14 changed files with 54 additions and 156 deletions
Download patch file Download diff file Expand all files Collapse all files

8
include/openmc/math_functions.h
View file

@ -22,7 +22,7 @@ namespace openmc {
//! @return The requested percentile
//==============================================================================
extern "C" double normal_percentile_c(double p);
extern "C" double normal_percentile(double p);
//==============================================================================
//! Calculate the percentile of the Student's t distribution with a specified
@ -58,7 +58,7 @@ extern "C" void calc_pn_c(int n, double x, double pnx[]);
//! evaluated at x
//==============================================================================
extern "C" double evaluate_legendre_c(int n, const double data[], double x);
extern "C" double evaluate_legendre(int n, const double data[], double x);
//==============================================================================
//! Calculate the n-th order real spherical harmonics for a given angle (in
@ -144,7 +144,7 @@ Direction rotate_angle(Direction u, double mu, double* phi);
//! @result The sampled outgoing energy
//==============================================================================
extern "C" double maxwell_spectrum_c(double T);
extern "C" double maxwell_spectrum(double T);
//==============================================================================
//! Samples an energy from a Watt energy-dependent fission distribution.
@ -159,7 +159,7 @@ extern "C" double maxwell_spectrum_c(double T);
//! @result The sampled outgoing energy
//==============================================================================
extern "C" double watt_spectrum_c(double a, double b);
extern "C" double watt_spectrum(double a, double b);
//==============================================================================
//! Doppler broadens the windowed multipole curvefit.

10
include/openmc/mgxs_interface.h
View file

@ -18,6 +18,9 @@ extern std::vector<Mgxs> nuclides_MG;
extern std::vector<Mgxs> macro_xs;
extern "C" int num_energy_groups;
//TODO: When more of the Fortran is converted (input_xml, tallies, etc, also
// bring over energy_bin_avg, energy_bins, etc, as vectors)
//==============================================================================
// Mgxs data loading interface methods
//==============================================================================
@ -44,13 +47,6 @@ extern "C" void
calculate_xs_c(int i_mat, int gin, double sqrtkT, const double uvw[3],
double& total_xs, double& abs_xs, double& nu_fiss_xs);
extern "C" void
sample_scatter_c(int i_mat, int gin, int& gout, double& mu, double& wgt,
double uvw[3]);
extern "C" void
sample_fission_energy_c(int i_mat, int gin, int& dg, int& gout);
extern "C" double
get_nuclide_xs_c(int index, int xstype, int gin, int* gout, double* mu, int* dg);

13
include/openmc/physics_mg.h
View file

@ -10,22 +10,21 @@
namespace openmc {
//==============================================================================
// SCATTER
//==============================================================================
//TODO: Remove energy_bin_avg and material_xs parameters when they reside on
// the C-side this should happen after materials, physics, input, and tallies
// are brought over
//! \brief samples particle behavior after a collision event.
extern "C" void
collision_mg(Particle* p, Bank* fission_bank, const int64_t fission_bank_size,
const double* energy_bin_avg, const MaterialMacroXS& material_xs);
collision_mg(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS& material_xs);
//! \brief samples a reaction type.
//!
//! Note that there is special logic when suvival biasing is turned on since
//! fission and disappearance are treated implicitly.
void
sample_reaction(Particle* p, Bank* fission_bank,
const int64_t fission_bank_size, const double* energy_bin_avg,
sample_reaction(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS& material_xs);
//! \brief Samples the scattering event

23
openmc/capi/math.py
View file

@ -12,8 +12,8 @@ _dll.t_percentile_c.argtypes = [c_double, c_int]
_dll.calc_pn_c.restype = None
_dll.calc_pn_c.argtypes = [c_int, c_double, ndpointer(c_double)]
_dll.evaluate_legendre_c.restype = c_double
_dll.evaluate_legendre_c.argtypes = [c_int, POINTER(c_double), c_double]
_dll.evaluate_legendre.restype = c_double
_dll.evaluate_legendre.argtypes = [c_int, POINTER(c_double), c_double]
_dll.calc_rn_c.restype = None
_dll.calc_rn_c.argtypes = [c_int, ndpointer(c_double), ndpointer(c_double)]
@ -27,11 +27,11 @@ _dll.calc_zn_rad_c.argtypes = [c_int, c_double, ndpointer(c_double)]
_dll.rotate_angle_c.restype = None
_dll.rotate_angle_c.argtypes = [ndpointer(c_double), c_double,
POINTER(c_double)]
_dll.maxwell_spectrum_c.restype = c_double
_dll.maxwell_spectrum_c.argtypes = [c_double]
_dll.maxwell_spectrum.restype = c_double
_dll.maxwell_spectrum.argtypes = [c_double]
_dll.watt_spectrum_c.restype = c_double
_dll.watt_spectrum_c.argtypes = [c_double, c_double]
_dll.watt_spectrum.restype = c_double
_dll.watt_spectrum.argtypes = [c_double, c_double]
_dll.broaden_wmp_polynomials_c.restype = None
_dll.broaden_wmp_polynomials_c.argtypes = [c_double, c_double, c_int,
@ -100,9 +100,8 @@ def evaluate_legendre(data, x):
"""
data_arr = np.array(data, dtype=np.float64)
return _dll.evaluate_legendre_c(len(data),
data_arr.ctypes.data_as(POINTER(c_double)),
x)
return _dll.evaluate_legendre(len(data),
data_arr.ctypes.data_as(POINTER(c_double)), x)
def calc_rn(n, uvw):
@ -182,7 +181,7 @@ def calc_zn_rad(n, rho):
zn_rad = np.zeros(num_bins, dtype=np.float64)
_dll.calc_zn_rad_c(n, rho, zn_rad)
return zn_rad
def rotate_angle(uvw0, mu, phi=None):
""" Rotates direction cosines through a polar angle whose cosine is
@ -231,7 +230,7 @@ def maxwell_spectrum(T):
"""
return _dll.maxwell_spectrum_c(T)
return _dll.maxwell_spectrum(T)
def watt_spectrum(a, b):
@ -251,7 +250,7 @@ def watt_spectrum(a, b):
"""
return _dll.watt_spectrum_c(a, b)
return _dll.watt_spectrum(a, b)
def broaden_wmp_polynomials(E, dopp, n):

4
src/distribution.cpp
View file

@ -90,7 +90,7 @@ Maxwell::Maxwell(pugi::xml_node node)
double Maxwell::sample() const
{
return maxwell_spectrum_c(theta_);
return maxwell_spectrum(theta_);
}
//==============================================================================
@ -110,7 +110,7 @@ Watt::Watt(pugi::xml_node node)
double Watt::sample() const
{
return watt_spectrum_c(a_, b_);
return watt_spectrum(a_, b_);
}
//==============================================================================

4
src/distribution_energy.cpp
View file

@ -279,7 +279,7 @@ double MaxwellEnergy::sample(double E) const
while (true) {
// Sample maxwell fission spectrum
double E_out = maxwell_spectrum_c(theta);
double E_out = maxwell_spectrum(theta);
// Accept energy based on restriction energy
if (E_out <= E - u_) return E_out;
@ -343,7 +343,7 @@ double WattEnergy::sample(double E) const
while (true) {
// Sample energy-dependent Watt fission spectrum
double E_out = watt_spectrum_c(a, b);
double E_out = watt_spectrum(a, b);
// Accept energy based on restriction energy
if (E_out <= E - u_) return E_out;

44
src/math.F90
View file

@ -12,10 +12,7 @@ module math
public :: calc_rn
public :: calc_zn
public :: calc_zn_rad
public :: evaluate_legendre
public :: rotate_angle
public :: maxwell_spectrum
public :: watt_spectrum
public :: faddeeva
public :: w_derivative
public :: broaden_wmp_polynomials
@ -42,16 +39,6 @@ module math
real(C_DOUBLE), intent(out) :: pnx(n + 1)
end subroutine calc_pn
pure function evaluate_legendre_c_intfc(n, data, x) &
bind(C, name='evaluate_legendre_c') result(val)
use ISO_C_BINDING
implicit none
integer(C_INT), value, intent(in) :: n
real(C_DOUBLE), intent(in) :: data(n)
real(C_DOUBLE), value, intent(in) :: x
real(C_DOUBLE) :: val
end function evaluate_legendre_c_intfc
pure subroutine calc_rn(n, uvw, rn) bind(C, name='calc_rn_c')
use ISO_C_BINDING
implicit none
@ -85,23 +72,6 @@ module math
real(C_DOUBLE), optional, intent(in) :: phi
end subroutine rotate_angle_c_intfc
function maxwell_spectrum(T) bind(C, name='maxwell_spectrum_c') &
result(E_out)
use ISO_C_BINDING
implicit none
real(C_DOUBLE), value, intent(in) :: T
real(C_DOUBLE) :: E_out
end function maxwell_spectrum
function watt_spectrum(a, b) bind(C, name='watt_spectrum_c') &
result(E_out)
use ISO_C_BINDING
implicit none
real(C_DOUBLE), value, intent(in) :: a
real(C_DOUBLE), value, intent(in) :: b
real(C_DOUBLE) :: E_out
end function watt_spectrum
subroutine broaden_wmp_polynomials(E, dopp, n, factors) &
bind(C, name='broaden_wmp_polynomials_c')
use ISO_C_BINDING
@ -157,20 +127,6 @@ module math
contains
!===============================================================================
! EVALUATE_LEGENDRE Find the value of f(x) given a set of Legendre coefficients
! and the value of x
!===============================================================================
pure function evaluate_legendre(data, x) result(val) bind(C)
real(C_DOUBLE), intent(in) :: data(:)
real(C_DOUBLE), intent(in) :: x
real(C_DOUBLE) :: val
val = evaluate_legendre_c_intfc(size(data) - 1, data, x)
end function evaluate_legendre
!===============================================================================
! ROTATE_ANGLE rotates direction cosines through a polar angle whose cosine is
! mu and through an azimuthal angle sampled uniformly. Note that this is done

12
src/math_functions.cpp
View file

@ -6,7 +6,7 @@ namespace openmc {
// Mathematical methods
//==============================================================================
double normal_percentile_c(double p) {
double normal_percentile(double p) {
constexpr double p_low = 0.02425;
constexpr double a[6] = {-3.969683028665376e1, 2.209460984245205e2,
-2.759285104469687e2, 1.383577518672690e2,
@ -79,7 +79,7 @@ double t_percentile_c(double p, int df){
// 16 (4), pp. 1123-1132 (1987).
double n = df;
double k = 1. / (n - 2.);
double z = normal_percentile_c(p);
double z = normal_percentile(p);
double z2 = z * z;
t = std::sqrt(n * k) * (z + (z2 - 3.) * z * k / 4. + ((5. * z2 - 56.) * z2 +
75.) * z * k * k / 96. + (((z2 - 27.) * 3. * z2 + 417.) * z2 - 315.) *
@ -103,7 +103,7 @@ void calc_pn_c(int n, double x, double pnx[]) {
}
double evaluate_legendre_c(int n, const double data[], double x) {
double evaluate_legendre(int n, const double data[], double x) {
double pnx[n + 1];
double val = 0.0;
calc_pn_c(n, x, pnx);
@ -658,7 +658,7 @@ Direction rotate_angle(Direction u, double mu, double* phi)
}
double maxwell_spectrum_c(double T) {
double maxwell_spectrum(double T) {
// Set the random numbers
double r1 = prn();
double r2 = prn();
@ -674,8 +674,8 @@ double maxwell_spectrum_c(double T) {
}
double watt_spectrum_c(double a, double b) {
double w = maxwell_spectrum_c(a);
double watt_spectrum(double a, double b) {
double w = maxwell_spectrum(a);
double E_out = w + 0.25 * a * a * b + (2. * prn() - 1.) * std::sqrt(a * a * b * w);
return E_out;

20
src/mgxs_interface.F90
View file

@ -62,26 +62,6 @@ module mgxs_interface
real(C_DOUBLE), intent(inout) :: nu_fiss_xs
end subroutine calculate_xs_c
subroutine sample_scatter_c(i_mat, gin, gout, mu, wgt, uvw) bind(C)
use ISO_C_BINDING
implicit none
integer(C_INT), value, intent(in) :: i_mat
integer(C_INT), value, intent(in) :: gin
integer(C_INT), intent(inout) :: gout
real(C_DOUBLE), intent(inout) :: mu
real(C_DOUBLE), intent(inout) :: wgt
real(C_DOUBLE), intent(inout) :: uvw(1:3)
end subroutine sample_scatter_c
subroutine sample_fission_energy_c(i_mat, gin, dg, gout) bind(C)
use ISO_C_BINDING
implicit none
integer(C_INT), value, intent(in) :: i_mat
integer(C_INT), value, intent(in) :: gin
integer(C_INT), intent(inout) :: dg
integer(C_INT), intent(inout) :: gout
end subroutine sample_fission_energy_c
subroutine get_name_c(index, name_len, name) bind(C)
use ISO_C_BINDING
implicit none

30
src/mgxs_interface.cpp
View file

@ -97,36 +97,6 @@ calculate_xs_c(int i_mat, int gin, double sqrtkT, const double uvw[3],
//==============================================================================
void
sample_scatter_c(int i_mat, int gin, int& gout, double& mu, double& wgt,
double uvw[3])
{
int gout_c = gout - 1;
macro_xs[i_mat - 1].sample_scatter(gin - 1, gout_c, mu, wgt);
// adjust return value for fortran indexing
gout = gout_c + 1;
// Rotate the angle
rotate_angle_c(uvw, mu, nullptr);
}
//==============================================================================
void
sample_fission_energy_c(int i_mat, int gin, int& dg, int& gout)
{
int dg_c = 0;
int gout_c = 0;
macro_xs[i_mat - 1].sample_fission_energy(gin - 1, dg_c, gout_c);
// adjust return values for fortran indexing
dg = dg_c + 1;
gout = gout_c + 1;
}
//==============================================================================
double
get_nuclide_xs_c(int index, int xstype, int gin, int* gout, double* mu, int* dg)
{

16
src/physics_mg.cpp
View file

@ -20,15 +20,14 @@
namespace openmc {
void
collision_mg(Particle* p, Bank* fission_bank, const int64_t fission_bank_size,
const double* energy_bin_avg, const MaterialMacroXS& material_xs)
collision_mg(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS& material_xs)
{
// Add to the collision counter for the particle
p->n_collision++;
// Sample the reaction type
sample_reaction(p, fission_bank, fission_bank_size, energy_bin_avg,
material_xs);
sample_reaction(p, energy_bin_avg, material_xs);
// Display information about collision
if ((settings::verbosity >= 10) || (openmc_trace)) {
@ -39,8 +38,7 @@ collision_mg(Particle* p, Bank* fission_bank, const int64_t fission_bank_size,
}
void
sample_reaction(Particle* p, Bank* fission_bank,
const int64_t fission_bank_size, const double* energy_bin_avg,
sample_reaction(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS& material_xs)
{
// Create fission bank sites. Note that while a fission reaction is sampled,
@ -50,7 +48,11 @@ sample_reaction(Particle* p, Bank* fission_bank,
if (materials[p->material - 1]->fissionable) {
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
create_fission_sites(p, fission_bank, n_bank, fission_bank_size,
Bank* result_bank;
int64_t result_bank_size;
// Get pointer to fission bank from Fortran side
openmc_fission_bank(&result_bank, &result_bank_size);
create_fission_sites(p, result_bank, n_bank, result_bank_size,
material_xs);
} else if ((settings::run_mode == RUN_MODE_FIXEDSOURCE) &&
(settings::create_fission_neutrons)) {

12
src/scattdata.cpp
View file

@ -316,8 +316,8 @@ ScattDataLegendre::update_max_val()
}
// Calculate probability
double f = evaluate_legendre_c(dist[gin][i_gout].size() - 1,
dist[gin][i_gout].data(), mu);
double f = evaluate_legendre(dist[gin][i_gout].size() - 1,
dist[gin][i_gout].data(), mu);
// if this is a new maximum, store it
if (f > max_val[gin][i_gout]) max_val[gin][i_gout] = f;
@ -339,8 +339,8 @@ ScattDataLegendre::calc_f(int gin, int gout, double mu)
f = 0.;
} else {
int i_gout = gout - gmin[gin];
f = evaluate_legendre_c(dist[gin][i_gout].size() - 1,
dist[gin][i_gout].data(), mu);
f = evaluate_legendre(dist[gin][i_gout].size() - 1,
dist[gin][i_gout].data(), mu);
}
return f;
}
@ -881,8 +881,8 @@ convert_legendre_to_tabular(ScattDataLegendre& leg, ScattDataTabular& tab,
tab.fmu[gin][i_gout].resize(n_mu);
for (int imu = 0; imu < n_mu; imu++) {
tab.fmu[gin][i_gout][imu] =
evaluate_legendre_c(leg.dist[gin][i_gout].size() - 1,
leg.dist[gin][i_gout].data(), tab.mu[imu]);
evaluate_legendre(leg.dist[gin][i_gout].size() - 1,
leg.dist[gin][i_gout].data(), tab.mu[imu]);
}
// Ensure positivity

4
src/secondary_nbody.cpp
View file

@ -32,13 +32,13 @@ void NBodyPhaseSpace::sample(double E_in, double& E_out, double& mu) const
double E_max = (Ap - 1.0)/Ap * (A_/(A_ + 1.0)*E_in + Q_);
// x is essentially a Maxwellian distribution
double x = maxwell_spectrum_c(1.0);
double x = maxwell_spectrum(1.0);
double y;
double r1, r2, r3, r4, r5, r6;
switch (n_bodies_) {
case 3:
y = maxwell_spectrum_c(1.0);
y = maxwell_spectrum(1.0);
break;
case 4:
r1 = prn();

10
src/tracking.F90
View file

@ -34,12 +34,9 @@ module tracking
implicit none
interface
subroutine collision_mg(p, fission_bank, fission_bank_size, &
energy_bin_avg, material_xs) bind(C)
import Particle, Bank, C_INT64_T, C_DOUBLE, MaterialMacroXS
subroutine collision_mg(p, energy_bin_avg, material_xs) bind(C)
import Particle, C_DOUBLE, MaterialMacroXS
type(Particle), intent(inout) :: p
type(Bank), intent(inout) :: fission_bank(*)
integer(C_INT64_T), value, intent(in) :: fission_bank_size
real(C_DOUBLE), intent(in) :: energy_bin_avg(*)
type(MaterialMacroXS), intent(in) :: material_xs
end subroutine collision_mg
@ -239,8 +236,7 @@ contains
if (run_CE) then
call collision(p)
else
call collision_mg(p, fission_bank, &
size(fission_bank, kind=C_INT64_T), energy_bin_avg, material_xs)
call collision_mg(p, energy_bin_avg, material_xs)
end if
! Score collision estimator tallies -- this is done after a collision