PB3D
[2.45]
Ideal linear high-n MHD stability in 3-D
|
Utilities pertaining to HDF5 and XDMF. More...
Functions/Subroutines | |
integer function, public | set_1d_vars (lim_tot, lim_loc, space_id, c_plist_id) |
Sets up the chunk property list and/or the 1D hyperslabs that correspond to a local subset of lim_tot given by the limits lim_loc . More... | |
integer function, public | probe_hdf5_group (HDF5_name, group_name, group_exists) |
Probe HDF5 file for group existence. More... | |
integer function, public | list_all_vars_in_group (group_id) |
Lists all variables in a HDF5 group. More... | |
Variables | |
logical, public | debug_set_1d_vars = .false. |
set to true to debug set_1D_vars More... | |
Utilities pertaining to HDF5 and XDMF.
integer function, public hdf5_utilities::list_all_vars_in_group | ( | integer(hid_t), intent(in) | group_id | ) |
Lists all variables in a HDF5 group.
[in] | group_id | group identifier |
Definition at line 319 of file HDF5_utilities.f90.
integer function, public hdf5_utilities::probe_hdf5_group | ( | character(len=*), intent(in) | HDF5_name, |
character(len=*), intent(in) | group_name, | ||
logical, intent(inout) | group_exists | ||
) |
Probe HDF5 file for group existence.
[in] | hdf5_name | name of HDF5 file |
[in] | group_name | name of group to probe for |
[in,out] | group_exists | whether group exists |
Definition at line 251 of file HDF5_utilities.f90.
integer function, public hdf5_utilities::set_1d_vars | ( | integer, dimension(:,:), intent(in) | lim_tot, |
integer, dimension(:,:), intent(in) | lim_loc, | ||
integer(hid_t), intent(in), optional | space_id, | ||
integer(hid_t), intent(inout), optional | c_plist_id | ||
) |
Sets up the chunk property list and/or the 1D hyperslabs that correspond to a local subset of lim_tot
given by the limits lim_loc
.
An example for the hyperslab is given in 3 dimensions with size n = 5,3,3:
The 1D equivalent can be represented as
Now, dimension 1 will only allow for the following elements (given by -
):
dimension 2 will only allow for the following elements:
dimension 3 will only allow for the following elements:
so that the total is given by:
In practice, it is easiest to work with a full selection, where for every dimension the ranges that do not correspond to the range of that dimension are taken away. Clearly, if a dimension has the full range, nothing will be taken out:
block_i
= \(n_1 n_2 \cdots n_{i-1} (b_i-a_{i+1}) \) ,stride_i
= \(n_1 n_2 \cdots n_{i-1} (B_i-A_{i+1}) \) ,offset_i
= \(n_1 n_2 \cdots n_{i-1} (a_i-A_i) \) ,count_i
= \(n_{i+1} n_{i+2} \cdots n_N \) ,where \(a_i\) and \(b_i\) represent the local limits of dimension \(i\), \(A_i\) and \(B_i\) the total ones and the number of dimensions is \(N\), as can be verified.
The chunk property list for creation can be set up so that its size is equal to the largest dimensions of full range, or an integer part of that if it exceeds 4GB. Since the variables don't need to be used more than once, either nbytes
or nslots
could be set to 0 in an access property list, but this is currently not done.
[in] | lim_tot | total limits |
[in] | lim_loc | local limits |
[in] | space_id | dataspace identifier |
[in,out] | c_plist_id | chunk creation property list identifier |
Definition at line 81 of file HDF5_utilities.f90.
logical, public hdf5_utilities::debug_set_1d_vars = .false. |
set to true to debug set_1D_vars
Definition at line 24 of file HDF5_utilities.f90.