From 42e1148fa68bce7442e5767281093ee54fa722e9 Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Thu, 10 May 2018 00:12:22 -0400 Subject: [PATCH] Move most lattice summary.h5 writing to C++ --- src/geometry_header.F90 | 2 - src/hdf5_interface.h | 10 ++++- src/input_xml.F90 | 22 ----------- src/lattice.cpp | 87 +++++++++++++++++++++++++++++++++++++++-- src/lattice.h | 1 - src/summary.F90 | 65 +----------------------------- 6 files changed, 93 insertions(+), 94 deletions(-) diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 index 025fd69c5..02f9b7e82 100644 --- a/src/geometry_header.F90 +++ b/src/geometry_header.F90 @@ -198,7 +198,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 end type RectLattice !=============================================================================== @@ -208,7 +207,6 @@ module geometry_header type, extends(Lattice) :: HexLattice 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 end type HexLattice !=============================================================================== diff --git a/src/hdf5_interface.h b/src/hdf5_interface.h index 66757df8e..232d2f6e7 100644 --- a/src/hdf5_interface.h +++ b/src/hdf5_interface.h @@ -98,8 +98,16 @@ write_int(hid_t group_id, char const *name, int32_t buffer) H5Dclose(dataset); } +template void +write_int(hid_t group_id, char const *name, + const std::array &buffer, bool indep) +{ + hsize_t dims[1] {array_len}; + write_dataset(group_id, 1, dims, name, H5T_NATIVE_INT, buffer.data(), indep); +} -template inline void + +template void write_double_1D(hid_t group_id, char const *name, const std::array &buffer) { diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 717e2cada..4d48e8d14 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -1284,15 +1284,6 @@ contains call fatal_error("Rectangular lattice must be two or three dimensions.") end if - ! Read lattice lower-left location - if (node_word_count(node_lat, "lower_left") /= n) then - call fatal_error("Number of entries on must be the same & - &as the number of entries on .") - end if - - allocate(lat % lower_left(n)) - call get_node_array(node_lat, "lower_left", lat % lower_left) - ! TODO: Remove deprecation warning in a future release. if (check_for_node(node_lat, "type")) then call warning("The use of 'type' is no longer needed. The utility & @@ -1365,19 +1356,6 @@ contains lat % is_3d = .false. end if - ! Read lattice lower-left location - n = node_word_count(node_lat, "center") - if (lat % is_3d .and. n /= 3) then - call fatal_error("A hexagonal lattice with must have & - &
specified by 3 numbers.") - else if ((.not. lat % is_3d) .and. n /= 2) then - call fatal_error("A hexagonal lattice without must have & - &
specified by 2 numbers.") - end if - - allocate(lat % center(n)) - call get_node_array(node_lat, "center", lat % center) - ! Copy number of dimensions n_rings = lat % n_rings n_z = lat % n_axial diff --git a/src/lattice.cpp b/src/lattice.cpp index a28d7efe7..5723b2280 100644 --- a/src/lattice.cpp +++ b/src/lattice.cpp @@ -272,11 +272,57 @@ RectLattice::get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const void RectLattice::to_hdf5_inner(hid_t lat_group) const { + // Write basic lattice information. + write_string(lat_group, "type", "rectangular", false); if (is_3d) { write_double_1D(lat_group, "pitch", pitch); + write_double_1D(lat_group, "lower_left", lower_left); + write_int(lat_group, "dimension", n_cells, false); } else { - std::array out {{pitch[0], pitch[1]}}; - write_double_1D(lat_group, "pitch", out); + std::array pitch_short {{pitch[0], pitch[1]}}; + write_double_1D(lat_group, "pitch", pitch_short); + std::array ll_short {{lower_left[0], lower_left[1]}}; + write_double_1D(lat_group, "lower_left", ll_short); + std::array nc_short {{n_cells[0], n_cells[1]}}; + write_int(lat_group, "dimension", nc_short, false); + } + + // Write the universe ids. The convention here is to switch the ordering on + // the y-axis to match the way universes are input in a text file. + if (is_3d) { + hsize_t nx {static_cast(n_cells[0])}; + hsize_t ny {static_cast(n_cells[1])}; + hsize_t nz {static_cast(n_cells[2])}; + int out[nx*ny*nz]; + + for (int m = 0; m < nz; m++) { + for (int k = 0; k < ny; k++) { + for (int j = 0; j < nx; j++) { + int indx1 = nx*ny*m + nx*k + j; + int indx2 = nx*ny*m + nx*(ny-k-1) + j; + out[indx2] = universes_c[universes[indx1]]->id; + } + } + } + + hsize_t dims[3] {nz, ny, nx}; + write_int(lat_group, 3, dims, "universes", out, false); + + } else { + hsize_t nx {static_cast(n_cells[0])}; + hsize_t ny {static_cast(n_cells[1])}; + int out[nx*ny]; + + for (int k = 0; k < ny; k++) { + for (int j = 0; j < nx; j++) { + int indx1 = nx*k + j; + int indx2 = nx*(ny-k-1) + j; + out[indx2] = universes_c[universes[indx1]]->id; + } + } + + hsize_t dims[3] {1, ny, nx}; + write_int(lat_group, 3, dims, "universes", out, false); } } @@ -670,12 +716,45 @@ HexLattice::get_local_xyz(const double global_xyz[3], const int i_xyz[3]) const void HexLattice::to_hdf5_inner(hid_t lat_group) const { + // Write basic lattice information. + write_string(lat_group, "type", "hexagonal", false); + write_int(lat_group, 0, nullptr, "n_rings", &n_rings, false); + write_int(lat_group, 0, nullptr, "n_axial", &n_axial, false); if (is_3d) { write_double_1D(lat_group, "pitch", pitch); + write_double_1D(lat_group, "center", center); } else { - std::array out {{pitch[0]}}; - write_double_1D(lat_group, "pitch", out); + std::array pitch_short {{pitch[0]}}; + write_double_1D(lat_group, "pitch", pitch_short); + std::array center_short {{center[0], center[1]}}; + write_double_1D(lat_group, "center", center_short); } + + // Write the universe ids. + hsize_t nx {static_cast(2*n_rings - 1)}; + hsize_t ny {static_cast(2*n_rings - 1)}; + hsize_t nz {static_cast(n_axial)}; + int out[nx*ny*nz]; + + for (int m = 0; m < nz; m++) { + for (int k = 0; k < ny; k++) { + for (int j = 0; j < nx; j++) { + int indx = nx*ny*m + nx*k + j; + if (j + k < n_rings - 1) { + // This array position is never used; put a -1 to indicate this. + out[indx] = -1; + } else if (j + k > 3*n_rings - 3) { + // This array position is never used; put a -1 to indicate this. + out[indx] = -1; + } else { + out[indx] = universes_c[universes[indx]]->id; + } + } + } + } + + hsize_t dims[3] {nz, ny, nx}; + write_int(lat_group, 3, dims, "universes", out, false); } //============================================================================== diff --git a/src/lattice.h b/src/lattice.h index 78d83a2cb..8bdf469c3 100644 --- a/src/lattice.h +++ b/src/lattice.h @@ -41,7 +41,6 @@ class Lattice public: int32_t id; //! Universe ID number std::string name; //! User-defined name - //std::vector pitch; //! Pitch along each basis std::vector universes; //! Universes filling each lattice tile int32_t outer{NO_OUTER_UNIVERSE}; //! Universe tiled outside the lattice //std::vector offset; //! Distribcell offsets diff --git a/src/summary.F90 b/src/summary.F90 index 4ff27963e..22de03e7b 100644 --- a/src/summary.F90 +++ b/src/summary.F90 @@ -160,8 +160,7 @@ contains subroutine write_geometry(file_id) integer(HID_T), intent(in) :: file_id - integer :: i, j, k, m - integer, allocatable :: lattice_universes(:,:,:) + integer :: i, j integer, allocatable :: cell_materials(:) integer, allocatable :: cell_ids(:) real(8), allocatable :: cell_temperatures(:) @@ -311,68 +310,6 @@ contains call write_dataset(lattice_group, "outer", lat % outer) end if - select type (lat) - type is (RectLattice) - ! Write lattice type. - call write_dataset(lattice_group, "type", "rectangular") - - ! Write lattice dimensions, lower left corner, and pitch - if (lat % is_3d) then - call write_dataset(lattice_group, "dimension", lat % n_cells) - call write_dataset(lattice_group, "lower_left", lat % lower_left) - else - call write_dataset(lattice_group, "dimension", lat % n_cells(1:2)) - call write_dataset(lattice_group, "lower_left", lat % lower_left) - end if - - ! Write lattice universes. - allocate(lattice_universes(lat%n_cells(1), lat%n_cells(2), & - &lat%n_cells(3))) - do j = 1, lat%n_cells(1) - do k = 0, lat%n_cells(2) - 1 - do m = 1, lat%n_cells(3) - lattice_universes(j, k+1, m) = & - universes(lat%get([j-1, lat%n_cells(2)-k-1, m-1])+1)%id - end do - end do - end do - - type is (HexLattice) - ! Write lattice type. - call write_dataset(lattice_group, "type", "hexagonal") - - ! Write number of lattice cells. - call write_dataset(lattice_group, "n_rings", lat%n_rings) - call write_dataset(lattice_group, "n_axial", lat%n_axial) - - ! Write lattice center - call write_dataset(lattice_group, "center", lat%center) - - ! Write lattice universes. - allocate(lattice_universes(2*lat%n_rings - 1, 2*lat%n_rings - 1, & - &lat%n_axial)) - do m = 1, lat%n_axial - do k = 1, 2*lat%n_rings - 1 - do j = 1, 2*lat%n_rings - 1 - if (j + k < lat%n_rings + 1) then - ! This array position is never used; put a -1 to indicate this - lattice_universes(j,k,m) = -1 - cycle - else if (j + k > 3*lat%n_rings - 1) then - ! This array position is never used; put a -1 to indicate this - lattice_universes(j,k,m) = -1 - cycle - end if - lattice_universes(j,k,m) = universes(lat%get([j-1,k-1,m-1])+1)%id - end do - end do - end do - end select - - ! Write lattice universes - call write_dataset(lattice_group, "universes", lattice_universes) - deallocate(lattice_universes) - call close_group(lattice_group) end do LATTICE_LOOP