Switched to consistent use of output_interface

This commit is contained in:
Sterling Harper 2013-07-09 11:04:04 -04:00
parent 52233d0173
commit 06c3bd3690

View file

@ -1,160 +1,95 @@
!===============================================================================
! PARTICLE_TRACK handles output of particle tracks to disk.
!
! TODO: This module writes binary output files via a stream access i.e. it
! writes the particle's location to disk everytime write_particle_track is
! called. But HDF5 does not allow for a simple implimentation of stream access.
! An attempt was made to write HDF5 files in a stream fashion using an
! extendable dataset, but it always failed at the write call because of an out-
! of-bounds error. The extendable dataset code is commented out by enclosing it
! within "#if 0"/"#endif" blocks, and has been replaced by code that writes
! particle coordinates to an ever-expanding array which is written to disk in
! its entirety when finalize_particle_track is called.
!
! In the interest of consistency, maybe this module should be changed so that it
! uses the extendable HDF5 dataset, or it writes binary output all in one call
! using the output_interface (like the current HDF5 implimentation).
! PARTICLE_TRACK handles output of particle tracks (the paths taken by particles
! as they are transported through the geometry).
!===============================================================================
module particle_track
use constants
use global
use string, only: to_str
#if 0
use hdf5
#endif
#ifdef HDF5
use output_interface, only: file_create, file_close, write_data
#endif
use string, only: to_str
implicit none
#if 0
integer, private :: hdf5_err ! HDF error code
integer(HID_T), private :: track_dset ! dataset handle
integer(HID_T), private :: track_dspace ! dataspace handle
integer(HID_T), private :: track_fh ! HDF file handle
integer(HID_T), private :: cparms ! chunk parameters
integer(HSIZE_T), private :: n_tracks ! total number of tracks
#endif
#ifdef HDF5
character(MAX_FILE_LEN), private :: fname ! file name
integer, private :: n_tracks ! total number of tracks
real(8), private, allocatable :: coords(:,:) ! track coordinates
#endif
contains
!===============================================================================
! INITIALIZE_PARTICLE_TRACK opens a particle track output file.
! INITIALIZE_PARTICLE_TRACK
!===============================================================================
subroutine initialize_particle_track()
character(MAX_FILE_LEN) :: filename
#if 0
integer(HSIZE_T) :: dims(2), max_dims(2), i, j(1)
integer(HID_T) :: cparms
filename = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
// '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p % id)) &
// '.h5'
n_tracks = 0
! Create file.
call h5fcreate_f(filename, H5F_ACC_TRUNC_F, track_fh, hdf5_err)
! Create a dataspace with an unlimited max dimension.
dims = (/3, 1/)
i = 3
max_dims = (/i, H5S_UNLIMITED_F/)
call h5screate_simple_f(2, dims, track_dspace, hdf5_err, max_dims)
! Set dataspace chunking.
call h5pcreate_f(H5P_DATASET_CREATE_F, cparms, hdf5_err)
call h5pset_chunk_f(cparms, 2, dims, hdf5_err)
! Set fill value.
call h5pset_fill_value_f(cparms, H5T_NATIVE_DOUBLE, 0, hdf5_err)
! Create dataset.
call h5dcreate_f(track_fh, 'coordinates', H5T_NATIVE_DOUBLE, track_dspace, &
track_dset, hdf5_err, cparms)
#endif
#ifdef HDF5
filename = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
// '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p % id)) &
// '.h5'
n_tracks = 0
fname = filename
allocate(coords(3,1))
#else
filename = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
// '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p % id)) &
// '.binary'
open(UNIT=UNIT_TRACK, FILE=filename, ACTION="write", &
STATUS='replace', ACCESS='stream')
#endif
allocate(coords(1,1))
end subroutine initialize_particle_track
!===============================================================================
! WRITE_PARTICLE_TRACK outputs particle position to a binary file.
! WRITE_PARTICLE_TRACK copies particle position to an array.
!===============================================================================
subroutine write_particle_track()
#if 0
integer(HSIZE_T) :: dims(2), max_dims(2), offset(2), i
integer(HSIZE_T), parameter :: cnt(2) = (/3, 1/), write_dims(2) = (/3, 1/)
real(8) :: coords(3)
! There is another set of track coordinates. Incriment the counter.
n_tracks = n_tracks + 1
! Make the dataset the right size to fit the track coordinates.
i = 3
dims = (/i, n_tracks/)
call h5dset_extent_f(track_dset, dims, hdf5_err)
! Get the dataspace with the updated dimensions.
call h5dget_space_f(track_dset, track_dspace, hdf5_err)
! Select the hyperslab where the latest coordinates will be written.
i = 1
offset = (/i, n_tracks/)
if (n_tracks < 20) then
call h5sselect_hyperslab_f(track_dspace, H5S_SELECT_SET_F, offset, cnt, &
hdf5_err)
endif
! Write the coordinates to the dataset.
coords = p % coord0 % xyz
call h5dwrite_f(track_dset, H5T_NATIVE_DOUBLE, coords, &
write_dims, hdf5_err, file_space_id=track_dspace)
#endif
#ifdef HDF5
real(8) :: old_coords(3, n_tracks)
! Save the coordinates gathered in previous calls.
old_coords = coords
! Add another column to coords.
n_tracks = n_tracks + 1
deallocate(coords)
allocate(coords(3, n_tracks))
! Put the old coordinates back into the array.
coords(:, 1:n_tracks-1) = old_coords
! Write current coordinates into the newest column.
coords(:, n_tracks) = p % coord0 % xyz
#else
write(UNIT_TRACK) p % coord0 % xyz
#endif
end subroutine write_particle_track
!===============================================================================
! FINALIZE_PARTICLE_TRACK closes the particle track file.
! FINALIZE_PARTICLE_TRACK writes the particle track array to disk.
!
! output_interface currently does not support writing 2D binary arrays so there
! are two different versions of this subroutine; an HDF version which simply
! writes the array and a binary version which flattens the array into a 1D shape
! before writing.
!===============================================================================
subroutine finalize_particle_track()
#if 0
call h5dclose_f(track_dset, hdf5_err)
call h5sclose_f(track_dspace, hdf5_err)
call h5pclose_f(cparms, hdf5_err)
call h5fclose_f(track_fh, hdf5_err)
#endif
#ifdef HDF5
subroutine finalize_particle_track()
character(MAX_FILE_LEN) :: fname
fname = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
// '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p % id)) &
// '.h5'
call file_create(fname, 'serial')
call write_data(coords, 'coordinates', length=(/3, n_tracks/))
call file_close('serial')
deallocate(coords)
#else
close(UNIT=UNIT_TRACK)
#endif
end subroutine finalize_particle_track
#else
subroutine finalize_particle_track()
character(MAX_FILE_LEN) :: fname
integer :: i
real(8) :: flat_coords(3*n_tracks)
fname = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
// '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p % id)) &
// '.binary'
call file_create(fname, 'serial')
do i=1, n_tracks
flat_coords(3*i-2 : 3*i) = coords(:,i)
end do
call write_data(flat_coords, 'coordinates', length=3*n_tracks)
call file_close('serial')
deallocate(coords)
end subroutine finalize_particle_track
#endif
end module particle_track