22 integer :: fissnuc_index
31 integer,
dimension(:),
allocatable :: fissparts
32 integer,
dimension(:,:),
allocatable :: cscf_ind
34 real(
r_kind),
dimension(:),
allocatable :: q_value
35 real(
r_kind),
dimension(:),
allocatable :: channelprob
36 real(
r_kind),
dimension(:),
allocatable :: ch_amount
54 integer,
dimension(:),
allocatable :: net_idx
55 integer,
dimension(:),
allocatable :: z
56 integer,
dimension(:),
allocatable :: a
57 real(
r_kind),
dimension(:),
allocatable :: y
103 integer :: alloc_stat
117 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate" failed.',&
118 "init_fission_rates",190001)
151 integer,
dimension(:),
allocatable :: indices
153 info_entry(
"reorder_fragments")
160 if (
allocated(indices))
deallocate(indices)
161 allocate(indices(
fissrate(i)%dimens-2))
171 info_exit(
"reorder_fragments")
186 character(len=7) :: tmp
188 if (verbose_level .ge. 1)
then
191 elseif (verbose_level .ge. 2)
then
192 if (
nfiss .gt. 0)
write(*,
"(A)")
""
193 if (
nfiss .gt. 0)
write(*,
"(A)")
" Fission rates: "
194 if (
nfiss .gt. 0)
write(*,
"(A)")
" |----------------------|"
196 if (
nfiss .gt. 0)
write(*,
"(A)")
" | Total :"//adjustr(tmp)//
" |"
198 if (
n_nf .gt. 0)
write(*,
"(A)")
" | n-induced :"//adjustr(tmp)//
" |"
200 if (
n_bf .gt. 0)
write(*,
"(A)")
" | b-delayed :"//adjustr(tmp)//
" |"
202 if (
n_sf .gt. 0)
write(*,
"(A)")
" | spontaneous :"//adjustr(tmp)//
" |"
203 if (
nfiss .gt. 0)
write(*,
"(A)")
" |----------------------|"
204 if (
nfiss .gt. 0)
write(*,
"(A)")
""
224 character(len=*),
intent(in) :: path
227 integer :: alloc_stat
228 logical :: is_allocated
243 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate" failed.',&
244 "read_binary_fission_reaction_data",190001)
248 read(file_id)
fissrate(i)%fissnuc_index
257 read(file_id) is_allocated
258 if (is_allocated)
then
261 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate(i)%fissparts" failed.',&
262 "read_binary_fission_reaction_data",190001)
266 read(file_id) is_allocated
267 if (is_allocated)
then
269 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate(i)%cscf_ind" failed.',&
270 "read_binary_fission_reaction_data",190001)
276 read(file_id) is_allocated
277 if (is_allocated)
then
279 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate(i)%q_value" failed.',&
280 "read_binary_fission_reaction_data",190001)
284 read(file_id) is_allocated
285 if (is_allocated)
then
287 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate(i)%channelprob" failed.',&
288 "read_binary_fission_reaction_data",190001)
289 read(file_id)
fissrate(i)%channelprob
292 read(file_id) is_allocated
293 if (is_allocated)
then
295 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate(i)%ch_amount" failed.',&
296 "read_binary_fission_reaction_data",190001)
319 character(len=*),
intent(in) :: path
336 write(file_id)
fissrate(i)%fissnuc_index
343 write(file_id)
fissrate(i)%reac_type
344 write(file_id)
allocated(
fissrate(i)%fissparts)
345 if (
allocated(
fissrate(i)%fissparts))
then
346 write(file_id)
fissrate(i)%fissparts
348 write(file_id)
allocated(
fissrate(i)%cscf_ind)
349 if (
allocated(
fissrate(i)%cscf_ind))
then
353 write(file_id)
allocated(
fissrate(i)%q_value)
354 if (
allocated(
fissrate(i)%q_value))
then
357 write(file_id)
allocated(
fissrate(i)%channelprob)
358 if (
allocated(
fissrate(i)%channelprob))
then
359 write(file_id)
fissrate(i)%channelprob
361 write(file_id)
allocated(
fissrate(i)%ch_amount)
362 if (
allocated(
fissrate(i)%ch_amount))
then
363 write(file_id)
fissrate(i)%ch_amount
393 integer,
intent(inout) :: rrate_length
394 integer,
intent(out) :: fiss_count
395 integer :: new_length
400 new_length = rrate_length
401 if (
nfiss .ne. 0)
then
408 new_length = rrate_length
411 rrate_length = new_length
451 integer,
intent(inout) :: rrate_length
452 type(
fissionrate_type),
dimension(:),
allocatable,
intent(inout) :: fissrate_in
453 integer,
intent(inout) :: nfiss_in
455 real(
r_kind),
dimension(net_size) :: lambdas
457 real(
r_kind) :: tot_fiss_prob
458 integer :: parent_idx
461 logical,
dimension(rrate_length) :: rrate_mask
462 logical,
dimension(nfiss_in) :: fissrate_mask
465 integer :: new_rrate_length
466 integer :: new_nfiss_in
467 integer :: alloc_stat
468 real(
r_kind) :: lambda_rate
470 info_entry(
"modify_halflifes")
474 fissrate_mask(:) = .true.
477 if (
allocated(rrate))
then
478 rrate_mask(:) = .true.
484 if ((rrate(i)%reac_src .eq.
rrs_twr) .or. &
485 (rrate(i)%reac_src .eq.
rrs_nu) .or. &
486 (rrate(i)%reac_src .eq.
rrs_fiss) .or. &
487 (rrate(i)%reac_src .eq.
rrs_aext) .or. &
488 (rrate(i)%reac_src .eq.
rrs_detb)) cycle
490 if (rrate(i)%reac_type.eq.
rrt_betm)
then
492 if ((rrate(i)%parts(j).ne.
ineu) .and. (rrate(i)%parts(j).ne.
ipro))
then
497 if ((rrate(i)%reac_src .eq.
rrs_reacl) .or. &
498 (rrate(i)%reac_src .eq.
rrs_wext))
then
499 lambda_rate = dexp(rrate(i)%param(1))
500 lambdas(rrate(i)%parts(j)) = lambdas(rrate(i)%parts(j))+lambda_rate
501 if ((1d0-tot_fiss_prob) .ne. 0d0)
then
502 rrate(i)%param(1) = rrate(i)%param(1)+dlog(1d0-tot_fiss_prob)
505 rrate_mask(i) = .false.
508 else if (rrate(i)%reac_src .eq.
rrs_tabl)
then
514 lambdas(rrate(i)%parts(j)) = lambdas(rrate(i)%parts(j))+lambda_rate
515 if ((1d0-tot_fiss_prob) .ne. 0d0)
then
519 rrate_mask(i) = .false.
523 "modify_halflifes",190014)
533 new_rrate_length = count(rrate_mask)
534 if (new_rrate_length .ne. rrate_length)
then
537 if (verbose_level .ge. 2)
then
538 write(*,*)
"Reducing the rate array to the correct size by removing "//&
540 " beta-decays and putting it into fission rates."
544 allocate(rrate_tmp(new_rrate_length),stat=alloc_stat)
545 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "rrate_tmp" failed.',&
546 "modify_halflifes",190001)
548 rrate_tmp(1:new_rrate_length) = pack(rrate,rrate_mask)
549 deallocate(rrate,stat=alloc_stat)
550 if ( alloc_stat /= 0)
call raise_exception(
'Deallocation of "rrate" failed.',&
551 "modify_halflifes",190002)
553 allocate(rrate(new_rrate_length),stat=alloc_stat)
554 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "rrate" failed.',&
555 "modify_halflifes",190001)
557 rrate(1:new_rrate_length) = rrate_tmp(1:new_rrate_length)
558 deallocate(rrate_tmp,stat=alloc_stat)
559 if ( alloc_stat /= 0)
call raise_exception(
'Deallocation of "rrate_tmp" failed.',&
560 "modify_halflifes",190002)
561 rrate_length = new_rrate_length
566 do i=1,
size(fissrate_in)
567 if (fissrate_in(i)%reac_type.eq.
rrt_bf)
then
569 parent_idx = fissrate_in(i)%fissnuc_index
570 px_idx = int(fissrate_in(i)%released_n+1d0)
572 if (lambdas(parent_idx) .eq. 0d0)
then
574 fissrate_mask(i) = .false.
577 fissrate_in(i)%param(:) = 0d0
578 fissrate_in(i)%param(1) = dlog(px * lambdas(parent_idx))
584 new_nfiss_in = count(fissrate_mask)
585 if (new_nfiss_in .ne. nfiss_in)
then
587 if (verbose_level .ge. 2)
then
588 write(*,*)
"Reducing the fission rate array to the correct size by removing "//&
589 int_to_str(nfiss_in-new_nfiss_in)//
" beta-delayed fission rates."
592 allocate(fissrate_tmp(new_nfiss_in),stat=alloc_stat)
593 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp" failed.',&
594 "modify_halflifes",190001)
597 if (fissrate_mask(i))
then
599 allocate(fissrate_tmp(j)%fissparts(fissrate_in(i)%dimens),stat=alloc_stat)
600 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp(j)%fissparts" failed.',&
601 "modify_halflifes",190001)
602 allocate(fissrate_tmp(j)%cscf_ind(fissrate_in(i)%dimens,fissrate_in(i)%dimens),stat=alloc_stat)
603 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp(j)%cscf_ind" failed.',&
604 "modify_halflifes",190001)
605 allocate(fissrate_tmp(j)%q_value(fissrate_in(i)%channels),stat=alloc_stat)
606 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp(j)%q_value" failed.',&
607 "modify_halflifes",190001)
608 allocate(fissrate_tmp(j)%channelprob(fissrate_in(i)%channels),stat=alloc_stat)
609 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp(j)%channelprob" failed.',&
610 "modify_halflifes",190001)
611 allocate(fissrate_tmp(j)%ch_amount(fissrate_in(i)%dimens),stat=alloc_stat)
612 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_tmp(j)%ch_amount" failed.',&
613 "modify_halflifes",190001)
615 fissrate_tmp(j) = fissrate_in(i)
621 deallocate(fissrate_in,stat=alloc_stat)
622 if ( alloc_stat /= 0)
call raise_exception(
'Deallocation of "fissrate_in" failed.',&
623 "modify_halflifes",190002)
624 allocate(fissrate_in(new_nfiss_in),stat=alloc_stat)
625 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in" or failed.',&
626 "modify_halflifes",190001)
630 allocate(fissrate_in(i)%fissparts(fissrate_tmp(i)%dimens),stat=alloc_stat)
631 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in(i)%fissparts" failed.',&
632 "modify_halflifes",190001)
633 allocate(fissrate_in(i)%cscf_ind(fissrate_tmp(i)%dimens,fissrate_tmp(i)%dimens),stat=alloc_stat)
634 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in(i)%cscf_ind" failed.',&
635 "modify_halflifes",190001)
636 allocate(fissrate_in(i)%q_value(fissrate_tmp(i)%channels),stat=alloc_stat)
637 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in(i)%q_value" failed.',&
638 "modify_halflifes",190001)
639 allocate(fissrate_in(i)%channelprob(fissrate_tmp(i)%channels),stat=alloc_stat)
640 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in(i)%channelprob" failed.',&
641 "modify_halflifes",190001)
642 allocate(fissrate_in(i)%ch_amount(fissrate_tmp(i)%dimens),stat=alloc_stat)
643 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "fissrate_in(i)%ch_amount" failed.',&
644 "modify_halflifes",190001)
645 fissrate_in(i) = fissrate_tmp(i)
650 deallocate(fissrate_tmp,stat=alloc_stat)
651 if ( alloc_stat /= 0)
call raise_exception(
'Deallocation of "fissrate_tmp" or failed.',&
652 "modify_halflifes",190002)
653 nfiss_in = new_nfiss_in
659 if ( alloc_stat /= 0)
call raise_exception(
'Deallocation of "beta_delayed_fiss_probs" failed.',&
660 "modify_halflifes",190002)
662 info_exit(
"modify_halflifes")
699 integer :: alloc_stat
701 info_entry(
"count_fission_rates")
712 call raise_exception(
"Fission format for spontaneous fission rates not implemented, got '"//&
714 "'fission_format_spontaneous' parameter.",&
715 "count_fission_rates",190010)
725 call raise_exception(
"Fission format for neutron-induced fission rates not implemented, got '"//&
727 "'fission_format_n_induced' parameter.",&
728 "count_fission_rates",190010)
743 if ( alloc_stat /= 0)
call raise_exception(
'Allocation of "beta_delayed_fiss_probs" failed.',&
744 "count_fission_rates",190001)
747 call raise_exception(
"Fission format for beta-delayed fission rates not implemented, got '"//&
749 "'fission_format_beta_delayed' parameter.",&
750 "count_fission_rates",190010)
756 info_exit(
"count_fission_rates")
775 character(len=*),
intent(in) :: fission_rate_file
776 integer,
intent(out) :: count_rates
779 integer :: fissionlib
780 character(5),
dimension(6) :: parts
781 integer,
dimension(6) :: parts_index
782 real(
r_kind),
dimension(9) :: params
798 read(fissionlib,
my_format(1), iostat = read_stat) &
799 grp, parts(1:6), src, res, rev, q
800 if (read_stat /= 0)
exit
801 read(fissionlib,
"(4E13.6)") params(1:4)
802 read(fissionlib,
"(5E13.6)") params(5:9)
803 if (grp.ne.0) cycle fission
806 if (parts(j) .eq.
' ')
exit indices
807 parts_index(j) =
benam(parts(j))
809 if (parts_index(j) .eq. 0) cycle fission
844 character(len=*),
intent(in) :: fission_rate_file
845 integer,
intent(out) :: count_rates
848 integer :: fissionlib
849 character(5),
dimension(6) :: parts
850 integer,
dimension(6) :: parts_index
852 real(
r_kind),
dimension(9) :: params
868 read(fissionlib, *, iostat = read_stat) &
869 parts(1), params(1), src
870 if (read_stat /= 0)
exit
874 if (parts(j) .eq.
' ')
exit indices_hl
875 parts_index(j) =
benam(parts(j))
877 if (parts_index(j) .eq. 0) cycle fission_hl
928 character(len=*),
intent(in) :: fission_rate_file
929 integer,
intent(out) :: count_rates
930 integer,
intent(out) :: nr_cols
931 integer,
parameter :: maxcols = 30
934 integer :: fissionlib
935 character(5) :: parent
936 integer :: idx_parent
938 real(
r_kind),
dimension(maxcols):: buff
939 character(500) :: dummy
941 integer :: fiss_probs
955 read(fissionlib, * , iostat = read_stat) parent, src
957 if (read_stat /= 0)
exit
960 if (ifiss .eq. 0)
then
961 read(fissionlib,
'(A)', iostat = read_stat) dummy
963 read(dummy,*,iostat=read_stat) buff(1:j)
964 if (read_stat .ne. 0)
exit
967 nr_cols = count(buff .ne. -1d0)
969 read(fissionlib,*, iostat = read_stat) buff(1:nr_cols)
970 if (read_stat .ne. 0)
call raise_exception(
'Reading the fission probabilities failed.',&
971 "count_fission_rates_probability_format",190012)
975 fiss_probs = count((buff .ne. -1d0) .and. (buff .ne. 0d0))
978 idx_parent =
benam(parent)
980 if (idx_parent .eq. 0) cycle fission_prob
984 ifiss = ifiss + fiss_probs
1017 integer :: neutfission=0
1018 integer :: betafission=0
1019 integer :: fissindex
1022 integer :: nfrac_distr
1031 if (astat /= 0)
call raise_exception(
'Allocation of "fissfrags_bdel" failed.',&
1032 "add_fission_fragments",190001)
1054 do fissindex=1,
nfiss
1074 //
") not implemented yet. "//&
1075 "Set it to a supported value!",&
1076 "add_fission_fragments",&
1084 //
") not implemented yet. "//&
1085 "Set it to a supported value!",&
1086 "add_fission_fragments",&
1107 //
") not implemented yet. "//&
1108 "Set it to a supported value!",&
1109 "add_fission_fragments",&
1117 //
") not implemented yet. "//&
1118 "Set it to a supported value!",&
1119 "add_fission_fragments",&
1140 //
") not implemented yet. "//&
1141 "Set it to a supported value!",&
1142 "add_fission_fragments",&
1150 //
") not implemented yet. "//&
1151 "Set it to a supported value!",&
1152 "add_fission_fragments",&
1167 if (astat.ne.0)
call raise_exception(
"Could not deallocate fissfrags_beta_delayed array.",&
1168 "add_fission_fragments",190002)
1172 if (astat.ne.0)
call raise_exception(
"Could not deallocate fissfrags_n_induced array.",&
1173 "add_fission_fragments",190002)
1177 if (astat.ne.0)
call raise_exception(
"Could not deallocate fissfrags_spontaneous array.",&
1178 "add_fission_fragments",190002)
1213 integer :: fisscount
1230 call raise_exception(
"Fission format for spontaneous fission rates not implemented, got '"//&
1232 "'fission_format_spontaneous' parameter.",&
1233 "read_fission_rates",190010)
1244 call raise_exception(
"Fission format for neutron-induced fission rates not implemented, got '"//&
1246 "'fission_format_n_induced' parameter.",&
1247 "read_fission_rates",190010)
1265 call raise_exception(
"Fission format for beta-delayed fission rates not implemented, got '"//&
1267 "'fission_format_beta_delayed' parameter.",&
1268 "read_fission_rates",190010)
1307 character(len=*),
intent(in) :: fission_path
1308 integer,
intent(in) :: reac_type
1309 integer,
intent(inout) :: start_idx
1310 integer :: read_stat
1311 integer :: fissionlib
1312 character(5),
dimension(6) :: parts
1313 integer,
dimension(6) :: parts_index
1314 real(
r_kind),
dimension(9) :: params
1320 integer :: group_index
1321 integer :: fissindex
1330 fissindex = start_idx
1332 read(fissionlib,
my_format(1), iostat = read_stat) &
1333 grp, parts(1:6), src, res, rev, q
1334 if (read_stat /= 0)
exit fission_loop
1335 read(fissionlib,
"(4E13.6)") params(1:4)
1336 read(fissionlib,
"(5E13.6)") params(5:9)
1344 inner_fission:
do j=1,6
1345 if (parts(j) .eq.
' ')
exit inner_fission
1346 parts_index(j) =
benam(parts(j))
1348 if (parts_index(j) .eq. 0) cycle fission_loop
1349 end do inner_fission
1360 if (reac_type.eq.
rrt_nf)
then
1361 fissrate(fissindex)%reac_type = reac_type
1364 if (parts_index(1) .eq.
ineu)
then
1365 fissrate(fissindex)%fissnuc_index = parts_index(2)
1367 fissrate(fissindex)%fissnuc_index = parts_index(1)
1370 else if (reac_type.eq.
rrt_sf)
then
1371 fissrate(fissindex)%reac_type = reac_type
1374 fissrate(fissindex)%fissnuc_index = parts_index(1)
1376 else if (reac_type.eq.
rrt_bf)
then
1378 fissrate(fissindex)%released_n = count(parts_index(:) .eq.
ineu)
1379 fissrate(fissindex)%reac_type = reac_type
1382 fissrate(fissindex)%fissnuc_index = parts_index(1)
1389 fissindex = fissindex + 1
1393 start_idx = fissindex
1440 character(len=*),
intent(in) :: fission_path
1441 integer,
intent(in) :: reac_type
1442 integer,
intent(inout) :: start_idx
1444 integer :: read_stat
1445 integer :: fissionlib
1446 character(5),
dimension(6) :: parts
1447 integer,
dimension(6) :: parts_index
1448 real(
r_kind),
dimension(9) :: params
1450 integer :: group_index
1451 integer :: fissindex
1453 integer :: idx_parent
1454 character(5) :: parent
1455 real(
r_kind),
dimension(amount_cols):: pnf
1458 fissindex = start_idx
1467 read(fissionlib,*, iostat = read_stat) parent, src
1468 if (read_stat /= 0)
exit
1469 read(fissionlib,*, iostat = read_stat) pnf
1471 if (read_stat /= 0)
call raise_exception(
"Could not read fission probabilities.",&
1472 "read_fission_rates_probability_format",&
1477 idx_parent =
benam(parent)
1478 if (idx_parent .eq. 0) cycle fission_prob
1482 parts_index(1) = idx_parent
1485 if (pnf(j) .eq. 0d0) cycle
1491 if (reac_type.eq.
rrt_bf)
then
1493 fissrate(fissindex)%reac_type = reac_type
1496 fissrate(fissindex)%fissnuc_index = parts_index(1)
1497 fissrate(fissindex)%released_n = j-1
1499 call raise_exception(
"Fission type not implemented for probability format!",&
1500 "read_fission_rates_probability_format",&
1513 fissindex = fissindex + 1
1518 " is larger than one.",
"read_fission_rates_probability_format",&
1553 character(len=*),
intent(in) :: fission_path
1554 integer,
intent(in) :: reac_type
1555 integer,
intent(inout) :: start_idx
1556 integer :: read_stat
1557 integer :: fissionlib
1558 character(5),
dimension(6) :: parts
1559 integer,
dimension(6) :: parts_index
1560 real(
r_kind),
dimension(9) :: params
1562 integer :: group_index
1563 integer :: fissindex
1572 fissindex = start_idx
1577 read(fissionlib, *, iostat = read_stat) &
1578 parts(1), params(1), src
1579 if (read_stat /= 0)
exit fission_loop
1581 params(1) = dlog(dlog(2d0)/params(1))
1583 inner_fission:
do j=1,6
1584 if (parts(j) .eq.
' ')
exit inner_fission
1585 parts_index(j) =
benam(parts(j))
1587 if (parts_index(j) .eq. 0) cycle fission_loop
1588 end do inner_fission
1597 if (reac_type .eq.
rrt_sf)
then
1598 fissrate(fissindex)%reac_type = reac_type
1601 fissrate(fissindex)%fissnuc_index = parts_index(1)
1603 else if (reac_type .eq.
rrt_bf)
then
1604 fissrate(fissindex)%reac_type = reac_type
1607 fissrate(fissindex)%fissnuc_index = parts_index(1)
1609 call raise_exception(
"Fission type not implemented for half life format!",&
1610 "read_fission_rates_halflife_format",&
1621 fissindex = fissindex + 1
1625 start_idx = fissindex
1657 integer :: fiss_mode
1659 integer :: alloc_stat
1661 info_entry(
"fiss_dist")
1668 mass =
isotope(fissrate_inout%fissnuc_index)%mass
1669 pnr =
isotope(fissrate_inout%fissnuc_index)%p_nr
1672 if (fissrate_inout%reac_type .eq.
rrt_nf)
then
1677 else if (fissrate_inout%reac_type .eq.
rrt_sf)
then
1682 else if (fissrate_inout%reac_type .eq.
rrt_bf)
then
1687 nemiss = fissrate_inout%released_n
1700 z2 = nint(52.01 - (zf - 80)/1.d1)
1707 if (a1.gt.
minmax(z1,2))
then
1708 nemiss = a1 -
minmax(z1,2)
1710 else if (a1.lt.
minmax(z1,1))
then
1711 if (verbose_level .ge. 2)
then
1712 print *,
'a1 error:', z1, a1, zf, af
1713 print *,
'adjusting...'
1719 if (a2.gt.
minmax(z2,2))
then
1720 nemiss = nemiss + a2 -
minmax(z2,2)
1723 else if (a2.lt.
minmax(z2,1))
then
1724 if (verbose_level .ge. 2)
then
1725 print *,
'a2 error:', z2, a2, zf, af
1726 print *,
'adjusting...'
1733 fissrate_inout%fissnuc_index =
findaz(mass,pnr)
1734 fissrate_inout%channels = 1
1736 allocate(fissrate_inout%channelprob(1),fissrate_inout%q_value(1),stat=alloc_stat)
1737 if (alloc_stat.ne.0)
call raise_exception(
"Could not allocate memory for channelprob and q_value.",&
1740 fissrate_inout%channelprob(1) = 1.d0
1741 fissrate_inout%mode = fiss_mode
1743 if (fiss_mode.eq.1)
then
1744 fissrate_inout%dimens = 4
1745 allocate(fissrate_inout%fissparts(fissrate_inout%dimens), &
1746 fissrate_inout%ch_amount(fissrate_inout%dimens),stat=alloc_stat)
1747 if (alloc_stat.ne.0)
call raise_exception(
"Could not allocate 'fissparts' and 'ch_amount'.",&
1749 fissrate_inout%fissparts(1) =
ineu
1750 fissrate_inout%ch_amount(1) = float(nemiss) - 1.d0
1751 fissrate_inout%fissparts(2) = fissrate_inout%fissnuc_index
1752 fissrate_inout%ch_amount(2) = -1.d0
1753 fissrate_inout%fissparts(3) =
findaz(a1,z1)
1754 fissrate_inout%ch_amount(3) = 1.d0
1755 fissrate_inout%fissparts(4) =
findaz(a2,z2)
1756 fissrate_inout%ch_amount(4) = 1.d0
1757 fissrate_inout%q_value(1) = (1-nemiss) *
isotope(
ineu)%mass_exc + &
1758 isotope(fissrate_inout%fissnuc_index)%mass_exc - &
1759 isotope(fissrate_inout%fissparts(3))%mass_exc - &
1760 isotope(fissrate_inout%fissparts(4))%mass_exc
1762 if (nemiss .eq. 0)
then
1763 fissrate_inout%dimens = 3
1764 allocate(fissrate_inout%fissparts(fissrate_inout%dimens),&
1765 fissrate_inout%ch_amount(fissrate_inout%dimens),stat=alloc_stat)
1766 if (alloc_stat.ne.0)
call raise_exception(
"Could not allocate 'fissparts' and 'ch_amount'.",&
1769 fissrate_inout%dimens = 4
1770 allocate(fissrate_inout%fissparts(fissrate_inout%dimens),&
1771 fissrate_inout%ch_amount(fissrate_inout%dimens),stat=alloc_stat)
1772 if (alloc_stat.ne.0)
call raise_exception(
"Could not allocate 'fissparts' and 'ch_amount'.",&
1774 fissrate_inout%fissparts(4) =
ineu
1775 fissrate_inout%ch_amount(4) = float(nemiss)
1777 fissrate_inout%fissparts(1) = fissrate_inout%fissnuc_index
1778 fissrate_inout%ch_amount(1) = -1.d0
1779 fissrate_inout%fissparts(2) =
findaz(a1,z1)
1780 fissrate_inout%ch_amount(2) = 1.d0
1781 fissrate_inout%fissparts(3) =
findaz(a2,z2)
1782 fissrate_inout%ch_amount(3) = 1.d0
1784 if (fissrate_inout%fissparts(2) .eq. -1)
then
1785 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'? "//&
1789 if (fissrate_inout%fissparts(3) .eq. -1)
then
1790 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'? "//&
1795 fissrate_inout%q_value(1) =
isotope(fissrate_inout%fissparts(1))%mass_exc - &
1796 isotope(fissrate_inout%fissparts(2))%mass_exc - &
1797 isotope(fissrate_inout%fissparts(3))%mass_exc - &
1801 qval = fissrate_inout%q_value(1)
1802 fissrate_inout%averageQ = qval
1803 info_exit(
"fiss_dist")
1837 integer :: fiss_mode
1840 integer :: a1,a2,z1,z2
1842 real(
r_kind) :: za,al,ah,cz,ca
1843 real(
r_kind),
parameter :: lim = 1.d-6
1845 integer :: alloc_stat
1848 integer :: neutronflag
1849 real(
r_kind) :: nemiss_total
1850 integer :: additional_neutrons
1853 info_entry(
"kodtakdist")
1857 additional_neutrons = 0d0
1860 mass =
isotope(fissrate_inout%fissnuc_index)%mass
1861 pnr =
isotope(fissrate_inout%fissnuc_index)%p_nr
1864 if (fissrate_inout%reac_type .eq.
rrt_nf)
then
1869 else if (fissrate_inout%reac_type .eq.
rrt_sf)
then
1874 else if (fissrate_inout%reac_type .eq.
rrt_bf)
then
1879 additional_neutrons = fissrate_inout%released_n
1880 af = mass-additional_neutrons
1884 al = 0.85d0*af - 104.98d0
1885 ah = 0.15d0*af + 103.87d0
1893 za = zf*(a+0.6d0)/af
1895 paz = 0.5d0*dexp(-((z-za)**2)/cz) * &
1896 ((dexp(-((a-al)**2)/ca))+(dexp(-((a-ah)**2)/ca)))/ &
1897 (
unit%pi*dsqrt(cz*ca))
1907 if (a1.gt.
minmax(z1,2))
then
1909 else if (a1.lt.
minmax(z1,1))
then
1910 if (verbose_level .ge. 2)
then
1911 print *,
'Warning in kodtakdist, fragment was lighter than included nuclei'
1912 print *,
'Fragment was: ', a1, z1
1916 if (a2.gt.
minmax(z2,2))
then
1918 else if (a2.lt.
minmax(z2,1))
then
1919 if (verbose_level .ge. 2)
then
1920 print *,
'Warning in kodtakdist, fragment was lighter than included nuclei'
1921 print *,
'Fragment was: ', a2, z2
1929 fissrate_inout%channels = nf
1930 fissrate_inout%fissnuc_index =
findaz(mass,pnr)
1931 fissrate_inout%mode = fiss_mode
1934 allocate(fissrate_inout%channelprob(nf),fissrate_inout%q_value(nf), &
1936 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of fission rate prob. arrays failed.",&
1937 "kodtakdist",190001)
1940 select case (fiss_mode)
1942 dimens = 2 * fissrate_inout%channels + 2
1943 fissrate_inout%dimens = dimens
1944 allocate(fissrate_inout%fissparts(dimens),fissrate_inout%ch_amount(dimens),&
1946 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of fission rate parts failed.",&
1947 "kodtakdist",190001)
1949 fissrate_inout%fissparts(1) =
ineu
1950 fissrate_inout%fissparts(2) = fissrate_inout%fissnuc_index
1951 fissrate_inout%ch_amount(2) = -1.d0
1953 dimens = 2 * fissrate_inout%channels + 1
1954 if (neutronflag.eq.1) dimens = dimens+1
1955 fissrate_inout%dimens = dimens
1956 allocate(fissrate_inout%fissparts(dimens),fissrate_inout%ch_amount(dimens),&
1958 if (alloc_stat .ne. 0)
call raise_exception(
"Allocation of fission rate parts failed.",&
1959 "kodtakdist",190001)
1960 fissrate_inout%ch_amount(1) = -1.d0
1961 fissrate_inout%fissparts(1) = fissrate_inout%fissnuc_index
1965 nemiss_total = float(additional_neutrons)
1969 massloop:
do a = 1,af
1970 za = zf*(a+0.6d0)/af
1971 paz = 0.5d0*dexp(-((z-za)**2)/cz) * &
1972 ((dexp(-((a-al)**2)/ca))+(dexp(-((a-ah)**2)/ca)))/ &
1973 (
unit%pi*dsqrt(cz*ca))
1974 if(paz.lt.lim) cycle massloop
1984 if (a1.gt.
minmax(z1,2))
then
1985 nemiss = a1 -
minmax(z1,2)
1988 else if (a1.lt.
minmax(z1,1))
then
1989 if (verbose_level .ge. 2)
then
1990 print *,
'Warning in kodtakdist, fragment was lighter than included nuclei'
1991 print *,
'Fragment was: ', a1, z1
1995 if (a2.gt.
minmax(z2,2))
then
1996 nemiss = nemiss + a2 -
minmax(z2,2)
1998 else if (a2.lt.
minmax(z2,1))
then
2000 if (verbose_level .ge. 2)
then
2001 print *,
'Warning in kodtakdist, fragment was lighter than included nuclei'
2002 print *,
'Fragment was: ', a2, z2
2007 nemiss_total = nemiss_total+ float(nemiss)*paz
2009 select case(fiss_mode)
2011 fissrate_inout%channelprob(i) = paz
2012 fissrate_inout%fissparts(i+2) =
findaz(a1,z1)
2014 if (fissrate_inout%fissparts(i+2) .eq. -1)
then
2015 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'? "//&
2017 "kodtakdist",190013)
2020 fissrate_inout%ch_amount(i+2) = fissrate_inout%channelprob(i)
2021 fissrate_inout%fissparts(i+2+fissrate_inout%channels) =
findaz(a2,z2)
2023 if (fissrate_inout%fissparts(i+2+fissrate_inout%channels) .eq. -1)
then
2024 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'? "//&
2026 "kodtakdist",190013)
2028 fissrate_inout%ch_amount(i+2+fissrate_inout%channels) = fissrate_inout%channelprob(i)
2029 fissrate_inout%q_value(i) =
isotope(fissrate_inout%fissnuc_index)%mass_exc - &
2031 isotope(fissrate_inout%fissparts(i+2))%mass_exc - &
2032 isotope(fissrate_inout%fissparts(i+2+fissrate_inout%channels))%mass_exc
2034 fissrate_inout%channelprob(i) = paz
2035 fissrate_inout%fissparts(i+1) =
findaz(a1,z1)
2037 if (fissrate_inout%fissparts(i+1) .eq. -1)
then
2038 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'?",&
2039 "kodtakdist",190013)
2042 fissrate_inout%ch_amount(i+1) = fissrate_inout%channelprob(i)
2043 fissrate_inout%fissparts(i+1+fissrate_inout%channels) =
findaz(a2,z2)
2045 if (fissrate_inout%fissparts(i+1+fissrate_inout%channels) .eq. -1)
then
2046 call raise_exception(
"Fragment was not found in nucleus library, does your library have 'holes'?",&
2047 "kodtakdist",190013)
2050 fissrate_inout%ch_amount(i+1+fissrate_inout%channels) = fissrate_inout%channelprob(i)
2051 if (neutronflag.eq.1)
then
2052 fissrate_inout%fissparts(dimens) =
ineu
2053 fissrate_inout%ch_amount(dimens) = real(nemiss)
2055 fissrate_inout%q_value(i) =
isotope(fissrate_inout%fissnuc_index)%mass_exc - &
2057 isotope(fissrate_inout%fissparts(i+1))%mass_exc - &
2058 isotope(fissrate_inout%fissparts(i+1+fissrate_inout%channels))%mass_exc
2064 select case(fiss_mode)
2066 fissrate_inout%ch_amount(1) = nemiss_total - 1.d0
2068 if (neutronflag.eq.1)
then
2069 fissrate_inout%fissparts(dimens) =
ineu
2070 fissrate_inout%ch_amount(dimens) = nemiss_total
2075 do i=1,fissrate_inout%dimens
2076 qval = qval -
isotope(fissrate_inout%fissparts(i))%mass_exc*fissrate_inout%ch_amount(i)
2078 fissrate_inout%averageQ = qval
2082 info_exit(
"kodtakdist")
2125 character(len=*),
intent(in) :: file
2126 type(
fragtype),
dimension(:),
allocatable,
intent(out) :: fragment_array
2129 integer :: frag_counter
2134 integer :: read_stat
2138 real(r_kind) :: Ynem
2139 real(r_kind) :: helper
2140 integer :: mina,maxa
2143 info_entry(
"read_fiss_fragment_file")
2152 read(file_id,*,iostat=read_stat) zf, af, nof
2153 if (read_stat .ne. 0)
exit
2160 if (((
findaz(af,zf) .ne. -1) .or. &
2161 (
findaz(af-1,zf) .ne. -1) .or. &
2162 (
findaz(af,zf-1) .ne. -1)))
then
2163 frag_counter = frag_counter + 1
2168 allocate(fragment_array(frag_counter),stat=astat)
2169 if (astat .ne. 0)
then
2171 "read_fiss_fragment_file",190001)
2177 fragloop:
do i = 1,frag_counter
2178 read(file_id,*,iostat=read_stat) zf, af, nof
2185 if (.not. ((
findaz(af,zf) .ne. -1) .or. &
2186 (
findaz(af-1,zf) .ne. -1) .or. &
2187 (
findaz(af,zf-1) .ne. -1)))
then
2196 allocate(fragment_array(i)%net_idx(nof),&
2197 fragment_array(i)%Z(nof),&
2198 fragment_array(i)%A(nof),&
2199 fragment_array(i)%Y(nof),&
2201 if (astat .ne. 0)
then
2203 "read_fiss_fragment_file",190001)
2207 fragment_array(i)%nr_frags = nof
2208 fragment_array(i)%Zp = zf
2209 fragment_array(i)%Ap = af
2217 read(file_id,*) fragment_array(i)%Z(j),fragment_array(i)%A(j),fragment_array(i)%Y(j)
2218 fragment_array(i)%net_idx(j) =
findaz(fragment_array(i)%A(j),fragment_array(i)%Z(j))
2222 if (fragment_array(i)%net_idx(j) .eq. -1)
then
2224 if (verbose_level .ge. 2)
then
2225 write(*,*)
"Fragment ",fragment_array(i)%Z(j),fragment_array(i)%A(j),&
2226 " not in network. Splitting into lighter isotope and neutrons."
2229 mina =
minmax(fragment_array(i)%Z(j),1)
2230 maxa =
minmax(fragment_array(i)%Z(j),2)
2231 if (fragment_array(i)%A(j) .gt. maxa)
then
2233 nem = nem + (fragment_array(i)%A(j)-maxa)
2237 ynem = ynem + fragment_array(i)%Y(j)*(fragment_array(i)%A(j)-maxa)
2238 fragment_array(i)%Y(j) = fragment_array(i)%Y(j)
2239 fragment_array(i)%A(j) = maxa
2240 fragment_array(i)%net_idx(j) =
findaz(maxa,fragment_array(i)%Z(j))
2248 if (fragment_array(i)%net_idx(j) .eq. -1)
then
2251 if (verbose_level .ge. 2)
then
2252 write(*,*)
"Fission fragment Z=",fragment_array(i)%Z(j), &
2253 "A=",fragment_array(i)%A(j), &
2254 "Y=",fragment_array(i)%Y(j), &
2255 "not in network, trying to remap."
2259 mina =
minmax(fragment_array(i)%Z(j),1)
2260 bn = (mina-fragment_array(i)%A(j))
2265 if ((fragment_array(i)%net_idx(m) .ne. -1) .and. &
2266 ((bn*fragment_array(i)%Y(j) / fragment_array(i)%A(m)) .lt. fragment_array(i)%Y(m)))
then
2274 "consider to include more nuclides in the network. "//&
2276 "Nucleus Z="//trim(int_to_str(fragment_array(i)%Z(j)))// &
2277 ", A="//trim(int_to_str(fragment_array(i)%A(j)))// &
2278 " was not included. "//new_line(
'A')//&
2279 "Need "//int_to_str(bn)//
" neutrons.",&
2280 "read_fiss_fragment_file",&
2285 fragment_array(i)%Y(m) = fragment_array(i)%Y(m) - (bn*fragment_array(i)%Y(j) / fragment_array(i)%A(m))
2286 fragment_array(i)%Y(j) = fragment_array(i)%Y(j)
2288 fragment_array(i)%A(j) = mina
2289 fragment_array(i)%net_idx(j) =
findaz(mina,fragment_array(i)%Z(j))
2296 helper = helper+fragment_array(i)%A(j)*fragment_array(i)%Y(j)
2298 if (abs(helper - float(fragment_array(i)%Ap)) .gt. 1e-2)
then
2299 call raise_exception(
"Yield of fission fragments is not conserved. "// &
2300 "This happened for fragment distribution of "// &
2301 "Nucleus Z="//trim(int_to_str(fragment_array(i)%Zp))// &
2302 ", A="//trim(int_to_str(fragment_array(i)%Ap))// &
2303 ". Ensure that the sum is 2.",
"read_fiss_fragment_file",&
2307 fragment_array(i)%neutrons = nem
2308 fragment_array(i)%Yn = ynem
2315 info_exit(
"read_fiss_fragment_file")
2339 integer,
optional,
intent(in) :: missing_frags
2343 real(r_kind) :: qval
2344 integer :: nfrac_distr
2347 integer :: ind_parent
2353 integer :: additional_neutrons
2354 integer :: missing_frags_internal
2356 info_entry(
"file_fiss_frag")
2359 if (.not.
present(missing_frags))
then
2360 missing_frags_internal = 0
2362 missing_frags_internal = missing_frags
2366 mass =
isotope(fissrate_inout%fissnuc_index)%mass
2367 pnr =
isotope(fissrate_inout%fissnuc_index)%p_nr
2370 if (fissrate_inout%reac_type .eq.
rrt_bf)
then
2372 else if (fissrate_inout%reac_type .eq.
rrt_nf)
then
2374 else if (fissrate_inout%reac_type .eq.
rrt_sf)
then
2379 additional_neutrons = fissrate_inout%released_n
2383 afiss = mass - additional_neutrons
2386 if (fissrate_inout%reac_type .eq.
rrt_bf)
then
2389 elseif (fissrate_inout%reac_type .eq.
rrt_nf)
then
2395 ind_parent =
findaz(mass,pnr)
2396 fissrate_inout%fissnuc_index = ind_parent
2400 if (fissrate_inout%reac_type .eq.
rrt_bf)
then
2402 else if (fissrate_inout%reac_type .eq.
rrt_nf)
then
2404 else if (fissrate_inout%reac_type .eq.
rrt_sf)
then
2408 if ((fissfrags_tmp%Zp .eq. zfiss) .and. (fissfrags_tmp%Ap .eq. afiss))
then
2415 if (.not. found)
then
2418 if (missing_frags_internal .eq. 0)
then
2420 "A="//
int_to_str(afiss)//
" in fragment file. Either add it to the fission fragment file or set "//&
2421 "'fission_frag_missing' to 1 (Panov et al. 2001) or 2 (Kodama & Takahashi 1975)",&
2422 "file_fiss_frag", 190016)
2426 if (verbose_level .ge.2)
then
2427 write(*,*)
"No fragment distribution found for Z=",zfiss,
", A=",afiss
2428 if (missing_frags_internal .eq. 1)
then
2429 write(*,*)
"Using Panov 2001 instead."
2430 else if (missing_frags_internal .eq. 2)
then
2431 write(*,*)
"Using Kodama & Takahashi 1975 instead."
2435 if (missing_frags_internal .eq. 1)
then
2438 else if (missing_frags_internal .eq. 2)
then
2443 "Got "//
int_to_str(missing_frags_internal)//
".",&
2444 "file_fiss_frag",190015)
2448 if (fissrate_inout%reac_type .eq.
rrt_bf)
then
2450 else if (fissrate_inout%reac_type .eq.
rrt_nf)
then
2452 else if (fissrate_inout%reac_type .eq.
rrt_sf)
then
2457 allocate(fissrate_inout%fissparts(fissfrags_tmp%nr_frags+2),&
2458 fissrate_inout%ch_amount(fissfrags_tmp%nr_frags+2),&
2461 if (astat .ne. 0)
then
2462 call raise_exception(
"Could not allocate memory for fission fragments.",&
2463 "file_fiss_frag",190001)
2468 fissrate_inout%dimens = fissfrags_tmp%nr_frags+2
2471 if (.not. (fissrate_inout%reac_type .eq.
rrt_nf))
then
2472 fissrate_inout%fissparts(1) = ind_parent
2473 fissrate_inout%fissparts(2) =
ineu
2474 fissrate_inout%ch_amount(1) = -1d0
2475 fissrate_inout%ch_amount(2) = fissfrags_tmp%Yn+float(additional_neutrons)
2477 fissrate_inout%fissparts(1) =
ineu
2478 fissrate_inout%fissparts(2) = ind_parent
2479 fissrate_inout%ch_amount(1) = fissfrags_tmp%Yn-1d0+float(additional_neutrons)
2480 fissrate_inout%ch_amount(2) = -1d0
2484 do i=1,fissfrags_tmp%nr_frags
2485 fissrate_inout%fissparts(i+2) = fissfrags_tmp%net_idx(i)
2486 fissrate_inout%fissparts(i+2) = fissfrags_tmp%net_idx(i)
2487 fissrate_inout%ch_amount(i+2) = fissfrags_tmp%Y(i)
2493 qvalue:
do i=1,fissrate_inout%dimens
2494 qval = qval -
isotope(fissrate_inout%fissparts(i))%mass_exc*fissrate_inout%ch_amount(i)
2496 fissrate_inout%averageQ = qval
2498 info_exit(
"file_fiss_frag")
2513 integer,
intent(in) :: mass,pnr
2514 integer,
intent(in) :: pos
2515 real(
r_kind),
intent(out) :: qval
2516 integer :: fmass,fpnr
2517 integer :: read_stat
2518 integer :: alloc_stat
2520 integer,
dimension(:),
allocatable :: a1,z1
2521 integer,
dimension(:),
allocatable :: a2,z2
2522 integer,
intent(in) :: neutfission
2523 integer :: i,nof, nemiss
2525 real(
r_kind) :: nrate,npre,nafter
2526 real(
r_kind),
dimension(:),
allocatable :: paz
2529 character(2) :: dummy
2531 info_entry(
"abla_nfiss")
2535 read(neutfission,
"(I3,1X,I3,5X,I3)",iostat=read_stat) zf, af, nof
2536 if (read_stat /= 0)
then
2541 read(neutfission,*,iostat=read_stat) nrate, npre, nafter
2542 if (read_stat /= 0)
then
2543 call raise_exception(
"Could not read fission file.",
"abla_nfiss",190004)
2546 if ((af.eq.mass).and.(zf.eq.pnr).and.(nof.ne.0))
then
2548 n_mult = npre + nafter
2549 nemiss = int(dnint(n_mult))
2552 allocate(a1(nof),stat=alloc_stat)
2553 allocate(a2(nof),stat=alloc_stat)
2554 allocate(z1(nof),stat=alloc_stat)
2555 allocate(z2(nof),stat=alloc_stat)
2556 allocate(paz(nof),stat=alloc_stat)
2562 allocate(
fissrate(pos)%channelprob(nof))
2563 allocate(
fissrate(pos)%q_value(nof))
2564 dimens = 2 *
fissrate(pos)%channels + 2
2566 allocate(
fissrate(pos)%fissparts(dimens))
2567 allocate(
fissrate(pos)%ch_amount(dimens))
2569 fissrate(pos)%ch_amount(1) = float(nemiss) - 1.d0
2576 read(neutfission,*,iostat=read_stat) &
2577 & a1(i),z1(i),paz(i)
2578 if (read_stat /= 0)
then
2579 call raise_exception(
"Could not read fission file.",
"abla_nfiss",190004)
2582 if (a1(i).gt.
minmax(z1(i),2))
then
2583 nemiss = nemiss + a1(i) -
minmax(z1(i),2)
2585 else if (a1(i).lt.
minmax(z1(i),1))
then
2586 print *,
'a1 error:'
2588 a2(i) = mass + 1 - a1(i) - nemiss
2590 psum = psum + paz(i)
2593 if (a2(i).gt.
minmax(z2(i),2))
then
2594 nemiss = nemiss + a2(i) -
minmax(z2(i),2)
2596 else if (a2(i).lt.
minmax(z2(i),1))
then
2597 print *,
'a2 error:', a2(i), z2(i), a1(i), z1(i), mass, pnr
2602 paz(i) = paz(i)/psum
2603 fissrate(pos)%channelprob(i) = paz(i)
2618 qvalue:
do i=1,
fissrate(pos)%channels
2630 rewind(neutfission,iostat=read_stat)
2631 if (read_stat /= 0)
then
2633 "abla_nfiss",190005)
2640 read(neutfission,*,iostat=read_stat) dummy
2641 if (read_stat /= 0)
then
2643 "abla_nfiss",190004)
2652 allocate(a1(1),stat=alloc_stat)
2653 allocate(a2(1),stat=alloc_stat)
2654 allocate(z1(1),stat=alloc_stat)
2655 allocate(z2(1),stat=alloc_stat)
2661 a2(1) = nint(fmass/2.d0)
2662 z2(1) = nint(fpnr/2.d0)
2663 a1(1) = fmass - a2(1)
2664 z1(1) = fpnr - z2(1)
2666 a2(1) = nint(fmass/2.d0)
2667 z2(1) = nint(fpnr/2.d0)
2668 a1(1) = fmass - a2(1)
2669 z1(1) = fpnr - z2(1)
2672 z2(1) = nint(52.01 - (fpnr - 80)/1.d1)
2675 a1(1) = fmass - a2(1)
2676 z1(1) = fpnr - z2(1)
2681 if (a1(1).gt.
minmax(z1(1),2))
then
2682 nemiss = a1(1) -
minmax(z1(1),2)
2684 else if (a1(1).lt.
minmax(z1(1),1))
then
2685 print *,
'a1 error:'
2688 if (a2(1).gt.
minmax(z2(1),2))
then
2689 nemiss = nemiss + a2(1) -
minmax(z2(1),2)
2691 else if (a2(1).lt.
minmax(z2(1),1))
then
2692 print *,
'a2 error:'
2699 allocate(
fissrate(pos)%channelprob(1))
2704 fissrate(pos)%channelprob(1) = 1.d0
2707 fissrate(pos)%ch_amount(1) = float(nemiss) - 1.d0
2727 rewind(neutfission,iostat=read_stat)
2728 if (read_stat /= 0)
then
2730 "abla_nfiss",190005)
2733 info_exit(
"abla_nfiss")
2746 integer,
intent(in) :: mass,pnr
2747 integer,
intent(in) :: pos
2748 real(
r_kind),
intent(out) :: qval
2749 integer :: nemiss,fmass,fpnr
2750 integer :: read_stat
2751 integer :: alloc_stat
2753 integer,
dimension(:),
allocatable :: a1,z1
2754 integer,
dimension(:),
allocatable :: a2,z2
2755 integer,
intent(in) :: betafission
2758 real(
r_kind) :: nrate,npre,nafter
2759 real(
r_kind),
dimension(:),
allocatable :: paz
2762 character(2) :: dummy
2764 info_entry(
"abla_betafiss")
2768 read(betafission,*, iostat = read_stat) zf, af, nof
2770 if (read_stat /= 0)
exit
2771 read(betafission,*,iostat=read_stat) nrate, npre, nafter
2772 if (read_stat /= 0)
then
2774 "abla_betafiss",190006)
2776 if ((af.eq.mass).and.(zf.eq.pnr).and.(nof.ne.0))
then
2778 n_mult = npre + nafter
2779 nemiss = int(dnint(n_mult))
2782 allocate(a1(nof),stat=alloc_stat)
2783 allocate(a2(nof),stat=alloc_stat)
2784 allocate(z1(nof),stat=alloc_stat)
2785 allocate(z2(nof),stat=alloc_stat)
2786 allocate(paz(nof),stat=alloc_stat)
2796 allocate(
fissrate(pos)%channelprob(nof))
2797 allocate(
fissrate(pos)%q_value(nof))
2798 dimens = 2 *
fissrate(pos)%channels + 1
2799 if (nemiss.gt.0) dimens = dimens + 1
2801 allocate(
fissrate(pos)%fissparts(dimens))
2802 allocate(
fissrate(pos)%ch_amount(dimens))
2806 if (nemiss.gt.0)
then
2808 fissrate(pos)%ch_amount(dimens) = float(nemiss)
2813 read(betafission,
"(I3,1X,I3,2X,F6.4)",iostat=read_stat) a1(i), z1(i), paz(i)
2814 if (read_stat /= 0)
then
2816 "abla_betafiss",190006)
2820 if (a1(i).gt.
minmax(z1(i),2))
then
2821 nemiss = nemiss + a1(i) -
minmax(z1(i),2)
2823 else if (a1(i).lt.
minmax(z1(i),1))
then
2824 print *,
'a1 error:'
2827 a2(i) = mass - a1(i) - nemiss
2829 if (
fissrate(pos)%mode.eq.3) z2(i) = z2(i) + 1
2830 psum = psum + paz(i)
2833 if (a2(i).gt.
minmax(z2(i),2))
then
2834 nemiss = nemiss + a2(i) -
minmax(z2(i),2)
2836 else if (a2(i).lt.
minmax(z2(i),1))
then
2837 print *,
'a2 error:',a2(i),
minmax(z2(i),1)
2842 paz(i) = paz(i)/psum
2844 fissrate(pos)%channelprob(i) = paz(i)
2859 qvalue:
do i=1,
fissrate(pos)%channels
2870 rewind(betafission,iostat=read_stat)
2871 if (read_stat /= 0)
then
2873 "abla_betafiss",190007)
2880 read(betafission,*,iostat=read_stat) dummy
2881 if (read_stat /= 0)
then
2883 "abla_betafiss",190007)
2892 allocate(a1(1),stat=alloc_stat)
2893 allocate(a2(1),stat=alloc_stat)
2894 allocate(z1(1),stat=alloc_stat)
2895 allocate(z2(1),stat=alloc_stat)
2908 a2(1) = nint(fmass/2.d0)
2909 z2(1) = nint(fpnr/2.d0)
2910 a1(1) = fmass - a2(1)
2911 z1(1) = fpnr - z2(1)
2913 a2(1) = nint(fmass/2.d0)
2914 z2(1) = nint(fpnr/2.d0)
2915 a1(1) = fmass - a2(1)
2916 z1(1) = fpnr - z2(1)
2919 z2(1) = nint(52.01 - (fpnr - 80)/1.d1)
2922 a1(1) = fmass - a2(1)
2923 z1(1) = fpnr - z2(1)
2928 if (a1(1).gt.
minmax(z1(1),2))
then
2929 nemiss = a1(1) -
minmax(z1(1),2)
2931 else if (a1(1).lt.
minmax(z1(1),1))
then
2932 print *,
'a1 error:'
2935 if (a2(1).gt.
minmax(z2(1),2))
then
2936 nemiss = nemiss + a2(1) -
minmax(z2(1),2)
2938 else if (a2(1).lt.
minmax(z2(1),1))
then
2939 print *,
'a2 error:',a2(1),
minmax(z2(1),1)
2943 if (nemiss .eq. 0)
then
2951 allocate(
fissrate(pos)%channelprob(1))
2955 fissrate(pos)%channelprob(1) = 1.d0
2964 if (nemiss.gt.0)
then
2966 fissrate(pos)%ch_amount(4) = float(nemiss)
2980 rewind(betafission,iostat=read_stat)
2981 if (read_stat /= 0)
then
2983 "abla_betafiss",190007)
2986 info_exit(
"abla_betafiss")