From f92bdee1b77a5c70731442057d013726c77ea0a6 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Mon, 23 Jan 2012 13:31:48 -0500 Subject: [PATCH] Expanded HDF5 output to include geometry and materials. --- src/DEPENDENCIES | 3 + src/hdf5_interface.F90 | 404 ++++++++++++++++++++++++++++++++++------- 2 files changed, 340 insertions(+), 67 deletions(-) diff --git a/src/DEPENDENCIES b/src/DEPENDENCIES index d13f68d9dd..ceb2ea4728 100644 --- a/src/DEPENDENCIES +++ b/src/DEPENDENCIES @@ -74,7 +74,10 @@ global.o: tally_header.o global.o: timing.o hdf5_interface.o: constants.o +hdf5_interface.o: geometry_header.o hdf5_interface.o: global.o +hdf5_interface.o: material_header.o +hdf5_interface.o: string.o initialize.o: ace.o initialize.o: bank_header.o diff --git a/src/hdf5_interface.F90 b/src/hdf5_interface.F90 index 9f5d8cec80..007f99bda3 100644 --- a/src/hdf5_interface.F90 +++ b/src/hdf5_interface.F90 @@ -1,7 +1,10 @@ module hdf5_interface use constants + use geometry_header, only: Cell, Surface, Universe, Lattice use global + use material_header, only: Material + use string, only: to_str #ifdef HDF5 use hdf5 @@ -36,30 +39,24 @@ contains subroutine hdf5_write_header() - integer :: rank = 1 - integer(HSIZE_T) :: dims(1) = (/1/) - character(8) :: date_ - character(10) :: time_ - character(19) :: current_time + character(8) :: date_ + character(10) :: time_ + character(19) :: current_time ! Write version information - call h5ltmake_dataset_int_f(hdf5_output_file, "version_major", & - rank, dims, (/ VERSION_MAJOR /), hdf5_err) - call h5ltmake_dataset_int_f(hdf5_output_file, "version_minor", & - rank, dims, (/ VERSION_MINOR /), hdf5_err) - call h5ltmake_dataset_int_f(hdf5_output_file, "version_release", & - rank, dims, (/ VERSION_RELEASE /), hdf5_err) + call hdf5_make_integer(hdf5_output_file, "version_major", VERSION_MAJOR) + call hdf5_make_integer(hdf5_output_file, "version_minor", VERSION_MINOR) + call hdf5_make_integer(hdf5_output_file, "version_release", VERSION_RELEASE) ! Write current date and time call date_and_time(DATE=date_, TIME=time_) current_time = date_(1:4) // "-" // date_(5:6) // "-" // date_(7:8) // & " " // time_(1:2) // ":" // time_(3:4) // ":" // time_(5:6) - call h5ltmake_dataset_string_f(hdf5_output_file, "/date_and_time", & + call h5ltmake_dataset_string_f(hdf5_output_file, "date_and_time", & current_time, hdf5_err) ! Write MPI information - call h5ltmake_dataset_int_f(hdf5_output_file, "n_procs", & - rank, dims, (/ n_procs /), hdf5_err) + call hdf5_make_integer(hdf5_output_file, "n_procs", n_procs) call h5ltset_attribute_string_f(hdf5_output_file, "n_procs", & "description", "Number of MPI processes", hdf5_err) @@ -71,23 +68,16 @@ contains subroutine hdf5_write_summary() - integer :: rank = 1 - integer(HSIZE_T) :: dims(1) = (/1/) - ! Write criticality information if (problem_type == PROB_CRITICALITY) then ! Need to write integer(8)'s using double instead since there is no H5LT ! call for making a dataset of type long - call h5ltmake_dataset_double_f(hdf5_output_file, "n_particles", & - rank, dims, (/ real(n_particles,8) /), hdf5_err) + call hdf5_make_double(hdf5_output_file, "n_particles", real(n_particles,8)) ! Use H5LT interface to write n_cycles, n_inactive, and n_active - call h5ltmake_dataset_int_f(hdf5_output_file, "n_cycles", & - rank, dims, (/ n_cycles /), hdf5_err) - call h5ltmake_dataset_int_f(hdf5_output_file, "n_inactive", & - rank, dims, (/ n_inactive /), hdf5_err) - call h5ltmake_dataset_int_f(hdf5_output_file, "n_active", & - rank, dims, (/ n_cycles - n_inactive /), hdf5_err) + call hdf5_make_integer(hdf5_output_file, "n_cycles", n_cycles) + call hdf5_make_integer(hdf5_output_file, "n_inactive", n_inactive) + call hdf5_make_integer(hdf5_output_file, "n_active", n_cycles - n_inactive) ! Add description of each variable call h5ltset_attribute_string_f(hdf5_output_file, "n_particles", & @@ -111,22 +101,232 @@ contains subroutine hdf5_write_geometry() - integer :: rank = 1 - integer(HSIZE_T) :: dims(1) = (/1/) + integer :: i, j, k + integer(HSIZE_T) :: dims(1) + integer(HSIZE_T) :: dims2(2) integer(HID_T) :: geometry_group + integer(HID_T) :: cell_group + integer(HID_T) :: surface_group + integer(HID_T) :: universe_group + integer(HID_T) :: lattice_group + integer(HID_T) :: temp_group + integer, allocatable :: lattice_universes(:,:) + type(Cell), pointer :: c => null() + type(Surface), pointer :: s => null() + type(Universe), pointer :: u => null() + type(Lattice), pointer :: l => null() ! Create group for geometry call h5gcreate_f(hdf5_output_file, "/geometry", geometry_group, hdf5_err) ! Use H5LT interface to write number of geometry objects - call h5ltmake_dataset_int_f(geometry_group, "n_cells", & - rank, dims, (/ n_cells /), hdf5_err) - call h5ltmake_dataset_int_f(geometry_group, "n_surfaces", & - rank, dims, (/ n_surfaces /), hdf5_err) - call h5ltmake_dataset_int_f(geometry_group, "n_universes", & - rank, dims, (/ n_universes /), hdf5_err) - call h5ltmake_dataset_int_f(geometry_group, "n_lattices", & - rank, dims, (/ n_lattices /), hdf5_err) + call hdf5_make_integer(geometry_group, "n_cells", n_cells) + call hdf5_make_integer(geometry_group, "n_surfaces", n_surfaces) + call hdf5_make_integer(geometry_group, "n_universes", n_universes) + call hdf5_make_integer(geometry_group, "n_lattices", n_lattices) + + ! ========================================================================== + ! WRITE INFORMATION ON CELLS + + call h5gcreate_f(geometry_group, "cells", cell_group, hdf5_err) + + ! Write information on each cell + do i = 1, n_cells + c => cells(i) + + ! Create group for i-th cell + call h5gcreate_f(cell_group, "cell " // trim(to_str(c % id)), & + temp_group, hdf5_err) + + ! Write universe for this cell + call hdf5_make_integer(temp_group, "universe", & + universes(c % universe) % id) + + ! Write information on what fills this cell + select case (c % type) + case (CELL_NORMAL) + call h5ltmake_dataset_string_f(temp_group, "fill_type", "normal", & + hdf5_err) + call hdf5_make_integer(temp_group, "material", & + materials(c % material) % id) + case (CELL_FILL) + call h5ltmake_dataset_string_f(temp_group, "fill_type", "universe", & + hdf5_err) + call hdf5_make_integer(temp_group, "material", & + universes(c % fill) % id) + case (CELL_LATTICE) + call h5ltmake_dataset_string_f(temp_group, "fill_type", "lattice", & + hdf5_err) + call hdf5_make_integer(temp_group, "lattice", & + lattices(c % fill) % id) + end select + + ! Write list of bounding surfaces + if (c % n_surfaces > 0) then + dims(1) = c % n_surfaces + call h5ltmake_dataset_int_f(temp_group, "surfaces", 1, & + dims, c % surfaces, hdf5_err) + end if + + ! Close group for i-th cell + call h5gclose_f(temp_group, hdf5_err) + end do + + call h5gclose_f(cell_group, hdf5_err) + + ! ========================================================================== + ! WRITE INFORMATION ON SURFACES + + call h5gcreate_f(geometry_group, "surfaces", surface_group, hdf5_err) + + ! Write information on each surface + do i = 1, n_surfaces + s => surfaces(i) + + ! Create group for i-th surface + call h5gcreate_f(surface_group, "surface " // trim(to_str(s % id)), & + temp_group, hdf5_err) + + ! Write surface type + select case (s % type) + case (SURF_PX) + call h5ltmake_dataset_string_f(temp_group, "type", "X Plane", hdf5_err) + case (SURF_PY) + call h5ltmake_dataset_string_f(temp_group, "type", "Y Plane", hdf5_err) + case (SURF_PZ) + call h5ltmake_dataset_string_f(temp_group, "type", "Z Plane", hdf5_err) + case (SURF_PLANE) + call h5ltmake_dataset_string_f(temp_group, "type", "Plane", hdf5_err) + case (SURF_CYL_X) + call h5ltmake_dataset_string_f(temp_group, "type", "X Cylinder", hdf5_err) + case (SURF_CYL_Y) + call h5ltmake_dataset_string_f(temp_group, "type", "Y Cylinder", hdf5_err) + case (SURF_CYL_Z) + call h5ltmake_dataset_string_f(temp_group, "type", "Z Cylinder", hdf5_err) + case (SURF_SPHERE) + call h5ltmake_dataset_string_f(temp_group, "type", "Sphere", hdf5_err) + case (SURF_BOX_X) + case (SURF_BOX_Y) + case (SURF_BOX_Z) + case (SURF_BOX) + case (SURF_GQ) + call h5ltmake_dataset_string_f(temp_group, "type", "General Quadratic", hdf5_err) + end select + + ! Write coefficients for surface + dims(1) = size(s % coeffs) + call h5ltmake_dataset_double_f(temp_group, "coefficients", 1, dims, & + s % coeffs, hdf5_err) + + ! Write positive neighbors + if (allocated(s % neighbor_pos)) then + dims(1) = size(s % neighbor_pos) + call h5ltmake_dataset_int_f(temp_group, "neighbors_positive", 1, dims, & + s % neighbor_pos, hdf5_err) + end if + + ! Write negative neighbors + if (allocated(s % neighbor_neg)) then + dims(1) = size(s % neighbor_neg) + call h5ltmake_dataset_int_f(temp_group, "neighbors_negative", 1, dims, & + s % neighbor_neg, hdf5_err) + end if + + ! Write boundary condition + select case (s % bc) + case (BC_TRANSMIT) + call h5ltmake_dataset_string_f(temp_group, "boundary_condition", & + "transmission", hdf5_err) + case (BC_VACUUM) + call h5ltmake_dataset_string_f(temp_group, "boundary_condition", & + "vacuum", hdf5_err) + case (BC_REFLECT) + call h5ltmake_dataset_string_f(temp_group, "boundary_condition", & + "reflective", hdf5_err) + case (BC_PERIODIC) + call h5ltmake_dataset_string_f(temp_group, "boundary_condition", & + "periodic", hdf5_err) + end select + + ! Close group for i-th surface + call h5gclose_f(temp_group, hdf5_err) + end do + + call h5gclose_f(surface_group, hdf5_err) + + ! ========================================================================== + ! WRITE INFORMATION ON UNIVERSES + + call h5gcreate_f(geometry_group, "universes", universe_group, hdf5_err) + + ! Write information on each universe + do i = 1, n_universes + u => universes(i) + + ! Create group for i-th universe + call h5gcreate_f(universe_group, "universe " // trim(to_str(u % id)), & + temp_group, hdf5_err) + + ! Write list of cells in this universe + if (u % n_cells > 0) then + dims(1) = u % n_cells + call h5ltmake_dataset_int_f(temp_group, "cells", 1, dims, & + u % cells, hdf5_err) + end if + + ! Close group for i-th universe + call h5gclose_f(temp_group, hdf5_err) + end do + + call h5gclose_f(universe_group, hdf5_err) + + ! ========================================================================== + ! WRITE INFORMATION ON LATTICES + + call h5gcreate_f(geometry_group, "lattices", lattice_group, hdf5_err) + + ! Write information on each lattice + do i = 1, n_lattices + l => lattices(i) + + ! Create group for i-th lattice + call h5gcreate_f(lattice_group, "lattice " // trim(to_str(l % id)), & + temp_group, hdf5_err) + + ! Write lattice type + select case(l % type) + case (LATTICE_RECT) + call h5ltmake_dataset_string_f(temp_group, "type", "rectangular", hdf5_err) + case (LATTICE_HEX) + call h5ltmake_dataset_string_f(temp_group, "type", "hexagonal", hdf5_err) + end select + + ! Write lattice dimensions, lower left corner, and width of element + dims(1) = 2 + call h5ltmake_dataset_int_f(temp_group, "n_elements", 1, dims, & + (/ l % n_x, l % n_y /), hdf5_err) + call h5ltmake_dataset_double_f(temp_group, "lower_left", 1, dims, & + (/ l % x0, l % y0 /), hdf5_err) + call h5ltmake_dataset_double_f(temp_group, "element_width", 1, dims, & + (/ l % width_x, l % width_y /), hdf5_err) + + ! Write lattice elements + allocate(lattice_universes(l % n_x, l % n_y)) + do j = 1, l % n_x + do k = 1, l % n_y + lattice_universes(j,k) = universes(l % element(j,k)) % id + end do + end do + dims2 = (/ l % n_x, l % n_y /) + call h5ltmake_dataset_int_f(temp_group, "elements", 2, dims2, & + lattice_universes, hdf5_err) + deallocate(lattice_universes) + + ! Close group for i-th lattice + call h5gclose_f(temp_group, hdf5_err) + end do + + call h5gclose_f(lattice_group, hdf5_err) ! Close geometry group call h5gclose_f(geometry_group, hdf5_err) @@ -139,16 +339,65 @@ contains subroutine hdf5_write_materials() - integer :: rank = 1 - integer(HSIZE_T) :: dims(1) = (/1/) + integer :: i + integer(HSIZE_T) :: dims(1) + integer(HSIZE_T) :: size_string = 12 integer(HID_T) :: materials_group + integer(HID_T) :: temp_group + integer(HID_T) :: string_type + integer(HID_T) :: dspace_id + integer(HID_T) :: dset_id + type(Material), pointer :: m => null() ! Create group for materials call h5gcreate_f(hdf5_output_file, "/materials", materials_group, hdf5_err) ! Use H5LT interface to write number of materials - call h5ltmake_dataset_int_f(materials_group, "n_materials", & - rank, dims, (/ n_materials /), hdf5_err) + call hdf5_make_integer(materials_group, "n_materials", n_materials) + + ! Write information on each material + do i = 1, n_materials + m => materials(i) + + ! Create group for i-th universe + call h5gcreate_f(materials_group, "material " // trim(to_str(m % id)), & + temp_group, hdf5_err) + + ! Write atom density with units + call hdf5_make_double(temp_group, "atom_density", m % density) + call h5ltset_attribute_string_f(temp_group, "atom_density", & + "units", "atom/barn-cm", hdf5_err) + + ! Create string type of length 12 + call h5tcopy_f(H5T_C_S1, string_type, hdf5_err) + call h5tset_size_f(string_type, size_string, hdf5_err) + + ! Create dataspace and dataset for writing nuclides + dims(1) = m % n_nuclides + call h5screate_simple_f(1, dims, dspace_id, hdf5_err) + call h5dcreate_f(temp_group, "nuclides", string_type, dspace_id, & + dset_id, hdf5_err) + + ! Write list of nuclides + call h5dwrite_f(dset_id, string_type, m % names, dims, hdf5_err) + + ! Close dataspace and dataset for nuclides + call h5dclose_f(dset_id, hdf5_err) + call h5sclose_f(dspace_id, hdf5_err) + + ! Write atom densities + call h5ltmake_dataset_double_f(temp_group, "nuclide_densities", 1, & + dims, m % atom_density, hdf5_err) + + ! Write S(a,b) information if present + if (m % has_sab_table) then + call h5ltmake_dataset_string_f(temp_group, "sab_table", & + m % sab_name, hdf5_err) + end if + + ! Close group for i-th material + call h5gclose_f(temp_group, hdf5_err) + end do ! Close materials group call h5gclose_f(materials_group, hdf5_err) @@ -161,8 +410,6 @@ contains subroutine hdf5_write_timing() - integer :: rank = 1 - integer(HSIZE_T) :: dims(1) = (/1/) integer(HID_T) :: timing_group integer(8) :: total_particles real(8) :: speed @@ -171,30 +418,18 @@ contains call h5gcreate_f(hdf5_output_file, "/timing", timing_group, hdf5_err) ! Write timing data - call h5ltmake_dataset_double_f(timing_group, "time_initialize", & - rank, dims, (/ time_initialize % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_read_xs", & - rank, dims, (/ time_read_xs % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_unionize", & - rank, dims, (/ time_unionize % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_compute", & - rank, dims, (/ time_compute % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_intercycle", & - rank, dims, (/ time_intercycle % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_tallies", & - rank, dims, (/ time_ic_tallies % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_sample", & - rank, dims, (/ time_ic_sample % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_sendrecv", & - rank, dims, (/ time_ic_sendrecv % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_rebuild", & - rank, dims, (/ time_ic_rebuild % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_inactive", & - rank, dims, (/ time_inactive % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_active", & - rank, dims, (/ time_active % elapsed /), hdf5_err) - call h5ltmake_dataset_double_f(timing_group, "time_total", & - rank, dims, (/ time_total % elapsed /), hdf5_err) + call hdf5_make_double(timing_group, "time_initialize", time_initialize % elapsed) + call hdf5_make_double(timing_group, "time_read_xs", time_read_xs % elapsed) + call hdf5_make_double(timing_group, "time_unionize", time_unionize % elapsed) + call hdf5_make_double(timing_group, "time_compute", time_compute % elapsed) + call hdf5_make_double(timing_group, "time_intercycle", time_intercycle % elapsed) + call hdf5_make_double(timing_group, "time_tallies", time_ic_tallies % elapsed) + call hdf5_make_double(timing_group, "time_sample", time_ic_sample % elapsed) + call hdf5_make_double(timing_group, "time_sendrecv", time_ic_sendrecv % elapsed) + call hdf5_make_double(timing_group, "time_rebuild", time_ic_rebuild % elapsed) + call hdf5_make_double(timing_group, "time_inactive", time_inactive % elapsed) + call hdf5_make_double(timing_group, "time_active", time_active % elapsed) + call hdf5_make_double(timing_group, "time_total", time_total % elapsed) ! Add descriptions to timing data call h5ltset_attribute_string_f(timing_group, "time_initialize", & @@ -225,8 +460,7 @@ contains ! Write calculation rate total_particles = n_particles * n_cycles speed = real(total_particles) / time_compute % elapsed - call h5ltmake_dataset_double_f(timing_group, "neutrons_per_second", & - rank, dims, (/ speed /), hdf5_err) + call hdf5_make_double(timing_group, "neutrons_per_second", speed) ! Close timing group call h5gclose_f(timing_group, hdf5_err) @@ -247,6 +481,42 @@ contains end subroutine hdf5_close_output +!=============================================================================== +! HDF5_MAKE_INTEGER +!=============================================================================== + + subroutine hdf5_make_integer(group, name, buffer) + + integer(HID_T), intent(in) :: group + character(*), intent(in) :: name + integer, intent(in) :: buffer + + integer :: rank = 1 + integer(HSIZE_T) :: dims(1) = (/1/) + + call h5ltmake_dataset_int_f(group, name, rank, dims, & + (/ buffer /), hdf5_err) + + end subroutine hdf5_make_integer + +!=============================================================================== +! HDF5_MAKE_DOUBLE +!=============================================================================== + + subroutine hdf5_make_double(group, name, buffer) + + integer(HID_T), intent(in) :: group + character(*), intent(in) :: name + real(8), intent(in) :: buffer + + integer :: rank = 1 + integer(HSIZE_T) :: dims(1) = (/1/) + + call h5ltmake_dataset_double_f(group, name, rank, dims, & + (/ buffer /), hdf5_err) + + end subroutine hdf5_make_double + #endif end module hdf5_interface