From ad4a5da040ef69325d990bbd4320fa0310a4f08c Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Fri, 16 Feb 2018 18:25:52 -0500 Subject: [PATCH] Move lattice get_local_xyz to C++ --- src/geometry_header.F90 | 95 +++++++++-------------------------------- src/lattice.cpp | 27 ++++++++++++ src/lattice.h | 11 ++++- 3 files changed, 56 insertions(+), 77 deletions(-) diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 index c58ab1d6c..f5d79d0a6 100644 --- a/src/geometry_header.F90 +++ b/src/geometry_header.F90 @@ -125,6 +125,16 @@ module geometry_header integer(C_INT), intent(out) :: i_xyz(3) end subroutine lattice_get_indices_c + subroutine lattice_get_local_xyz_c(lat_ptr, global_xyz, i_xyz, local_xyz) & + bind(C, name='lattice_get_local_xyz') + use ISO_C_BINDING + implicit none + type(C_PTR), intent(in), value :: lat_ptr + real(C_DOUBLE), intent(in) :: global_xyz(3) + integer(C_INT), intent(in) :: i_xyz(3) + real(C_DOUBLE), intent(out) :: local_xyz(3) + end subroutine lattice_get_local_xyz_c + subroutine lattice_to_hdf5_c(lat_ptr, group) bind(C, name='lattice_to_hdf5') use ISO_C_BINDING use hdf5 @@ -167,25 +177,10 @@ module geometry_header procedure :: are_valid_indices => lattice_are_valid_indices procedure :: distance => lattice_distance procedure :: get_indices => lattice_get_indices + procedure :: get_local_xyz => lattice_get_local_xyz procedure :: to_hdf5 => lattice_to_hdf5 - - procedure(lattice_get_local_xyz_), deferred :: get_local_xyz end type Lattice - abstract interface - -!=============================================================================== -! GET_LOCAL_XYZ returns the translated local version of the given global xyz. - - function lattice_get_local_xyz_(this, global_xyz, i_xyz) result(local_xyz) - import Lattice - class(Lattice), intent(in) :: this - real(8), intent(in) :: global_xyz(3) - integer, intent(in) :: i_xyz(3) - real(8) :: local_xyz(3) - end function lattice_get_local_xyz_ - end interface - !=============================================================================== ! RECTLATTICE extends LATTICE for rectilinear arrays. !=============================================================================== @@ -193,10 +188,6 @@ module geometry_header type, extends(Lattice) :: RectLattice integer :: n_cells(3) ! Number of cells along each axis real(8), allocatable :: lower_left(:) ! Global lower-left corner of lat - - contains - - procedure :: get_local_xyz => get_local_rect end type RectLattice !=============================================================================== @@ -207,10 +198,6 @@ module geometry_header integer :: n_rings ! Number of radial ring cell positoins integer :: n_axial ! Number of axial cell positions real(8), allocatable :: center(:) ! Global center of lattice - - contains - - procedure :: get_local_xyz => get_local_hex end type HexLattice !=============================================================================== @@ -311,63 +298,21 @@ contains call lattice_get_indices_c(this % ptr, xyz, i_xyz) end function lattice_get_indices + function lattice_get_local_xyz(this, global_xyz, i_xyz) & + result(local_xyz) + class(Lattice), intent(in) :: this + real(C_DOUBLE), intent(in) :: global_xyz(3) + integer(C_INT), intent(in) :: i_xyz(3) + real(C_DOUBLE) :: local_xyz(3) + call lattice_get_local_xyz_c(this % ptr, global_xyz, i_xyz, local_xyz) + end function lattice_get_local_xyz + subroutine lattice_to_hdf5(this, group) class(Lattice), intent(in) :: this integer(HID_T), intent(in) :: group call lattice_to_hdf5_c(this % ptr, group) end subroutine lattice_to_hdf5 -!=============================================================================== - - function get_local_rect(this, global_xyz, i_xyz) result(local_xyz) - class(RectLattice), intent(in) :: this - real(8), intent(in) :: global_xyz(3) - integer, intent(in) :: i_xyz(3) - real(8) :: local_xyz(3) - - real(8) :: xyz(3) ! global_xyz alias - - xyz = global_xyz - - local_xyz(1) = xyz(1) - (this % lower_left(1) + & - (i_xyz(1) - HALF)*this % pitch(1)) - local_xyz(2) = xyz(2) - (this % lower_left(2) + & - (i_xyz(2) - HALF)*this % pitch(2)) - if (this % is_3d) then - local_xyz(3) = xyz(3) - (this % lower_left(3) + & - (i_xyz(3) - HALF)*this % pitch(3)) - else - local_xyz(3) = xyz(3) - end if - end function get_local_rect - -!=============================================================================== - - function get_local_hex(this, global_xyz, i_xyz) result(local_xyz) - class(HexLattice), intent(in) :: this - real(8), intent(in) :: global_xyz(3) - integer, intent(in) :: i_xyz(3) - real(8) :: local_xyz(3) - - real(8) :: xyz(3) ! global_xyz alias - - xyz = global_xyz - - ! x_l = x_g - (center + pitch_x*cos(30)*index_x) - local_xyz(1) = xyz(1) - (this % center(1) + & - sqrt(THREE) / TWO * (i_xyz(1) - this % n_rings) * this % pitch(1)) - ! y_l = y_g - (center + pitch_x*index_x + pitch_y*sin(30)*index_y) - local_xyz(2) = xyz(2) - (this % center(2) + & - (i_xyz(2) - this % n_rings) * this % pitch(1) + & - (i_xyz(1) - this % n_rings) * this % pitch(1) / TWO) - if (this % is_3d) then - local_xyz(3) = xyz(3) - this % center(3) & - + (HALF*this % n_axial - i_xyz(3) + HALF) * this % pitch(2) - else - local_xyz(3) = xyz(3) - end if - end function get_local_hex - !=============================================================================== function cell_id(this) result(id) diff --git a/src/lattice.cpp b/src/lattice.cpp index 82f15401d..fd89cbf1e 100644 --- a/src/lattice.cpp +++ b/src/lattice.cpp @@ -219,6 +219,22 @@ RectLattice::get_indices(const double xyz[3]) const return {ix, iy, iz}; } +//============================================================================== + +std::array +RectLattice::get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const +{ + std::array local_xyz; + local_xyz[0] = global_xyz[0] - (lower_left[0] + (i_xyz[0] - 0.5)*pitch[0]); + local_xyz[1] = global_xyz[1] - (lower_left[1] + (i_xyz[1] - 0.5)*pitch[1]); + if (is_3d) { + local_xyz[2] = global_xyz[2] - (lower_left[2] + (i_xyz[2] - 0.5)*pitch[2]); + } else { + local_xyz[2] = global_xyz[2]; + } + return local_xyz; +} + //============================================================================== // HexLattice implementation //============================================================================== @@ -573,6 +589,8 @@ HexLattice::get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const return local_xyz; } +//============================================================================== +// Non-method functions //============================================================================== extern "C" void @@ -629,6 +647,15 @@ extern "C" { i_xyz[2] = inds[2]; } + void lattice_get_local_xyz(Lattice *lat, const double global_xyz[3], + const int i_xyz[3], double local_xyz[3]) + { + std::array xyz {lat->get_local_xyz(global_xyz, i_xyz)}; + local_xyz[0] = xyz[0]; + local_xyz[1] = xyz[1]; + local_xyz[2] = xyz[2]; + } + void lattice_to_hdf5(Lattice *lat, hid_t group) {lat->to_hdf5(group);} } diff --git a/src/lattice.h b/src/lattice.h index 81b0c2bdf..9ade739c0 100644 --- a/src/lattice.h +++ b/src/lattice.h @@ -79,8 +79,12 @@ public: //! @return An array containing the indices of a lattice tile. virtual std::array get_indices(const double xyz[3]) const = 0; - //virtual void get_local_xyz(const double global_xyz[3], const int i_xyz[3], - // double local_xyz[3]) const = 0; + //! Get coordinates local to a lattice tile. + //! @param global_xyz[3] A 3D Cartesian coordinate. + //! @param i_xyz[3] The indices for a lattice tile. + //! @return Local 3D Cartesian coordinates. + virtual std::array + get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const = 0; //! Write all information needed to reconstruct the lattice to an HDF5 group. //! @param group_id An HDF5 group id. @@ -107,6 +111,9 @@ public: std::array get_indices(const double xyz[3]) const; + std::array + get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const; + protected: std::array n_cells; //! Number of cells along each axis std::array lower_left; //! Global lower-left corner of the lattice