Go to the documentation of this file.
43 character(LEN=*),
parameter,
private ::
dsetname_u =
"Units"
46 integer(HSIZE_T),
dimension(1),
private ::
dims_1d=(/1/)
49 integer(HSIZE_T),
dimension(1:2),
private ::
maxdims
69 integer(HSIZE_T),
dimension(1:2),
private ::
dims_y
165 character(len=9),
dimension(17),
private ::
ts_names = &
166 (/
"tau_ga ",
"tau_ag ",
"tau_ng ",
"tau_gn ",
"tau_pg ",
"tau_gp ",
"tau_np ",&
167 "tau_pn ",
"tau_an ",
"tau_na ",
"tau_ap ",
"tau_pa ",
"tau_beta ",
"tau_alpha", &
168 "tau_nfiss",
"tau_sfiss",
"tau_bfiss"/)
182 real(
r_kind),
private,
dimension(:,:),
allocatable :: & !< Storage for the chunk which is later written to the file.
192 character(len=11),
dimension(10),
private ::
en_names = &
193 (/
"engen_tot ",
"S_src ",
"engen_ng_gn",
"engen_pg_gp",
"engen_ag_ga",&
194 "engen_np_pn",
"engen_an_na",
"engen_ap_pa",
"engen_beta ", &
199 character(len=14),
private :: & !< Name of the decay energy splitted by parent
205 integer(HID_T),
private :: & !< Dataset identifier
314 call h5open_f(i_stat)
316 if (i_stat .ne. 0)
call raise_exception(
"Could not initialize fortran predefined datatypes",&
317 "init_hdf5_module",240003)
321 if (i_stat .ne. 0)
then
323 "Close the file if you have opened it.",&
324 "init_hdf5_module",240004)
373 if (i_stat .ne. 0)
then
389 if (i_stat .ne. 0)
then
391 "init_nu_loss_gain",240005)
417 real(
r_kind),
intent(in) :: time
418 real(
r_kind),
intent(in) :: temp
419 real(
r_kind),
intent(in) :: dens
420 real(
r_kind),
intent(in) :: rad
421 real(
r_kind),
intent(in) :: the
422 real(
r_kind),
intent(in) :: heat
423 real(
r_kind),
intent(in) :: bet
426 nut = the + heat + bet
480 logical,
intent(in) :: include_decay
485 integer :: filter_info,filter_info_both
486 integer,
dimension(net_size) :: helper
491 if (i_stat .ne. 0)
then
503 if (include_decay .eqv. .true.)
then
507 maxdims = (/h5s_unlimited_f, h5s_unlimited_f/)
509 if (i_stat .ne. 0)
then
521 CALL h5zfilter_avail_f(h5z_filter_deflate_f, avail, i_stat)
523 if (.not. avail)
then
526 call h5zget_filter_info_f(h5z_filter_deflate_f, filter_info, i_stat)
528 filter_info_both=ior(h5z_filter_encode_enabled_f,h5z_filter_decode_enabled_f)
529 if (filter_info .ne. filter_info_both)
then
530 call raise_exception(
"gzip filter not available for encoding and decoding.",&
535 call h5pcreate_f(h5p_dataset_create_f,
crp_list, i_stat)
536 if (i_stat .ne. 0)
then
541 call h5pset_deflate_f(
crp_list, 9, i_stat)
542 if (i_stat .ne. 0)
then
548 if (i_stat .ne. 0)
then
556 if (i_stat .ne. 0)
then
562 if (i_stat .ne. 0)
then
568 if (i_stat .ne. 0)
then
574 if (i_stat .ne. 0)
then
580 if (i_stat .ne. 0)
then
586 if (i_stat .ne. 0)
then
592 if (i_stat .ne. 0)
then
598 if (i_stat .ne. 0)
then
604 if (i_stat .ne. 0)
then
610 if (i_stat .ne. 0)
then
615 if (i_stat .ne. 0)
then
645 if (i_stat /= 0)
call raise_exception(
'Allocation failed.',
"init_engen",&
657 subroutine extend_engen(time,engen_tot, heat,engen_ng_gn, engen_pg_gp, &
658 engen_ag_ga, engen_np_pn, engen_an_na, &
659 engen_ap_pa, engen_beta, engen_fiss, &
660 det_bet_engen, det_ag_engen, det_pg_engen, &
661 det_ng_engen, det_pn_engen, det_ap_engen, &
662 det_an_engen, det_other_engen, det_fiss_engen, &
667 real(
r_kind),
intent(in) :: time
668 real(
r_kind),
intent(in) :: engen_tot, heat, engen_ng_gn, engen_pg_gp
669 real(
r_kind),
intent(in) :: engen_ag_ga, engen_np_pn, engen_an_na
670 real(
r_kind),
intent(in) :: engen_ap_pa, engen_beta, engen_fiss
671 real(
r_kind),
dimension(net_size),
intent(in) :: det_bet_engen
672 real(
r_kind),
dimension(net_size),
intent(in) :: det_ag_engen
673 real(
r_kind),
dimension(net_size),
intent(in) :: det_pg_engen
674 real(
r_kind),
dimension(net_size),
intent(in) :: det_ng_engen
675 real(
r_kind),
dimension(net_size),
intent(in) :: det_pn_engen
676 real(
r_kind),
dimension(net_size),
intent(in) :: det_ap_engen
677 real(
r_kind),
dimension(net_size),
intent(in) :: det_an_engen
678 real(
r_kind),
dimension(net_size),
intent(in) :: det_other_engen
679 real(
r_kind),
dimension(net_size),
intent(in) :: det_fiss_engen
680 real(
r_kind),
dimension(net_size),
intent(in) :: det_tot_engen
681 logical,
intent(in) :: include_decay
689 engen_ag_ga, engen_np_pn, engen_an_na, &
690 engen_ap_pa, engen_beta, engen_fiss/)
692 if (include_decay .eqv. .true.)
then
724 logical,
intent(in) :: include_decay
736 if (include_decay .eqv. .true.)
then
764 real(
r_kind),
intent(in) :: time
765 real(
r_kind),
intent(in) :: temp
766 real(
r_kind),
intent(in) :: dens
767 integer,
intent(in) :: n_flows
769 real(
r_kind),
dimension(net_size),
intent(in) :: y
773 integer :: alloc_stat
775 integer,
dimension(:),
allocatable :: n_in_nr,p_in_nr
776 integer,
dimension(:),
allocatable :: n_out_nr,p_out_nr
777 real(
r_kind),
dimension(:),
allocatable :: flow_in,flow_out
778 real(
r_kind),
dimension(:),
allocatable :: y_in,y_out
779 integer(HID_T) :: tmp_group_id
780 real(
r_kind),
dimension(1) :: helper_time
781 real(
r_kind),
dimension(1) :: helper_temp
782 real(
r_kind),
dimension(1) :: helper_dens
791 flow_diff =
flow(i)%fwd -
flow(i)%bwd
793 if (flow_diff .ne. 0) count = count + 1
797 allocate(n_out_nr(count),stat=alloc_stat)
798 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'n_out_nr' failed!",
"extend_flow",240001)
799 allocate(p_out_nr(count),stat=alloc_stat)
800 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'p_out_nr' failed!",
"extend_flow",240001)
801 allocate(flow_out(count),stat=alloc_stat)
802 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'flow_out' failed!",
"extend_flow",240001)
803 allocate( n_in_nr(count),stat=alloc_stat)
804 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'n_in_nr' failed!",
"extend_flow",240001)
805 allocate( p_in_nr(count),stat=alloc_stat)
806 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'p_in_nr' failed!",
"extend_flow",240001)
807 allocate( flow_in(count),stat=alloc_stat)
808 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'flow_in' failed!",
"extend_flow",240001)
809 allocate( y_in(count) ,stat=alloc_stat)
810 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'Y_in' failed!",
"extend_flow",240001)
811 allocate( y_out(count) ,stat=alloc_stat)
812 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of 'Y_out' failed!",
"extend_flow",240001)
817 flow_diff =
flow(i)%fwd -
flow(i)%bwd
819 if (flow_diff .eq. 0) cycle
826 flow_out(count) =
flow(i)%bwd
827 flow_in(count) =
flow(i)%fwd
828 y_in(count) = y(
flow(i)%iin)
829 y_out(count) = y(
flow(i)%iout)
832 helper_time(1) = time
833 helper_temp(1) = temp
834 helper_dens(1) = dens
839 flow_in(:) = flow_in(:) - flow_out(:)
852 call h5gclose_f(tmp_group_id , i_stat)
856 deallocate(n_out_nr,stat=alloc_stat)
857 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'n_out_nr' failed!",
"extend_flow",240002)
858 deallocate(p_out_nr,stat=alloc_stat)
859 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'p_out_nr' failed!",
"extend_flow",240002)
860 deallocate(flow_out,stat=alloc_stat)
861 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'flow_out' failed!",
"extend_flow",240002)
862 deallocate( n_in_nr,stat=alloc_stat)
863 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'n_in_nr' failed!",
"extend_flow",240002)
864 deallocate( p_in_nr,stat=alloc_stat)
865 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'p_in_nr' failed!",
"extend_flow",240002)
866 deallocate( flow_in,stat=alloc_stat)
867 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'flow_in' failed!",
"extend_flow",240002)
868 deallocate(y_in ,stat=alloc_stat)
869 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'Y_in' failed!",
"extend_flow",240002)
870 deallocate(y_out,stat=alloc_stat)
871 if (alloc_stat .ne. 0)
call raise_exception(
"Deallocation of 'Y_out' failed!",
"extend_flow",240002)
888 integer,
dimension(net_size) :: helper
890 integer :: filter_info,filter_info_both
894 if (i_stat .ne. 0)
then
902 maxdims = (/h5s_unlimited_f, h5s_unlimited_f/)
904 if (i_stat .ne. 0)
then
911 call h5zfilter_avail_f(h5z_filter_deflate_f, avail, i_stat)
913 if (.not. avail)
then
916 call h5zget_filter_info_f(h5z_filter_deflate_f, filter_info, i_stat)
918 filter_info_both=ior(h5z_filter_encode_enabled_f,h5z_filter_decode_enabled_f)
919 if (filter_info .ne. filter_info_both)
then
920 call raise_exception(
"gzip filter not available for encoding and decoding.",
"init_snaps",240009)
924 call h5pcreate_f(h5p_dataset_create_f,
crp_list, i_stat)
925 if (i_stat .ne. 0)
then
931 call h5pset_deflate_f(
crp_list, 9, i_stat)
932 if (i_stat .ne. 0)
then
938 if (i_stat .ne. 0)
then
946 if (i_stat .ne. 0)
then
951 if (i_stat .ne. 0)
then
993 if (i_stat .ne. 0)
then
995 "init_mainout",240005)
1040 integer,
dimension(track_nuclei_nr) :: helper
1045 if (i_stat .ne. 0)
then
1047 "init_track_nuclei",240005)
1054 maxdims = (/h5s_unlimited_f, h5s_unlimited_f/)
1056 if (i_stat .ne. 0)
then
1058 "init_snaps",240006)
1062 call h5pcreate_f(h5p_dataset_create_f,
crp_list, i_stat)
1063 if (i_stat .ne. 0)
then
1065 "init_snaps",240010)
1069 if (i_stat .ne. 0)
then
1071 "init_snaps",240012)
1077 if (i_stat .ne. 0)
then
1079 "init_snaps",240006)
1082 if (i_stat .ne. 0)
then
1084 "init_snaps",240007)
1089 if (i_stat .ne. 0)
then
1091 "init_track_nuclei",240001)
1119 if (i_stat .ne. 0)
then
1121 "init_track_nuclei",240001)
1137 if (i_stat .ne. 0)
then
1139 "init_av_timescales",240005)
1164 integer,
parameter :: length=10
1165 character(len=20),
dimension(length,2) :: quantity
1166 integer(HID_T) :: dset_id
1167 integer(HSIZE_T),
dimension(2) :: dims_static = (/length,2/)
1169 INTEGER(HID_T) :: strtype
1170 integer(SIZE_T) :: typesize
1172 quantity(1 ,1) =
"Abundances" ; quantity(1 ,2) =
"-"
1173 quantity(2 ,1) =
"Atomic number" ; quantity(2 ,2) =
"-"
1174 quantity(3 ,1) =
"Density" ; quantity(3 ,2) =
"g/ccm"
1175 quantity(4 ,1) =
"Electron fraction" ; quantity(4 ,2) =
"-"
1176 quantity(5 ,1) =
"Entropy" ; quantity(5 ,2) =
"kB/nuc"
1177 quantity(6 ,1) =
"Mass number" ; quantity(6 ,2) =
"-"
1178 quantity(7 ,1) =
"Neutron number" ; quantity(7 ,2) =
"-"
1179 quantity(8 ,1) =
"Radius" ; quantity(8 ,2) =
"Km"
1180 quantity(9 ,1) =
"Temperature" ; quantity(9 ,2) =
"GK"
1181 quantity(10,1) =
"Time" ; quantity(10,2) =
"s"
1185 call h5tcopy_f(h5t_native_character, strtype, i_stat)
1187 call h5tset_size_f(strtype, typesize, i_stat)
1190 call h5screate_simple_f(2, dims_static,
dataspace, i_stat)
1198 call h5dwrite_f(dset_id, strtype, quantity, dims_static, i_stat)
1201 call h5dclose_f(dset_id, i_stat)
1220 real(
r_kind),
dimension(net_size),
intent(in) :: y
1221 real(
r_kind),
dimension(:),
allocatable :: y_tmpa,x_tmpa
1222 real(
r_kind),
dimension(:),
allocatable :: y_tmp,x_tmp
1223 integer,
dimension(:),
allocatable :: a_tmp,z_tmp
1231 if (i_stat .ne. 0)
then
1233 "write_finab",240005)
1237 if (i_stat .ne. 0)
then
1239 "write_finab",240005)
1247 if (y(i) .gt. 1d-25)
then
1253 allocate(y_tmp(count))
1254 allocate(a_tmp(count))
1255 allocate(z_tmp(count))
1260 if (y(i) .gt. 1d-25)
then
1263 a_tmp(count) =
isotope(i)%mass
1264 z_tmp(count) =
isotope(i)%p_nr
1284 if (i_stat .ne. 0)
then
1286 "write_finab",240005)
1290 allocate(y_tmpa(mval))
1300 if (y_tmpa(i) .gt. 1d-20)
then
1306 allocate(y_tmp(count))
1307 allocate(a_tmp(count))
1312 if (y_tmpa(i) .gt. 1d-20)
then
1315 y_tmp(count) = y_tmpa(i)
1334 if (i_stat .ne. 0)
then
1336 "write_finab",240005)
1341 allocate(y_tmpa(mval))
1342 allocate(x_tmpa(mval))
1354 if (y_tmpa(i) .gt. 1d-20)
then
1360 allocate(y_tmp(count))
1361 allocate(x_tmp(count))
1362 allocate(z_tmp(count))
1370 if (y_tmpa(i) .gt. 1d-20)
then
1372 y_tmp(count) = y_tmpa(i)
1373 x_tmp(count) = x_tmpa(i)
1399 character(len=*),
intent(in) :: dsetname
1400 integer(HID_T),
intent(in) :: group_id
1402 integer(HID_T) :: dset_id
1403 integer(HSIZE_T),
dimension(1) :: maxdims_1d
1407 maxdims_1d = h5s_unlimited_f
1411 if (i_stat .ne. 0)
then
1413 "create_1d_dataset",240006)
1417 call h5pcreate_f(h5p_dataset_create_f,
crp_list, i_stat)
1418 if (i_stat .ne. 0)
then
1420 "create_1d_dataset",240010)
1424 if (i_stat .ne. 0)
then
1426 "create_1d_dataset",240012)
1430 call h5dcreate_f(group_id, dsetname, h5t_native_double,
dataspace, &
1432 if (i_stat .ne. 0)
then
1434 "create_1d_dataset",240006)
1438 if (i_stat .ne. 0)
then
1440 "create_1d_dataset",240023)
1461 character(len=*),
intent(in) :: dsetname
1462 integer(HID_T),
intent(in) :: group_id
1464 integer(HID_T) :: dset_id
1465 integer(HSIZE_T),
dimension(1) :: maxdims_1d
1469 maxdims_1d = h5s_unlimited_f
1473 if (i_stat .ne. 0)
then
1475 "create_1d_int_dataset",240006)
1479 call h5pcreate_f(h5p_dataset_create_f,
crp_list, i_stat)
1480 if (i_stat .ne. 0)
then
1482 "create_1d_int_dataset",240010)
1486 if (i_stat .ne. 0)
then
1488 "create_1d_int_dataset",240012)
1493 call h5dcreate_f(group_id, dsetname, h5t_native_integer,
dataspace, &
1495 if (i_stat .ne. 0)
then
1497 "create_1d_int_dataset",240006)
1501 if (i_stat .ne. 0)
then
1503 "create_1d_int_dataset",240007)
1521 integer,
intent(in) :: length
1522 integer,
dimension(length),
intent(in) :: data
1523 character(len=*),
intent(in) :: dsetname
1524 integer(HID_T) :: dset_id
1525 integer(HID_T) :: group_id
1526 integer(HSIZE_T),
dimension(1) :: maxdims_1d
1527 integer(HSIZE_T),
dimension(1) :: dims_static
1531 dims_static = length
1534 call h5screate_simple_f(1, dims_static,
dataspace, i_stat)
1535 if (i_stat .ne. 0)
then
1537 "create_constant_1d_int_arrays",240006)
1541 call h5dcreate_f(group_id, dsetname, h5t_native_integer,
dataspace, &
1543 if (i_stat .ne. 0)
then
1545 "create_constant_1d_int_arrays",240010)
1549 if (i_stat .ne. 0)
then
1551 "create_constant_1d_int_arrays",240007)
1554 call h5dwrite_f(dset_id, h5t_native_integer,
data, dims_static, i_stat)
1555 if (i_stat .ne. 0)
then
1557 "create_constant_1d_int_arrays",240013)
1561 call h5dclose_f(dset_id, i_stat)
1562 if (i_stat .ne. 0)
then
1564 "create_constant_1d_int_arrays",240007)
1578 integer,
intent(in) :: length
1579 real(
r_kind),
dimension(length),
intent(in) :: data
1580 character(len=*),
intent(in) :: dsetname
1581 integer(HID_T) :: dset_id
1582 integer(HID_T) :: group_id
1583 integer(HSIZE_T),
dimension(1) :: maxdims_1d
1584 integer(HSIZE_T),
dimension(1) :: dims_static
1588 dims_static = length
1591 call h5screate_simple_f(1, dims_static,
dataspace, i_stat)
1592 if (i_stat .ne. 0)
then
1594 "create_constant_1d_int_arrays",240006)
1598 call h5dcreate_f(group_id, dsetname, h5t_native_double,
dataspace, &
1600 if (i_stat .ne. 0)
then
1602 "create_constant_1d_int_arrays",240010)
1606 if (i_stat .ne. 0)
then
1608 "create_constant_1d_int_arrays",240007)
1611 call h5dwrite_f(dset_id, h5t_native_double,
data, dims_static, i_stat)
1612 if (i_stat .ne. 0)
then
1614 "create_constant_1d_int_arrays",240013)
1618 call h5dclose_f(dset_id, i_stat)
1619 if (i_stat .ne. 0)
then
1621 "create_constant_1d_int_arrays",240007)
1637 integer,
intent(in) :: chunksize
1638 real(
r_kind),
dimension(chunksize),
intent(in) ::
value
1639 integer(HID_T),
intent(in) :: dset_id
1640 integer,
intent(in) :: in_size
1641 integer(HID_T) :: memspace
1643 integer(HSIZE_T),
dimension(1) :: offset
1644 integer(HSIZE_T),
dimension(1) :: count
1645 integer(HSIZE_T),
dimension(1) :: new_size
1648 new_size(1) = in_size+chunksize
1650 call h5dset_extent_f(dset_id, new_size, i_stat)
1651 if (i_stat .ne. 0)
then
1653 "extend_1d_dataset",240014)
1656 call h5screate_simple_f (1,
dims_1d, memspace, i_stat)
1657 if (i_stat .ne. 0)
then
1659 "extend_1d_dataset",240015)
1665 count(1) = chunksize
1668 call h5dget_space_f(dset_id,
dataspace, i_stat)
1669 call h5sselect_hyperslab_f(
dataspace, h5s_select_set_f, &
1670 offset, count, i_stat)
1671 if (i_stat .ne. 0)
then
1673 "extend_1d_dataset",240016)
1676 call h5dwrite_f(dset_id, h5t_native_double,
value,
dims_1d, i_stat, &
1678 if (i_stat .ne. 0)
then
1680 "extend_1d_dataset",240013)
1695 integer,
intent(in) :: chunksize
1696 integer,
dimension(chunksize),
intent(in) ::
value
1697 integer(HID_T),
intent(in) :: dset_id
1698 integer,
intent(in) :: in_size
1699 integer(HID_T) :: memspace
1701 integer(HSIZE_T),
dimension(1) :: offset
1702 integer(HSIZE_T),
dimension(1) :: count
1703 integer(HSIZE_T),
dimension(1) :: new_size
1706 new_size(1) = in_size+chunksize
1708 call h5dset_extent_f(dset_id, new_size, i_stat)
1709 if (i_stat .ne. 0)
then
1711 "extend_1d_int_dataset",240014)
1714 call h5screate_simple_f (1,
dims_1d, memspace, i_stat)
1715 if (i_stat .ne. 0)
then
1717 "extend_1d_int_dataset",240015)
1723 count(1) = chunksize
1726 call h5dget_space_f(dset_id,
dataspace, i_stat)
1727 call h5sselect_hyperslab_f(
dataspace, h5s_select_set_f, &
1728 offset, count, i_stat)
1729 if (i_stat .ne. 0)
then
1731 "extend_1d_int_dataset",240016)
1735 call h5dwrite_f(dset_id, h5t_native_integer,
value,
dims_1d, i_stat, &
1737 if (i_stat .ne. 0)
then
1739 "extend_1d_int_dataset",240013)
1751 tau_pn, tau_an, tau_na, tau_ap, tau_pa, tau_beta, tau_alpha, &
1752 tau_nfiss, tau_sfiss, tau_bfiss)
1754 real(
r_kind),
intent(in) :: time
1755 real(
r_kind),
intent(in) :: tau_ng, tau_gn, tau_pg, tau_gp, tau_np
1756 real(
r_kind),
intent(in) :: tau_ga, tau_ag, tau_ap, tau_pa
1757 real(
r_kind),
intent(in) :: tau_pn, tau_an, tau_na, tau_beta, tau_alpha
1758 real(
r_kind),
intent(in) :: tau_nfiss, tau_sfiss, tau_bfiss
1766 tau_pn, tau_an, tau_na, tau_ap, tau_pa, tau_beta, tau_alpha, &
1767 tau_nfiss, tau_sfiss, tau_bfiss/)
1809 real(
r_kind),
dimension(net_size),
intent(in) :: y_array
1810 real(
r_kind),
intent(in) :: time
1868 real(
r_kind),
dimension(net_size),
intent(in) :: y_array
1869 real(
r_kind),
intent(in) :: time
1919 integer,
intent(in) :: cnt
1920 real(
r_kind),
intent(in) :: time
1921 real(
r_kind),
intent(in) :: temp
1922 real(
r_kind),
intent(in) :: dens
1923 real(
r_kind),
intent(in) :: entr
1924 real(
r_kind),
intent(in) :: rad
1925 real(
r_kind),
dimension(net_size),
intent(in) :: y
1927 real(
r_kind) :: ye,yneu,ypro,yhe4,ylight,yheavies,ysum,yasum,yzsum,abar,zbar
1934 ylight= max(1.0d-99,sum(y(
ipro+1:
ihe4-1)))
1945 if (
ineu.ge.1) yneu= max(1.0d-99,y(
ineu))
1948 if (
ipro.ge.1) ypro= max(1.0d-99,y(
ipro))
1951 if (
ihe4.ge.1) yhe4= max(1.0d-99,y(
ihe4))
2048 integer(HID_T),
intent(in) :: dset_id
2049 integer,
intent(in) :: nr_nuc
2050 integer,
intent(in) :: old_size
2051 integer,
intent(in) :: chunksize
2052 real(
r_kind),
dimension(nr_nuc,chunksize),
intent(in) :: y_array
2053 integer(HSIZE_T),
dimension(2) :: new_size
2054 integer(HSIZE_T),
dimension(1:2) :: offset
2055 integer(HSIZE_T),
dimension(1:2) :: count
2057 integer(HID_T) :: memspace
2060 new_size(1:2) = (/nr_nuc,old_size+chunksize/)
2062 dims_y = (/nr_nuc,chunksize/)
2064 call h5dset_extent_f(dset_id, new_size, i_stat)
2065 call h5screate_simple_f (2,
dims_y, memspace, i_stat)
2068 offset(1:2) = (/0,old_size/)
2070 count(1:2) = (/nr_nuc,chunksize/)
2073 call h5dget_space_f(dset_id,
dataspace, i_stat)
2074 call h5sselect_hyperslab_f(
dataspace, h5s_select_set_f, &
2075 offset, count, i_stat)
2077 call h5dwrite_f(dset_id, h5t_native_double, y_array,
dims_y, i_stat, &
2108 if (i_stat .ne. 0)
call raise_exception(
"Unable to close snapshot group.",&
2109 "hdf5_module_finalize",240017)
2112 if (i_stat .ne. 0)
call raise_exception(
"Unable to close abundance dataset.",&
2113 "hdf5_module_finalize",240007)
2115 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2116 "hdf5_module_finalize",240007)
2123 if (i_stat .ne. 0)
call raise_exception(
"Unable to close mainout group.",&
2124 "hdf5_module_finalize",240017)
2126 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2127 "hdf5_module_finalize",240007)
2129 if (i_stat .ne. 0)
call raise_exception(
"Unable to close temperature dataset.",&
2130 "hdf5_module_finalize",240007)
2132 if (i_stat .ne. 0)
call raise_exception(
"Unable to close density dataset.",&
2133 "hdf5_module_finalize",240007)
2135 if (i_stat .ne. 0)
call raise_exception(
"Unable to close entropy dataset.",&
2136 "hdf5_module_finalize",240007)
2138 if (i_stat .ne. 0)
call raise_exception(
"Unable to close radius dataset.",&
2139 "hdf5_module_finalize",240007)
2141 if (i_stat .ne. 0)
call raise_exception(
"Unable to close electron fraction dataset.",&
2142 "hdf5_module_finalize",240007)
2144 if (i_stat .ne. 0)
call raise_exception(
"Unable to close iteration dataset.",&
2145 "hdf5_module_finalize",240007)
2147 if (i_stat .ne. 0)
call raise_exception(
"Unable to close ylight dataset.",&
2148 "hdf5_module_finalize",240007)
2150 if (i_stat .ne. 0)
call raise_exception(
"Unable to close yheavy dataset.",&
2151 "hdf5_module_finalize",240007)
2153 if (i_stat .ne. 0)
call raise_exception(
"Unable to close yn dataset.",&
2154 "hdf5_module_finalize",240007)
2156 if (i_stat .ne. 0)
call raise_exception(
"Unable to close yp dataset.",&
2157 "hdf5_module_finalize",240007)
2159 if (i_stat .ne. 0)
call raise_exception(
"Unable to close ya dataset.",&
2160 "hdf5_module_finalize",240007)
2167 if (i_stat .ne. 0)
call raise_exception(
"Unable to close track_nuclei group.",&
2168 "hdf5_module_finalize",240017)
2170 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2171 "hdf5_module_finalize",240007)
2173 if (i_stat .ne. 0)
call raise_exception(
"Unable to close tracked abundance dataset.",&
2174 "hdf5_module_finalize",240007)
2181 if (i_stat .ne. 0)
call raise_exception(
"Unable to close nuloss group.",&
2182 "hdf5_module_finalize",240017)
2184 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2185 "hdf5_module_finalize",240007)
2187 if (i_stat .ne. 0)
call raise_exception(
"Unable to close beta dataset.",&
2188 "hdf5_module_finalize",240007)
2190 if (i_stat .ne. 0)
call raise_exception(
"Unable to close density dataset.",&
2191 "hdf5_module_finalize",240007)
2193 if (i_stat .ne. 0)
call raise_exception(
"Unable to close heat dataset.",&
2194 "hdf5_module_finalize",240007)
2196 if (i_stat .ne. 0)
call raise_exception(
"Unable to close nu total dataset.",&
2197 "hdf5_module_finalize",240007)
2199 if (i_stat .ne. 0)
call raise_exception(
"Unable to close radius dataset.",&
2200 "hdf5_module_finalize",240007)
2202 if (i_stat .ne. 0)
call raise_exception(
"Unable to close temperature dataset.",&
2203 "hdf5_module_finalize",240007)
2205 if (i_stat .ne. 0)
call raise_exception(
"Unable to close thermal dataset.",&
2206 "hdf5_module_finalize",240007)
2213 if (i_stat .ne. 0)
call raise_exception(
"Unable to close timescale group.",&
2214 "hdf5_module_finalize",240017)
2216 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2217 "hdf5_module_finalize",240007)
2220 if (i_stat .ne. 0)
call raise_exception(
"Unable to close timescale dataset.",&
2221 "hdf5_module_finalize",240007)
2229 if (i_stat .ne. 0)
call raise_exception(
"Unable to close energy group.",&
2230 "hdf5_module_finalize",240017)
2232 if (i_stat .ne. 0)
call raise_exception(
"Unable to close time dataset.",&
2233 "hdf5_module_finalize",240007)
2236 if (i_stat .ne. 0)
call raise_exception(
"Unable to close engen dataset.",&
2237 "hdf5_module_finalize",240007)
2241 if (i_stat .ne. 0)
call raise_exception(
"Unable to close decay energy dataset.",&
2242 "hdf5_module_finalize",240007)
2244 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (a,g) energy dataset.",&
2245 "hdf5_module_finalize",240007)
2247 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (p,g) energy dataset.",&
2248 "hdf5_module_finalize",240007)
2250 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (n,g) energy dataset.",&
2251 "hdf5_module_finalize",240007)
2253 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (p,n) energy dataset.",&
2254 "hdf5_module_finalize",240007)
2256 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (a,p) energy dataset.",&
2257 "hdf5_module_finalize",240007)
2259 if (i_stat .ne. 0)
call raise_exception(
"Unable to close (a,n) energy dataset.",&
2260 "hdf5_module_finalize",240007)
2262 if (i_stat .ne. 0)
call raise_exception(
"Unable to close other energy dataset.",&
2263 "hdf5_module_finalize",240007)
2265 if (i_stat .ne. 0)
call raise_exception(
"Unable to close fiss energy dataset.",&
2266 "hdf5_module_finalize",240007)
2268 if (i_stat .ne. 0)
call raise_exception(
"Unable to close tot energy dataset.",&
2269 "hdf5_module_finalize",240007)
2277 if (i_stat .ne. 0)
call raise_exception(
"Unable to close flow group.",&
2278 "hdf5_module_finalize",240017)
2282 call h5fclose_f(
file_id , i_stat)
2283 if (i_stat .ne. 0)
call raise_exception(
"Unable to close hdf5 file.",&
2284 "hdf5_module_finalize",240018)
2287 if (i_stat .ne. 0)
call raise_exception(
"Unable to close property list.",&
2288 "hdf5_module_finalize",240019)
character(len= *), parameter, private mainout_dsetname_temp
Name of the dataset.
integer h_track_nuclei_every
frequency of track nuclei output in hdf5 format
real(r_kind), dimension(chunk_size_mainout, 13), private chunk_store_mainout
Storage for the chunk which is later written to the file.
integer(hid_t), private nuloss_group_id
The ID of the group.
subroutine, private write_nu_loss_gain
Write the neutrino loss/gain into the hdf5 file.
integer, parameter, private chunk_size_ts
Chunk size.
integer(hid_t), private en_det_pn_id
integer(hid_t), private snaps_dset_id_t
Dataset identifier.
Module that contains necessary subroutines to write abundances to an hdf5 file.
subroutine, private write_engen(include_decay)
Write energy generations and time into the hdf5 file.
integer(hid_t), private nuloss_dset_id_t
Dataset identifier.
character(len= *), parameter, private mainout_dsetname_ye
Name of the dataset.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_bet
character *5, dimension(:), allocatable, public net_names
list of isotopes contained in the network
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_ng
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_pn
integer, public ipro
index of alphas, neutrons and protons
subroutine, private init_track_nuclei
Initialize the tracked_nuclei hdf5 group.
integer(hid_t), private en_det_bet_id
integer(hid_t), dimension(17), private ts_reac_ids
Dataset identifier.
character(len= *), parameter, private nuloss_dsetname_nut
Name of the dataset.
subroutine, private init_flow
Initialize the flow output.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_nut
character(len=14), private en_det_ap_name
character(len= *), parameter, private flow_dsetname_t
Name of the dataset.
logical, private hdf5_output
subroutine, public extend_snaps(time, Y_array)
Write abundances and time into the hdf5 file.
character(len= *), parameter, private mainout_dsetname_abar
Name of the dataset.
character(len= *), parameter, private mainout_dsetname_rad
Name of the dataset.
integer h_timescales_every
timescales output frequency in hdf5 format
subroutine, public extend_engen(time, engen_tot, heat, engen_ng_gn, engen_pg_gp, engen_ag_ga, engen_np_pn, engen_an_na, engen_ap_pa, engen_beta, engen_fiss, det_bet_engen, det_ag_engen, det_pg_engen, det_ng_engen, det_pn_engen, det_ap_engen, det_an_engen, det_other_engen, det_fiss_engen, det_tot_engen, include_decay)
Subroutine to write the energy generation per reaction type.
character(len= *), parameter, private mainout_dsetname_ylight
Name of the dataset.
integer, private chunk_counter_mainout
How full is the chunk already?
integer(hsize_t) function, private create_1d_dataset(dsetname, group_id)
Gets a dataset name and creates it.
integer, private chunk_counter_snaps
How full is the chunk already?
integer(hid_t), private nuloss_dset_id_the
Dataset identifier.
integer h_flow_every
flow output frequency in hdf5 format
integer(hid_t), private mainout_dset_id_yp
Dataset identifier.
character(len=14), private en_det_ag_name
integer, private iter_tracked
Iteration count, how often was already extended to the datasets?
character(len= *), parameter, private nuloss_dsetname_rad
Name of the dataset.
real(r_kind), dimension(:,:), allocatable, private chunk_store_track
Temporary storage for tracked nuclei.
integer(hid_t), private flow_group_id
The ID of the group.
integer(hsize_t) function, private create_1d_int_dataset(dsetname, group_id)
Gets a dataset name and creates an integer dataset.
Flow type to store flows.
integer(hid_t), private dset_id_y
Dataset identifier.
integer, parameter, private chunk_size_track
Chunk size.
integer(hid_t), private en_det_ng_id
character(len= *), parameter, private nuloss_dsetname_the
Name of the dataset.
character(len= *), parameter, private tracked_dsetname_t
Name of the dataset.
integer(hsize_t), dimension(1:2), private maxdims
Maximum dimensions (later set to unlimited)
integer, parameter, private chunk_size_nuloss
Chunk size.
character(:) function, allocatable, public int_to_str(num)
Converts a given integer to a string.
real(r_kind), dimension(:,:), allocatable, private chunk_store_snaps_y
Storage for the chunk which is later written to the file.
integer(hid_t), dimension(10), private en_reac_ids
Dataset identifier.
integer, private iter_snaps
Iteration count, how often was already extended to the datasets?
integer(hid_t), private nuloss_dset_id_heat
Dataset identifier.
integer(hid_t), private nuloss_dset_id_bet
Dataset identifier.
subroutine, public raise_exception(msg, sub, error_code)
Raise a exception with a given error message.
integer(hid_t), private en_det_tot_id
character(len= *), parameter, private ts_dsetname_t
Name of the dataset.
integer(hid_t), private en_det_an_id
subroutine, public init_hdf5_module()
Initialize the hdf5 module.
integer(hid_t), private mainout_group_id
The ID of the group.
integer, private iter_mainout
Iteration count, how often was already extended to the datasets?
character(len= *), parameter, private mainout_dsetname_t
Name of the dataset.
character(len= *), parameter, private mainout_dsetname_dens
Name of the dataset.
subroutine, private init_av_timescales
Initialize the timescales hdf5 group.
integer(hsize_t), dimension(1:2), private dims_track
Array dimensions.
character(len=14), private en_det_pn_name
integer, parameter, private chunk_size_en
Chunk size.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_ag
subroutine, private write_track_data()
Write the abundances of track nuclei into a file.
integer, private chunk_counter_ts
How full is the chunk already?
integer(hid_t), private en_dset_id_t
Dataset identifier.
character(len= *), parameter, private mainout_dsetname_ya
Name of the dataset.
logical h_finab
Store the finab in hdf5 format rather than in ascii format.
integer(hid_t), private en_det_ag_id
integer, dimension(:), allocatable, public track_nuclei_indices
Variables related to tracking individual nuclei.
integer(hsize_t), dimension(1:2), private dims_y
Array dimensions (net_size)
type(isotope_type), dimension(:), allocatable, public isotope
all nuclides used in the network
integer(hid_t), private mainout_dset_id_temp
Dataset identifier.
subroutine, private create_constant_1d_arrays(length, data, dsetname, group_id)
Create a hdf5 entry with constant values.
character(len= *), parameter, private nuloss_dsetname_bet
Name of the dataset.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_temp
character(len=14), private en_det_fiss_name
character(len=9), dimension(17), private ts_names
Names of the timescales.
integer(hid_t), private mainout_dset_id_rad
Dataset identifier.
integer, private iter_flow
Iteration count, how often was already extended to the datasets?
character(len= *), parameter, private dsetname_u
Name of the units dataset.
character(len=11), dimension(10), private en_names
Name of the energy generations.
subroutine, private write_snaps()
Write the abundances of all nuclei into a file.
subroutine, private init_nu_loss_gain
Initialize neutrino loss/gain output.
integer(hid_t), private mainout_dset_id_ye
Dataset identifier.
integer(hid_t), private mainout_dset_id_t
Dataset identifier.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_an
character(len= *), parameter, private nuloss_dsetname_dens
Name of the dataset.
integer(hid_t), private snaps_group_id
The ID of the group.
character(len= *), parameter, private mainout_dsetname_yheavy
Name of the dataset.
subroutine, public write_finab(Y)
Write final abundances into hdf5 file.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_dens
subroutine, public extend_flow(time, temp, dens, flow, n_flows, Y)
Extend the flow in the file With every call a new group is created and the flows are written in this ...
integer h_nu_loss_every
Output neutrino loss and gain in hdf5 format.
integer, private chunk_counter_en
How full is the chunk already?
character(len=14), private en_det_pg_name
logical h_engen_detailed
Output the energy per parent nucleus and reaction type.
character(len= *), parameter, private en_dsetname_t
Name of the dataset.
character(len=14), private en_det_ng_name
integer(hid_t), private en_det_ap_id
character(len=20), private hdf5_filename
Name of the hdf5 file.
integer(hid_t), private track_group_id
The ID of the group.
integer(hid_t), private nuloss_dset_id_temp
Dataset identifier.
character(len= *), parameter, private nuloss_dsetname_t
Name of the dataset.
character(len=14), private en_det_an_name
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_heat
real(r_kind), dimension(chunk_size_en, 11), private chunk_store_en
Temporary storage for energy generation.
subroutine, private extend_1d_int_dataset(dset_id, value, in_size, chunksize)
Extend a previously created 1D integer - dataset.
subroutine, public extend_nu_loss_gain(time, temp, dens, rad, the, heat, bet)
Write neutrino gain and loss.
integer, public net_size
total number of isotopes (network size)
integer, private iter_ts
Iteration count, how often was already extended to the datasets?
integer(hsize_t), dimension(1), private dims_1d
Helper variable to specify 1 dimension.
integer(hid_t), private mainout_dset_id_yheavy
Dataset identifier.
subroutine, private init_snaps
Initialize the snapshot hdf5 group.
integer(hid_t), private file_id
File identifier.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_t
character(len= *), parameter, private tracked_dsetname_y
Name of the dataset.
integer(hid_t), private crp_list
Dataset creation property identifier.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_tot
integer(hid_t), private nuloss_dset_id_rad
Dataset identifier.
integer(hid_t), private en_det_other_id
subroutine, private init_engen(include_decay)
Initialize energy generation output.
character(len= *), parameter, private nuloss_dsetname_temp
Name of the dataset.
integer(hid_t), private mainout_dset_id_dens
Dataset identifier.
subroutine, private extend_1d_dataset(dset_id, value, in_size, chunksize)
Extend a previously created 1D - dataset.
integer(hid_t), private dataspace
Dataspace identifier.
character(len=14), private en_det_other_name
integer, public track_nuclei_nr
amount of tracked nuclei
integer(hid_t), private mainout_dset_id_abar
Dataset identifier.
integer(hid_t), private en_det_fiss_id
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_ap
integer h_mainout_every
HDF5 output frequency of the mainout.
integer(hid_t), private en_group_id
The ID of the group.
integer(hid_t), private finab_group_id
The ID of the group.
integer, dimension(chunk_size_mainout), private chunk_store_int_mainout
Storage for integers in chunk which is later written to the file.
logical h_custom_snapshots
Same, but in hdf5 format.
subroutine, private write_av_timescales()
Write average timescales and time into the hdf5 file.
subroutine, private init_mainout
Initialize the mainout hdf5 group.
integer(hid_t), private mainout_dset_id_entr
Dataset identifier.
real(r_kind), dimension(chunk_size_ts, 18), private chunk_store_ts
Temporary storage for average timescales.
subroutine, private create_constant_1d_int_arrays(length, data, dsetname, group_id)
Create a hdf5 entry with constant values.
character(len= *), parameter, private mainout_dsetname_iter
Name of the dataset.
integer(hid_t), private en_det_pg_id
integer, parameter, private chunk_size_mainout
Chunk size.
integer, private iter_en
Iteration count, how often was already extended to the datasets?
integer, parameter, private chunk_size_snaps
Chunk size.
character(len= *), parameter, private mainout_dsetname_yn
Name of the dataset.
integer(hid_t), private mainout_dset_id_ya
Dataset identifier.
character(len= *), parameter, private dsetname_y
Name of the dataset.
integer(hid_t), private nuloss_dset_id_dens
Dataset identifier.
subroutine, private write_mainout_data()
Write the mainout data into a file.
character(len= *), parameter, private snaps_dsetname_t
Name of the dataset.
integer, private chunk_counter_track
How full is the chunk already?
character(len= *), parameter, private nuloss_dsetname_heat
Name of the dataset.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_fiss
integer, private chunk_counter_nuloss
How full is the chunk already?
subroutine, public extend_track_nuclei(time, Y_array)
Write tracked abundances and time into the hdf5 file.
character(len= *), parameter, private mainout_dsetname_yp
Name of the dataset.
real(r_kind), dimension(chunk_size_snaps), private chunk_store_snaps_time
Storage for the chunk which is later written to the file.
subroutine, public hdf5_module_finalize()
Finalize the hdf5 module.
subroutine, public extend_av_timescales(time, tau_ga, tau_ag, tau_ng, tau_gn, tau_pg, tau_gp, tau_np, tau_pn, tau_an, tau_na, tau_ap, tau_pa, tau_beta, tau_alpha, tau_nfiss, tau_sfiss, tau_bfiss)
Subroutine to write the average timescales (1/s) into the hdf5.
integer(hid_t), private mainout_dset_id_iter
Dataset identifier.
Contains types and objects shared between multiple modules.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_other
integer(hid_t), private ts_dset_id_t
Dataset identifier.
integer h_snapshot_every
snapshot output frequency in hdf5 format
integer, private iter_nuloss
Iteration count, how often was already extended to the datasets?
character(len= *), parameter, private mainout_dsetname_zbar
Name of the dataset.
character(len=5), dimension(:), allocatable, private tracked_names
Names of the tracked nuclei.
character(len=14), private en_det_bet_name
subroutine, public extend_mainout(cnt, time, temp, dens, entr, rad, Y)
Write the mainout data to the hdf5.
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_rad
integer(hid_t), private ts_group_id
The ID of the group.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_bet
character(len= *), parameter, private mainout_dsetname_entr
Name of the dataset.
integer(hid_t), private flow_dset_id_t
Dataset identifier.
subroutine, private write_units()
Write the units of the individual entries into the hdf5 file.
integer(hid_t), private mainout_dset_id_zbar
Dataset identifier.
integer(hid_t), private nuloss_dset_id_nut
Dataset identifier.
integer h_engen_every
Energy generation output frequency in hdf5 format.
integer(hid_t), private tracked_dset_id_t
Dataset identifier.
integer(hid_t), private tracked_abu_id
Dataset identifier.
integer(hid_t), private mainout_dset_id_yn
Dataset identifier.
character(len=14), private en_det_tot_name
Contains all runtime parameters as well as phys and math constants.
integer(hid_t), private mainout_dset_id_ylight
Dataset identifier.
subroutine, private extend_hdf5_y(dset_id, Y_array, nr_nuc, old_size, chunksize)
Extend the Y dataset for the net timestep.
real(r_kind), dimension(:,:), allocatable, private chunk_store_en_pg
real(r_kind), dimension(chunk_size_nuloss), private chunk_store_nuloss_the