From a6e75c35098f843cbf6973847749085d7842b6f2 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Sun, 13 May 2018 06:44:48 -0400 Subject: [PATCH] Cleaned up the optional parameters interface of rotate_angle --- openmc/capi/math.py | 9 ++++++--- src/math.F90 | 17 ++++++++++------- src/math_functions.cpp | 6 +++--- src/math_functions.h | 2 +- tests/unit_tests/test_math.py | 12 +++++++++++- 5 files changed, 31 insertions(+), 15 deletions(-) diff --git a/openmc/capi/math.py b/openmc/capi/math.py index c46390bced..95e30eff3f 100644 --- a/openmc/capi/math.py +++ b/openmc/capi/math.py @@ -1,4 +1,4 @@ -from ctypes import (c_int, c_double, POINTER) +from ctypes import (c_int, c_double, POINTER, CFUNCTYPE) import numpy as np from numpy.ctypeslib import ndpointer @@ -201,8 +201,11 @@ def rotate_angle(uvw0, mu, phi=None): uvw = np.zeros(3, dtype=np.float64) uvw0_arr = np.array(uvw0, dtype=np.float64) - _dll.rotate_angle(uvw0_arr, c_double(mu), - uvw, c_double(phi)) + if phi is None: + _dll.rotate_angle(uvw0_arr, c_double(mu), uvw, None) + else: + _dll.rotate_angle(uvw0_arr, c_double(mu), uvw, c_double(phi)) + return uvw diff --git a/src/math.F90 b/src/math.F90 index 7be3e5a6dc..057ebb213f 100644 --- a/src/math.F90 +++ b/src/math.F90 @@ -75,7 +75,7 @@ module math implicit none real(C_DOUBLE), intent(inout) :: uvw(3) real(C_DOUBLE), value, intent(in) :: mu - real(C_DOUBLE), value, intent(in) :: phi + real(C_DOUBLE), optional, intent(in) :: phi end subroutine rotate_angle_c_intfc function maxwell_spectrum_c_intfc(T) bind(C, name='maxwell_spectrum_c') & @@ -225,16 +225,19 @@ contains !=============================================================================== subroutine rotate_angle(uvw0, mu, uvw, phi) bind(C) - real(C_DOUBLE), intent(in) :: uvw0(3) ! directional cosine - real(C_DOUBLE), intent(in) :: mu ! cosine of angle in lab or CM - real(C_DOUBLE), intent(out) :: uvw(3) ! rotated directional cosine - real(C_DOUBLE), optional :: phi ! azimuthal angle + real(C_DOUBLE), intent(in) :: uvw0(3) ! directional cosine + real(C_DOUBLE), intent(in) :: mu ! cosine of angle in lab or CM + real(C_DOUBLE), intent(out) :: uvw(3) ! rotated directional cosine + real(C_DOUBLE), target, optional :: phi ! azimuthal angle + + real(C_DOUBLE), pointer :: phi_ptr uvw = uvw0 if (present(phi)) then - call rotate_angle_c_intfc(uvw, mu, phi) + phi_ptr => phi + call rotate_angle_c_intfc(uvw, mu, phi_ptr) else - call rotate_angle_c_intfc(uvw, mu, -10._8) + call rotate_angle_c_intfc(uvw, mu) end if end subroutine rotate_angle diff --git a/src/math_functions.cpp b/src/math_functions.cpp index c9d6fc3590..78ff08b055 100644 --- a/src/math_functions.cpp +++ b/src/math_functions.cpp @@ -659,7 +659,7 @@ void calc_zn_c(const int n, const double rho, const double phi, double zn[]) { // and SERPENT. //============================================================================== -void rotate_angle_c(double uvw[3], const double mu, double phi) { +void rotate_angle_c(double uvw[3], const double mu, double* phi) { double phi_; // azimuthal angle double sinphi; // std::sine of azimuthal angle double cosphi; // cosine of azimuthal angle @@ -675,8 +675,8 @@ void rotate_angle_c(double uvw[3], const double mu, double phi) { w0 = uvw[2]; // Sample azimuthal angle in [0,2pi) if none provided - if (phi != -10.) { - phi_ = phi; + if (phi != NULL) { + phi_ = (*phi); } else { phi_ = 2. * PI * prn(); } diff --git a/src/math_functions.h b/src/math_functions.h index 8df72f56cb..5417281ea4 100644 --- a/src/math_functions.h +++ b/src/math_functions.h @@ -72,7 +72,7 @@ extern "C" void calc_zn_c(const int n, const double rho, const double phi, //============================================================================== extern "C" void rotate_angle_c(double uvw[3], const double mu, - double phi = -10.); + double* phi=NULL); //============================================================================== // MAXWELL_SPECTRUM samples an energy from the Maxwell fission distribution diff --git a/tests/unit_tests/test_math.py b/tests/unit_tests/test_math.py index d48119a372..487d02b107 100644 --- a/tests/unit_tests/test_math.py +++ b/tests/unit_tests/test_math.py @@ -170,7 +170,17 @@ def test_rotate_angle(): assert np.array_equal(ref_uvw, test_uvw) - # Need to test phi=None somehow... + # Now to test phi is None + mu = 0.9 + settings = openmc.capi.settings + settings.seed = 1 + + # When seed = 1, phi will be sampled as 1.9116495709698769 + # The resultant reference is from hand-calculations given the above + ref_uvw = [0.9, 0.410813051297112, 0.1457142302040] + test_uvw = openmc.capi.math.rotate_angle(uvw0, mu) + + assert np.allclose(ref_uvw, test_uvw) def test_maxwell_spectrum():