ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_distBand.F90
Go to the documentation of this file.
1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3!!!! !!!!
4!!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5!!!! !!!!
6!!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7!!!! !!!!
8!!!! This file is part of the ParaMonte library. !!!!
9!!!! !!!!
10!!!! LICENSE !!!!
11!!!! !!!!
12!!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13!!!! !!!!
14!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16
94
95!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
96
98
99 use pm_kind, only: SK, IK, LK, RKB
100
101 implicit none
102
103 character(*, SK), parameter :: MODULE_NAME = "@pm_distBand"
104
118 real(RKB), parameter :: MEAN_ALPHA = -1.1_RKB
119
130 real(RKB), parameter :: MEAN_BETA = -2.3_RKB
131
132!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
133
168 end type
169
170!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
171
237 interface getBandEpeak
238
239 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
240
241#if RK5_ENABLED
242 PURE elemental module function getBandEpeak_RK5(alpha, beta, ebreak) result(epeak)
243#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
244 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEpeak_RK5
245#endif
246 use pm_kind, only: RKG => RK5
247 real(RKG) , intent(in) :: alpha, beta, ebreak
248 real(RKG) :: epeak
249 end function
250#endif
251
252#if RK4_ENABLED
253 PURE elemental module function getBandEpeak_RK4(alpha, beta, ebreak) result(epeak)
254#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
255 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEpeak_RK4
256#endif
257 use pm_kind, only: RKG => RK4
258 real(RKG) , intent(in) :: alpha, beta, ebreak
259 real(RKG) :: epeak
260 end function
261#endif
262
263#if RK3_ENABLED
264 PURE elemental module function getBandEpeak_RK3(alpha, beta, ebreak) result(epeak)
265#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
266 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEpeak_RK3
267#endif
268 use pm_kind, only: RKG => RK3
269 real(RKG) , intent(in) :: alpha, beta, ebreak
270 real(RKG) :: epeak
271 end function
272#endif
273
274#if RK2_ENABLED
275 PURE elemental module function getBandEpeak_RK2(alpha, beta, ebreak) result(epeak)
276#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
277 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEpeak_RK2
278#endif
279 use pm_kind, only: RKG => RK2
280 real(RKG) , intent(in) :: alpha, beta, ebreak
281 real(RKG) :: epeak
282 end function
283#endif
284
285#if RK1_ENABLED
286 PURE elemental module function getBandEpeak_RK1(alpha, beta, ebreak) result(epeak)
287#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
288 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEpeak_RK1
289#endif
290 use pm_kind, only: RKG => RK1
291 real(RKG) , intent(in) :: alpha, beta, ebreak
292 real(RKG) :: epeak
293 end function
294#endif
295
296 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
297
298 end interface
299
300!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
301
365
366 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
367
368#if RK5_ENABLED
369 PURE elemental module function getBandEbreak_RK5(alpha, beta, epeak) result(ebreak)
370#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
371 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEbreak_RK5
372#endif
373 use pm_kind, only: RKG => RK5
374 real(RKG) , intent(in) :: alpha, beta, epeak
375 real(RKG) :: ebreak
376 end function
377#endif
378
379#if RK4_ENABLED
380 PURE elemental module function getBandEbreak_RK4(alpha, beta, epeak) result(ebreak)
381#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
382 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEbreak_RK4
383#endif
384 use pm_kind, only: RKG => RK4
385 real(RKG) , intent(in) :: alpha, beta, epeak
386 real(RKG) :: ebreak
387 end function
388#endif
389
390#if RK3_ENABLED
391 PURE elemental module function getBandEbreak_RK3(alpha, beta, epeak) result(ebreak)
392#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
393 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEbreak_RK3
394#endif
395 use pm_kind, only: RKG => RK3
396 real(RKG) , intent(in) :: alpha, beta, epeak
397 real(RKG) :: ebreak
398 end function
399#endif
400
401#if RK2_ENABLED
402 PURE elemental module function getBandEbreak_RK2(alpha, beta, epeak) result(ebreak)
403#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
404 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEbreak_RK2
405#endif
406 use pm_kind, only: RKG => RK2
407 real(RKG) , intent(in) :: alpha, beta, epeak
408 real(RKG) :: ebreak
409 end function
410#endif
411
412#if RK1_ENABLED
413 PURE elemental module function getBandEbreak_RK1(alpha, beta, epeak) result(ebreak)
414#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
415 !DEC$ ATTRIBUTES DLLEXPORT :: getBandEbreak_RK1
416#endif
417 use pm_kind, only: RKG => RK1
418 real(RKG) , intent(in) :: alpha, beta, epeak
419 real(RKG) :: ebreak
420 end function
421#endif
422
423 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
424
425 end interface
426
427!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
428
492 interface getBandZeta
493
494 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
495
496#if RK5_ENABLED
497 PURE elemental module function getBandZeta_RK5(alpha, beta, ebreak) result(zeta)
498#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
499 !DEC$ ATTRIBUTES DLLEXPORT :: getBandZeta_RK5
500#endif
501 use pm_kind, only: RKG => RK5
502 real(RKG) , intent(in) :: alpha, beta, ebreak
503 real(RKG) :: zeta
504 end function
505#endif
506
507#if RK4_ENABLED
508 PURE elemental module function getBandZeta_RK4(alpha, beta, ebreak) result(zeta)
509#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
510 !DEC$ ATTRIBUTES DLLEXPORT :: getBandZeta_RK4
511#endif
512 use pm_kind, only: RKG => RK4
513 real(RKG) , intent(in) :: alpha, beta, ebreak
514 real(RKG) :: zeta
515 end function
516#endif
517
518#if RK3_ENABLED
519 PURE elemental module function getBandZeta_RK3(alpha, beta, ebreak) result(zeta)
520#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
521 !DEC$ ATTRIBUTES DLLEXPORT :: getBandZeta_RK3
522#endif
523 use pm_kind, only: RKG => RK3
524 real(RKG) , intent(in) :: alpha, beta, ebreak
525 real(RKG) :: zeta
526 end function
527#endif
528
529#if RK2_ENABLED
530 PURE elemental module function getBandZeta_RK2(alpha, beta, ebreak) result(zeta)
531#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
532 !DEC$ ATTRIBUTES DLLEXPORT :: getBandZeta_RK2
533#endif
534 use pm_kind, only: RKG => RK2
535 real(RKG) , intent(in) :: alpha, beta, ebreak
536 real(RKG) :: zeta
537 end function
538#endif
539
540#if RK1_ENABLED
541 PURE elemental module function getBandZeta_RK1(alpha, beta, ebreak) result(zeta)
542#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
543 !DEC$ ATTRIBUTES DLLEXPORT :: getBandZeta_RK1
544#endif
545 use pm_kind, only: RKG => RK1
546 real(RKG) , intent(in) :: alpha, beta, ebreak
547 real(RKG) :: zeta
548 end function
549#endif
550
551 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
552
553 end interface
554
555!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
556
656 interface getBandUDF
657
658 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
659
660#if RK5_ENABLED
661 PURE elemental module function getBandUDF_RK5(energy, alpha, beta, ebreak, zeta, invEfold) result(udf)
662#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
663 !DEC$ ATTRIBUTES DLLEXPORT :: getBandUDF_RK5
664#endif
665 use pm_kind, only: RKG => RK5
666 real(RKG) , intent(in) :: energy, alpha, beta, ebreak
667 real(RKG) , intent(in) , optional :: zeta, invEfold
668 real(RKG) :: udf
669 end function
670#endif
671
672#if RK4_ENABLED
673 PURE elemental module function getBandUDF_RK4(energy, alpha, beta, ebreak, zeta, invEfold) result(udf)
674#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
675 !DEC$ ATTRIBUTES DLLEXPORT :: getBandUDF_RK4
676#endif
677 use pm_kind, only: RKG => RK4
678 real(RKG) , intent(in) :: energy, alpha, beta, ebreak
679 real(RKG) , intent(in) , optional :: zeta, invEfold
680 real(RKG) :: udf
681 end function
682#endif
683
684#if RK3_ENABLED
685 PURE elemental module function getBandUDF_RK3(energy, alpha, beta, ebreak, zeta, invEfold) result(udf)
686#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
687 !DEC$ ATTRIBUTES DLLEXPORT :: getBandUDF_RK3
688#endif
689 use pm_kind, only: RKG => RK3
690 real(RKG) , intent(in) :: energy, alpha, beta, ebreak
691 real(RKG) , intent(in) , optional :: zeta, invEfold
692 real(RKG) :: udf
693 end function
694#endif
695
696#if RK2_ENABLED
697 PURE elemental module function getBandUDF_RK2(energy, alpha, beta, ebreak, zeta, invEfold) result(udf)
698#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
699 !DEC$ ATTRIBUTES DLLEXPORT :: getBandUDF_RK2
700#endif
701 use pm_kind, only: RKG => RK2
702 real(RKG) , intent(in) :: energy, alpha, beta, ebreak
703 real(RKG) , intent(in) , optional :: zeta, invEfold
704 real(RKG) :: udf
705 end function
706#endif
707
708#if RK1_ENABLED
709 PURE elemental module function getBandUDF_RK1(energy, alpha, beta, ebreak, zeta, invEfold) result(udf)
710#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
711 !DEC$ ATTRIBUTES DLLEXPORT :: getBandUDF_RK1
712#endif
713 use pm_kind, only: RKG => RK1
714 real(RKG) , intent(in) :: energy, alpha, beta, ebreak
715 real(RKG) , intent(in) , optional :: zeta, invEfold
716 real(RKG) :: udf
717 end function
718#endif
719
720 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
721
722 end interface
723
724!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
725
815 interface setBandUCDF
816
817 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
818
819#if RK5_ENABLED
820 impure elemental module subroutine setBandUCDF_RK5(ucdf, lb, ub, alpha, beta, ebreak, info)
821#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
822 !DEC$ ATTRIBUTES DLLEXPORT :: setBandUCDF_RK5
823#endif
824 use pm_kind, only: RKG => RK5
825 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
826 real(RKG) , intent(out) :: ucdf
827 integer(IK) , intent(out) :: info
828 end subroutine
829#endif
830
831#if RK4_ENABLED
832 impure elemental module subroutine setBandUCDF_RK4(ucdf, lb, ub, alpha, beta, ebreak, info)
833#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
834 !DEC$ ATTRIBUTES DLLEXPORT :: setBandUCDF_RK4
835#endif
836 use pm_kind, only: RKG => RK4
837 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
838 real(RKG) , intent(out) :: ucdf
839 integer(IK) , intent(out) :: info
840 end subroutine
841#endif
842
843#if RK3_ENABLED
844 impure elemental module subroutine setBandUCDF_RK3(ucdf, lb, ub, alpha, beta, ebreak, info)
845#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
846 !DEC$ ATTRIBUTES DLLEXPORT :: setBandUCDF_RK3
847#endif
848 use pm_kind, only: RKG => RK3
849 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
850 real(RKG) , intent(out) :: ucdf
851 integer(IK) , intent(out) :: info
852 end subroutine
853#endif
854
855#if RK2_ENABLED
856 impure elemental module subroutine setBandUCDF_RK2(ucdf, lb, ub, alpha, beta, ebreak, info)
857#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
858 !DEC$ ATTRIBUTES DLLEXPORT :: setBandUCDF_RK2
859#endif
860 use pm_kind, only: RKG => RK2
861 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
862 real(RKG) , intent(out) :: ucdf
863 integer(IK) , intent(out) :: info
864 end subroutine
865#endif
866
867#if RK1_ENABLED
868 impure elemental module subroutine setBandUCDF_RK1(ucdf, lb, ub, alpha, beta, ebreak, info)
869#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
870 !DEC$ ATTRIBUTES DLLEXPORT :: setBandUCDF_RK1
871#endif
872 use pm_kind, only: RKG => RK1
873 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
874 real(RKG) , intent(out) :: ucdf
875 integer(IK) , intent(out) :: info
876 end subroutine
877#endif
878
879 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
880
881 end interface
882
883!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
884
1008 interface setBandMean
1009
1010 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1011
1012#if RK5_ENABLED
1013 impure elemental module subroutine setBandMeanDef_RK5(mean, lb, ub, alpha, beta, ebreak, info)
1014#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1015 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanDef_RK5
1016#endif
1017 use pm_kind, only: RKG => RK5
1018 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
1019 real(RKG) , intent(out) :: mean
1020 integer(IK) , intent(out) :: info
1021 end subroutine
1022#endif
1023
1024#if RK4_ENABLED
1025 impure elemental module subroutine setBandMeanDef_RK4(mean, lb, ub, alpha, beta, ebreak, info)
1026#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1027 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanDef_RK4
1028#endif
1029 use pm_kind, only: RKG => RK4
1030 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
1031 real(RKG) , intent(out) :: mean
1032 integer(IK) , intent(out) :: info
1033 end subroutine
1034#endif
1035
1036#if RK3_ENABLED
1037 impure elemental module subroutine setBandMeanDef_RK3(mean, lb, ub, alpha, beta, ebreak, info)
1038#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1039 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanDef_RK3
1040#endif
1041 use pm_kind, only: RKG => RK3
1042 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
1043 real(RKG) , intent(out) :: mean
1044 integer(IK) , intent(out) :: info
1045 end subroutine
1046#endif
1047
1048#if RK2_ENABLED
1049 impure elemental module subroutine setBandMeanDef_RK2(mean, lb, ub, alpha, beta, ebreak, info)
1050#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1051 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanDef_RK2
1052#endif
1053 use pm_kind, only: RKG => RK2
1054 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
1055 real(RKG) , intent(out) :: mean
1056 integer(IK) , intent(out) :: info
1057 end subroutine
1058#endif
1059
1060#if RK1_ENABLED
1061 impure elemental module subroutine setBandMeanDef_RK1(mean, lb, ub, alpha, beta, ebreak, info)
1062#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1063 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanDef_RK1
1064#endif
1065 use pm_kind, only: RKG => RK1
1066 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak
1067 real(RKG) , intent(out) :: mean
1068 integer(IK) , intent(out) :: info
1069 end subroutine
1070#endif
1071
1072 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1073
1074#if RK5_ENABLED
1075 impure elemental module subroutine setBandMeanNew_RK5(mean, lb, ub, alpha, beta, ebreak, lbnew, ubnew, info)
1076#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1077 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanNew_RK5
1078#endif
1079 use pm_kind, only: RKG => RK5
1080 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak, lbnew, ubnew
1081 real(RKG) , intent(out) :: mean
1082 integer(IK) , intent(out) :: info
1083 end subroutine
1084#endif
1085
1086#if RK4_ENABLED
1087 impure elemental module subroutine setBandMeanNew_RK4(mean, lb, ub, alpha, beta, ebreak, lbnew, ubnew, info)
1088#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1089 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanNew_RK4
1090#endif
1091 use pm_kind, only: RKG => RK4
1092 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak, lbnew, ubnew
1093 real(RKG) , intent(out) :: mean
1094 integer(IK) , intent(out) :: info
1095 end subroutine
1096#endif
1097
1098#if RK3_ENABLED
1099 impure elemental module subroutine setBandMeanNew_RK3(mean, lb, ub, alpha, beta, ebreak, lbnew, ubnew, info)
1100#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1101 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanNew_RK3
1102#endif
1103 use pm_kind, only: RKG => RK3
1104 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak, lbnew, ubnew
1105 real(RKG) , intent(out) :: mean
1106 integer(IK) , intent(out) :: info
1107 end subroutine
1108#endif
1109
1110#if RK2_ENABLED
1111 impure elemental module subroutine setBandMeanNew_RK2(mean, lb, ub, alpha, beta, ebreak, lbnew, ubnew, info)
1112#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1113 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanNew_RK2
1114#endif
1115 use pm_kind, only: RKG => RK2
1116 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak, lbnew, ubnew
1117 real(RKG) , intent(out) :: mean
1118 integer(IK) , intent(out) :: info
1119 end subroutine
1120#endif
1121
1122#if RK1_ENABLED
1123 impure elemental module subroutine setBandMeanNew_RK1(mean, lb, ub, alpha, beta, ebreak, lbnew, ubnew, info)
1124#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1125 !DEC$ ATTRIBUTES DLLEXPORT :: setBandMeanNew_RK1
1126#endif
1127 use pm_kind, only: RKG => RK1
1128 real(RKG) , intent(in) :: lb, ub, alpha, beta, ebreak, lbnew, ubnew
1129 real(RKG) , intent(out) :: mean
1130 integer(IK) , intent(out) :: info
1131 end subroutine
1132#endif
1133
1134 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1135
1136 end interface
1137
1138!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1139
1277
1278 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1279
1280#if RK5_ENABLED
1281 impure elemental module subroutine setBandPhotonFromEnergyOldB_RK5(photon, energy, lb, ub, alpha, beta, ebreak, info)
1282#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1283 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyOldB_RK5
1284#endif
1285 use pm_kind, only: RKG => RK5
1286 real(RKG) , intent(in) :: energy, lb, ub, alpha, beta, ebreak
1287 real(RKG) , intent(out) :: photon
1288 integer(IK) , intent(out) :: info
1289 end subroutine
1290#endif
1291
1292#if RK4_ENABLED
1293 impure elemental module subroutine setBandPhotonFromEnergyOldB_RK4(photon, energy, lb, ub, alpha, beta, ebreak, info)
1294#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1295 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyOldB_RK4
1296#endif
1297 use pm_kind, only: RKG => RK4
1298 real(RKG) , intent(in) :: energy, lb, ub, alpha, beta, ebreak
1299 real(RKG) , intent(out) :: photon
1300 integer(IK) , intent(out) :: info
1301 end subroutine
1302#endif
1303
1304#if RK3_ENABLED
1305 impure elemental module subroutine setBandPhotonFromEnergyOldB_RK3(photon, energy, lb, ub, alpha, beta, ebreak, info)
1306#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1307 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyOldB_RK3
1308#endif
1309 use pm_kind, only: RKG => RK3
1310 real(RKG) , intent(in) :: energy, lb, ub, alpha, beta, ebreak
1311 real(RKG) , intent(out) :: photon
1312 integer(IK) , intent(out) :: info
1313 end subroutine
1314#endif
1315
1316#if RK2_ENABLED
1317 impure elemental module subroutine setBandPhotonFromEnergyOldB_RK2(photon, energy, lb, ub, alpha, beta, ebreak, info)
1318#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1319 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyOldB_RK2
1320#endif
1321 use pm_kind, only: RKG => RK2
1322 real(RKG) , intent(in) :: energy, lb, ub, alpha, beta, ebreak
1323 real(RKG) , intent(out) :: photon
1324 integer(IK) , intent(out) :: info
1325 end subroutine
1326#endif
1327
1328#if RK1_ENABLED
1329 impure elemental module subroutine setBandPhotonFromEnergyOldB_RK1(photon, energy, lb, ub, alpha, beta, ebreak, info)
1330#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1331 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyOldB_RK1
1332#endif
1333 use pm_kind, only: RKG => RK1
1334 real(RKG) , intent(in) :: energy, lb, ub, alpha, beta, ebreak
1335 real(RKG) , intent(out) :: photon
1336 integer(IK) , intent(out) :: info
1337 end subroutine
1338#endif
1339
1340 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1341
1342#if RK5_ENABLED
1343 impure elemental module subroutine setBandPhotonFromPhotonNewB_RK5(photon, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1344#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1345 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromPhotonNewB_RK5
1346#endif
1347 use pm_kind, only: RKG => RK5
1348 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1349 real(RKG) , intent(inout) :: photon
1350 integer(IK) , intent(out) :: info
1351 end subroutine
1352#endif
1353
1354#if RK4_ENABLED
1355 impure elemental module subroutine setBandPhotonFromPhotonNewB_RK4(photon, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1356#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1357 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromPhotonNewB_RK4
1358#endif
1359 use pm_kind, only: RKG => RK4
1360 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1361 real(RKG) , intent(inout) :: photon
1362 integer(IK) , intent(out) :: info
1363 end subroutine
1364#endif
1365
1366#if RK3_ENABLED
1367 impure elemental module subroutine setBandPhotonFromPhotonNewB_RK3(photon, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1368#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1369 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromPhotonNewB_RK3
1370#endif
1371 use pm_kind, only: RKG => RK3
1372 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1373 real(RKG) , intent(inout) :: photon
1374 integer(IK) , intent(out) :: info
1375 end subroutine
1376#endif
1377
1378#if RK2_ENABLED
1379 impure elemental module subroutine setBandPhotonFromPhotonNewB_RK2(photon, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1380#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1381 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromPhotonNewB_RK2
1382#endif
1383 use pm_kind, only: RKG => RK2
1384 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1385 real(RKG) , intent(inout) :: photon
1386 integer(IK) , intent(out) :: info
1387 end subroutine
1388#endif
1389
1390#if RK1_ENABLED
1391 impure elemental module subroutine setBandPhotonFromPhotonNewB_RK1(photon, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1392#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1393 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromPhotonNewB_RK1
1394#endif
1395 use pm_kind, only: RKG => RK1
1396 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1397 real(RKG) , intent(inout) :: photon
1398 integer(IK) , intent(out) :: info
1399 end subroutine
1400#endif
1401
1402 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1403
1404#if RK5_ENABLED
1405 impure elemental module subroutine setBandPhotonFromEnergyNewB_RK5(photon, lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak, info)
1406#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1407 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyNewB_RK5
1408#endif
1409 use pm_kind, only: RKG => RK5
1410 real(RKG) , intent(in) :: lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak
1411 real(RKG) , intent(out) :: photon
1412 integer(IK) , intent(out) :: info
1413 end subroutine
1414#endif
1415
1416#if RK4_ENABLED
1417 impure elemental module subroutine setBandPhotonFromEnergyNewB_RK4(photon, lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak, info)
1418#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1419 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyNewB_RK4
1420#endif
1421 use pm_kind, only: RKG => RK4
1422 real(RKG) , intent(in) :: lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak
1423 real(RKG) , intent(out) :: photon
1424 integer(IK) , intent(out) :: info
1425 end subroutine
1426#endif
1427
1428#if RK3_ENABLED
1429 impure elemental module subroutine setBandPhotonFromEnergyNewB_RK3(photon, lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak, info)
1430#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1431 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyNewB_RK3
1432#endif
1433 use pm_kind, only: RKG => RK3
1434 real(RKG) , intent(in) :: lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak
1435 real(RKG) , intent(out) :: photon
1436 integer(IK) , intent(out) :: info
1437 end subroutine
1438#endif
1439
1440#if RK2_ENABLED
1441 impure elemental module subroutine setBandPhotonFromEnergyNewB_RK2(photon, lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak, info)
1442#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1443 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyNewB_RK2
1444#endif
1445 use pm_kind, only: RKG => RK2
1446 real(RKG) , intent(in) :: lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak
1447 real(RKG) , intent(out) :: photon
1448 integer(IK) , intent(out) :: info
1449 end subroutine
1450#endif
1451
1452#if RK1_ENABLED
1453 impure elemental module subroutine setBandPhotonFromEnergyNewB_RK1(photon, lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak, info)
1454#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1455 !DEC$ ATTRIBUTES DLLEXPORT :: setBandPhotonFromEnergyNewB_RK1
1456#endif
1457 use pm_kind, only: RKG => RK1
1458 real(RKG) , intent(in) :: lbnew, ubnew, energy, lb, ub, alpha, beta, ebreak
1459 real(RKG) , intent(out) :: photon
1460 integer(IK) , intent(out) :: info
1461 end subroutine
1462#endif
1463
1464 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1465
1466 end interface
1467
1468!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1469
1606
1607 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1608
1609#if RK5_ENABLED
1610 impure elemental module subroutine setBandEnergyFromPhotonOldB_RK5(energy, photon, lb, ub, alpha, beta, ebreak, info)
1611#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1612 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonOldB_RK5
1613#endif
1614 use pm_kind, only: RKG => RK5
1615 real(RKG) , intent(in) :: photon, lb, ub, alpha, beta, ebreak
1616 real(RKG) , intent(out) :: energy
1617 integer(IK) , intent(out) :: info
1618 end subroutine
1619#endif
1620
1621#if RK4_ENABLED
1622 impure elemental module subroutine setBandEnergyFromPhotonOldB_RK4(energy, photon, lb, ub, alpha, beta, ebreak, info)
1623#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1624 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonOldB_RK4
1625#endif
1626 use pm_kind, only: RKG => RK4
1627 real(RKG) , intent(in) :: photon, lb, ub, alpha, beta, ebreak
1628 real(RKG) , intent(out) :: energy
1629 integer(IK) , intent(out) :: info
1630 end subroutine
1631#endif
1632
1633#if RK3_ENABLED
1634 impure elemental module subroutine setBandEnergyFromPhotonOldB_RK3(energy, photon, lb, ub, alpha, beta, ebreak, info)
1635#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1636 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonOldB_RK3
1637#endif
1638 use pm_kind, only: RKG => RK3
1639 real(RKG) , intent(in) :: photon, lb, ub, alpha, beta, ebreak
1640 real(RKG) , intent(out) :: energy
1641 integer(IK) , intent(out) :: info
1642 end subroutine
1643#endif
1644
1645#if RK2_ENABLED
1646 impure elemental module subroutine setBandEnergyFromPhotonOldB_RK2(energy, photon, lb, ub, alpha, beta, ebreak, info)
1647#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1648 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonOldB_RK2
1649#endif
1650 use pm_kind, only: RKG => RK2
1651 real(RKG) , intent(in) :: photon, lb, ub, alpha, beta, ebreak
1652 real(RKG) , intent(out) :: energy
1653 integer(IK) , intent(out) :: info
1654 end subroutine
1655#endif
1656
1657#if RK1_ENABLED
1658 impure elemental module subroutine setBandEnergyFromPhotonOldB_RK1(energy, photon, lb, ub, alpha, beta, ebreak, info)
1659#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1660 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonOldB_RK1
1661#endif
1662 use pm_kind, only: RKG => RK1
1663 real(RKG) , intent(in) :: photon, lb, ub, alpha, beta, ebreak
1664 real(RKG) , intent(out) :: energy
1665 integer(IK) , intent(out) :: info
1666 end subroutine
1667#endif
1668
1669 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1670
1671#if RK5_ENABLED
1672 impure elemental module subroutine setBandEnergyFromPhotonNewB_RK5(energy, lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak, info)
1673#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1674 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonNewB_RK5
1675#endif
1676 use pm_kind, only: RKG => RK5
1677 real(RKG) , intent(in) :: lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak
1678 real(RKG) , intent(out) :: energy
1679 integer(IK) , intent(out) :: info
1680 end subroutine
1681#endif
1682
1683#if RK4_ENABLED
1684 impure elemental module subroutine setBandEnergyFromPhotonNewB_RK4(energy, lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak, info)
1685#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1686 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonNewB_RK4
1687#endif
1688 use pm_kind, only: RKG => RK4
1689 real(RKG) , intent(in) :: lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak
1690 real(RKG) , intent(out) :: energy
1691 integer(IK) , intent(out) :: info
1692 end subroutine
1693#endif
1694
1695#if RK3_ENABLED
1696 impure elemental module subroutine setBandEnergyFromPhotonNewB_RK3(energy, lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak, info)
1697#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1698 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonNewB_RK3
1699#endif
1700 use pm_kind, only: RKG => RK3
1701 real(RKG) , intent(in) :: lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak
1702 real(RKG) , intent(out) :: energy
1703 integer(IK) , intent(out) :: info
1704 end subroutine
1705#endif
1706
1707#if RK2_ENABLED
1708 impure elemental module subroutine setBandEnergyFromPhotonNewB_RK2(energy, lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak, info)
1709#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1710 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonNewB_RK2
1711#endif
1712 use pm_kind, only: RKG => RK2
1713 real(RKG) , intent(in) :: lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak
1714 real(RKG) , intent(out) :: energy
1715 integer(IK) , intent(out) :: info
1716 end subroutine
1717#endif
1718
1719#if RK1_ENABLED
1720 impure elemental module subroutine setBandEnergyFromPhotonNewB_RK1(energy, lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak, info)
1721#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1722 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromPhotonNewB_RK1
1723#endif
1724 use pm_kind, only: RKG => RK1
1725 real(RKG) , intent(in) :: lbnew, ubnew, photon, lb, ub, alpha, beta, ebreak
1726 real(RKG) , intent(out) :: energy
1727 integer(IK) , intent(out) :: info
1728 end subroutine
1729#endif
1730
1731 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1732
1733#if RK5_ENABLED
1734 impure elemental module subroutine setBandEnergyFromEnergyNewB_RK5(energy, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1735#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1736 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromEnergyNewB_RK5
1737#endif
1738 use pm_kind, only: RKG => RK5
1739 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1740 real(RKG) , intent(inout) :: energy
1741 integer(IK) , intent(out) :: info
1742 end subroutine
1743#endif
1744
1745#if RK4_ENABLED
1746 impure elemental module subroutine setBandEnergyFromEnergyNewB_RK4(energy, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1747#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1748 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromEnergyNewB_RK4
1749#endif
1750 use pm_kind, only: RKG => RK4
1751 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1752 real(RKG) , intent(inout) :: energy
1753 integer(IK) , intent(out) :: info
1754 end subroutine
1755#endif
1756
1757#if RK3_ENABLED
1758 impure elemental module subroutine setBandEnergyFromEnergyNewB_RK3(energy, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1759#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1760 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromEnergyNewB_RK3
1761#endif
1762 use pm_kind, only: RKG => RK3
1763 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1764 real(RKG) , intent(inout) :: energy
1765 integer(IK) , intent(out) :: info
1766 end subroutine
1767#endif
1768
1769#if RK2_ENABLED
1770 impure elemental module subroutine setBandEnergyFromEnergyNewB_RK2(energy, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1771#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1772 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromEnergyNewB_RK2
1773#endif
1774 use pm_kind, only: RKG => RK2
1775 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1776 real(RKG) , intent(inout) :: energy
1777 integer(IK) , intent(out) :: info
1778 end subroutine
1779#endif
1780
1781#if RK1_ENABLED
1782 impure elemental module subroutine setBandEnergyFromEnergyNewB_RK1(energy, lbnew, ubnew, lb, ub, alpha, beta, ebreak, info)
1783#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1784 !DEC$ ATTRIBUTES DLLEXPORT :: setBandEnergyFromEnergyNewB_RK1
1785#endif
1786 use pm_kind, only: RKG => RK1
1787 real(RKG) , intent(in) :: lbnew, ubnew, lb, ub, alpha, beta, ebreak
1788 real(RKG) , intent(inout) :: energy
1789 integer(IK) , intent(out) :: info
1790 end subroutine
1791#endif
1792
1793 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1794
1795 end interface
1796
1797!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1798
1799! Bizarrely and frustratingly, Microsoft Windows Subsystem for Linux Ubuntu with GFortran yields Segmentation faults with internal procedure calls.
1800! This is so unfortunate. To bypass this issue for now, the following subroutine is implemented as separate submodule
1801! so that the internal shared parameters can be safely passed as submodule parameters.
1802
1803!WSL_GFORTRAN_BUG !> \brief
1804!WSL_GFORTRAN_BUG !> Integrate the Band differential spectrum over the input energy range.
1805!WSL_GFORTRAN_BUG !>
1806!WSL_GFORTRAN_BUG !> \param[in] lowerLim : The lower limit energy (in units of [keV]) of the integration.
1807!WSL_GFORTRAN_BUG !> \param[in] upperLim : The upper limit energy (in units of [keV]) of the integration.
1808!WSL_GFORTRAN_BUG !> \param[in] epk : The spectral peak energy in units of [keV].
1809!WSL_GFORTRAN_BUG !> \param[in] alpha : The lower spectral exponent of the Band model.
1810!WSL_GFORTRAN_BUG !> \param[in] beta : The upper spectral exponent of the Band model.
1811!WSL_GFORTRAN_BUG !> \param[in] tolerance : The relative accuracy tolerance of the integration.
1812!WSL_GFORTRAN_BUG !> \param[out] photonFluence : The fluence in units of photon counts within the input energy range.
1813!WSL_GFORTRAN_BUG !> \param[out] Err : An object of class [err_type](@ref pm_err::err_type) containing error-handling information.
1814!WSL_GFORTRAN_BUG subroutine getBandPhoton(lowerLim,upperLim,epk,alpha,beta,tolerance,photonFluence,Err)
1815!WSL_GFORTRAN_BUG #if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1816!WSL_GFORTRAN_BUG !DEC$ ATTRIBUTES DLLEXPORT :: getBandPhoton
1817!WSL_GFORTRAN_BUG #endif
1818!WSL_GFORTRAN_BUG
1819!WSL_GFORTRAN_BUG !use pm_quadRomb, only: getQuadRomb
1820!WSL_GFORTRAN_BUG !WSL_GFORTRAN_BUG use QuadPackRK1_pmod, only: qag
1821!WSL_GFORTRAN_BUG use pm_err, only: err_type
1822!WSL_GFORTRAN_BUG implicit none
1823!WSL_GFORTRAN_BUG real(RK) , intent(in) :: lowerLim, upperLim, epk, alpha, beta, tolerance
1824!WSL_GFORTRAN_BUG real(RK) , intent(out) :: photonFluence
1825!WSL_GFORTRAN_BUG type(err_type) , intent(out) :: Err
1826!WSL_GFORTRAN_BUG character(*, SK), parameter :: PROCEDURE_NAME = "@getBandPhoton()"
1827!WSL_GFORTRAN_BUG real(RK) :: ebrk, alphaPlusTwo
1828!WSL_GFORTRAN_BUG real(RK) :: thisUpperLim, alphaPlusTwoOverEpk, betaPlusOne
1829!WSL_GFORTRAN_BUG real(RK) :: alphaMinusBand, coef
1830!WSL_GFORTRAN_BUG real(RK) :: abserr
1831!WSL_GFORTRAN_BUG integer(IK) :: neval
1832!WSL_GFORTRAN_BUG integer(IK) :: ierr
1833!WSL_GFORTRAN_BUG !real(RK) :: alphaPlusOne, logGammaAlphaPlusOne
1834!WSL_GFORTRAN_BUG
1835!WSL_GFORTRAN_BUG err%occurred = .false._LK
1836!WSL_GFORTRAN_BUG
1837!WSL_GFORTRAN_BUG if (lowerLim>=upperLim) then
1838!WSL_GFORTRAN_BUG photonFluence = 0._RK
1839!WSL_GFORTRAN_BUG return
1840!WSL_GFORTRAN_BUG end if
1841!WSL_GFORTRAN_BUG
1842!WSL_GFORTRAN_BUG ! check if the photon indices are consistent with the mathematical rules
1843!WSL_GFORTRAN_BUG if (alpha<beta .or. alpha<-2._RK) then
1844!WSL_GFORTRAN_BUG photonFluence = -HUGE_RK
1845!WSL_GFORTRAN_BUG err%occurred = .true._LK
1846!WSL_GFORTRAN_BUG err%msg = MODULE_NAME//PROCEDURE_NAME//SK_": Error occurred: alpha<beta .or. alpha<-2._RK"
1847!WSL_GFORTRAN_BUG return
1848!WSL_GFORTRAN_BUG end if
1849!WSL_GFORTRAN_BUG
1850!WSL_GFORTRAN_BUG ! integrate the spectrum
1851!WSL_GFORTRAN_BUG alphaPlusTwo = alpha + 2._RK
1852!WSL_GFORTRAN_BUG alphaMinusBand = alpha - beta
1853!WSL_GFORTRAN_BUG ebrk = epk*alphaMinusBand/alphaPlusTwo
1854!WSL_GFORTRAN_BUG
1855!WSL_GFORTRAN_BUG if (lowerLim>ebrk) then
1856!WSL_GFORTRAN_BUG
1857!WSL_GFORTRAN_BUG ! there is only the high energy component in the photonFluence
1858!WSL_GFORTRAN_BUG betaPlusOne = beta + 1._RK
1859!WSL_GFORTRAN_BUG coef = ebrk**alphaMinusBand * exp(-alphaMinusBand);
1860!WSL_GFORTRAN_BUG photonFluence = coef * ( upperLim**betaPlusOne - lowerLim**betaPlusOne ) / betaPlusOne
1861!WSL_GFORTRAN_BUG return
1862!WSL_GFORTRAN_BUG
1863!WSL_GFORTRAN_BUG !#if WSL_ENABLED && CODECOV_ENABLED
1864!WSL_GFORTRAN_BUG !! LCOV_EXCL_START
1865!WSL_GFORTRAN_BUG !#endif
1866!WSL_GFORTRAN_BUG
1867!WSL_GFORTRAN_BUG elseif (lowerLim<ebrk) then
1868!WSL_GFORTRAN_BUG
1869!WSL_GFORTRAN_BUG alphaPlusTwoOverEpk = alphaPlusTwo / epk
1870!WSL_GFORTRAN_BUG thisUpperLim = min(upperLim,ebrk)
1871!WSL_GFORTRAN_BUG !alphaPlusOne = alpha + 1._RK
1872!WSL_GFORTRAN_BUG !if (alpha>-1._RK) then
1873!WSL_GFORTRAN_BUG ! logGammaAlphaPlusOne = log_gamma( alphaPlusOne )
1874!WSL_GFORTRAN_BUG ! ! use the analytical approach to compute the photonFluence:
1875!WSL_GFORTRAN_BUG ! ! https://www.wolframalpha.com/input/?i=integrate+x%5Ea+*+exp(-b*x)
1876!WSL_GFORTRAN_BUG ! photonFluence = getGamUpper( exponent = alphaPlusOne &
1877!WSL_GFORTRAN_BUG ! , logGammaExponent = logGammaAlphaPlusOne &
1878!WSL_GFORTRAN_BUG ! , lowerLim = alphaPlusTwoOverEpk * lowerLim &
1879!WSL_GFORTRAN_BUG ! , tolerance = tolerance &
1880!WSL_GFORTRAN_BUG ! ) &
1881!WSL_GFORTRAN_BUG ! - getGamUpper( exponent = alphaPlusOne &
1882!WSL_GFORTRAN_BUG ! , logGammaExponent = logGammaAlphaPlusOne &
1883!WSL_GFORTRAN_BUG ! , lowerLim = alphaPlusTwoOverEpk * thisUpperLim &
1884!WSL_GFORTRAN_BUG ! , tolerance = tolerance &
1885!WSL_GFORTRAN_BUG ! )
1886!WSL_GFORTRAN_BUG ! photonFluence = photonFluence / alphaPlusTwoOverEpk**alphaPlusOne
1887!WSL_GFORTRAN_BUG !else
1888!WSL_GFORTRAN_BUG ! use brute-force integration
1889!WSL_GFORTRAN_BUG call qag( f = getBandCompLowPhoton &
1890!WSL_GFORTRAN_BUG , a = lowerLim &
1891!WSL_GFORTRAN_BUG , b = thisUpperLim &
1892!WSL_GFORTRAN_BUG , epsabs = 0._RK &
1893!WSL_GFORTRAN_BUG , epsrel = tolerance &
1894!WSL_GFORTRAN_BUG , key = 1_IK &
1895!WSL_GFORTRAN_BUG , result = photonFluence &
1896!WSL_GFORTRAN_BUG , abserr = abserr &
1897!WSL_GFORTRAN_BUG , neval = neval &
1898!WSL_GFORTRAN_BUG , ier = ierr &
1899!WSL_GFORTRAN_BUG )
1900!WSL_GFORTRAN_BUG !write(*,*) neval
1901!WSL_GFORTRAN_BUG !call getQuadRomb ( getFunc = getBandCompLowPhoton &
1902!WSL_GFORTRAN_BUG ! , xmin = lowerLim &
1903!WSL_GFORTRAN_BUG ! , xmax = thisUpperLim &
1904!WSL_GFORTRAN_BUG ! , tolerance = 1.e-7_RK &
1905!WSL_GFORTRAN_BUG ! , nRefinement = 10_IK &
1906!WSL_GFORTRAN_BUG ! , photonFluence = photonFluence &
1907!WSL_GFORTRAN_BUG ! , ierr = ierr &
1908!WSL_GFORTRAN_BUG ! )
1909!WSL_GFORTRAN_BUG if (ierr/=0_IK) then
1910!WSL_GFORTRAN_BUG photonFluence = -HUGE_RK
1911!WSL_GFORTRAN_BUG err%occurred = .true._LK
1912!WSL_GFORTRAN_BUG err%stat = ierr
1913!WSL_GFORTRAN_BUG err%msg = MODULE_NAME//PROCEDURE_NAME//SK_": Error occurred at QuadPack routine. Check the error code to identify the root cause."
1914!WSL_GFORTRAN_BUG return
1915!WSL_GFORTRAN_BUG end if
1916!WSL_GFORTRAN_BUG !end if
1917!WSL_GFORTRAN_BUG
1918!WSL_GFORTRAN_BUG if (upperLim>ebrk) then
1919!WSL_GFORTRAN_BUG ! add the remaining part of the photonFluence from the high-energy component
1920!WSL_GFORTRAN_BUG betaPlusOne = beta + 1._RK
1921!WSL_GFORTRAN_BUG alphaMinusBand = alpha - beta
1922!WSL_GFORTRAN_BUG coef = ebrk**alphaMinusBand * exp(-alphaMinusBand);
1923!WSL_GFORTRAN_BUG photonFluence = photonFluence + coef * ( upperLim**betaPlusOne - ebrk**betaPlusOne ) / betaPlusOne
1924!WSL_GFORTRAN_BUG return
1925!WSL_GFORTRAN_BUG end if
1926!WSL_GFORTRAN_BUG
1927!WSL_GFORTRAN_BUG end if
1928!WSL_GFORTRAN_BUG
1929!WSL_GFORTRAN_BUG contains
1930!WSL_GFORTRAN_BUG
1931!WSL_GFORTRAN_BUG pure function getBandCompLowPhoton(energy) result(bandCompLow)
1932!WSL_GFORTRAN_BUG implicit none
1933!WSL_GFORTRAN_BUG real(RK) , intent(in) :: energy
1934!WSL_GFORTRAN_BUG real(RK) :: bandCompLow
1935!WSL_GFORTRAN_BUG bandCompLow = energy**alpha * exp(-alphaPlusTwoOverEpk*energy)
1936!WSL_GFORTRAN_BUG end function getBandCompLowPhoton
1937!WSL_GFORTRAN_BUG
1938!WSL_GFORTRAN_BUG !#if WSL_ENABLED && CODECOV_ENABLED
1939!WSL_GFORTRAN_BUG !! LCOV_EXCL_STOP
1940!WSL_GFORTRAN_BUG !#endif
1941!WSL_GFORTRAN_BUG
1942!WSL_GFORTRAN_BUG end subroutine getBandPhoton
1943
1944!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1945
1946! Bizarrely and frustratingly, Microsoft Windows Subsystem for Linux Ubuntu with GFortran yields Segmentation faults with internal procedure calls.
1947! This is so unfortunate. To bypass this issue for now, the following subroutine is implemented as separate submodule
1948! so that the internal shared parameters can be safely passed as submodule parameters.
1949
1950!WSL_GFORTRAN_BUG !> \brief
1951!WSL_GFORTRAN_BUG !> Integrate the Band differential spectrum over the input energy range in units of the input energy.
1952!WSL_GFORTRAN_BUG !>
1953!WSL_GFORTRAN_BUG !> \param[in] lowerLim : The lower limit energy (in units of [keV]) of the integration.
1954!WSL_GFORTRAN_BUG !> \param[in] upperLim : The upper limit energy (in units of [keV]) of the integration.
1955!WSL_GFORTRAN_BUG !> \param[in] epk : The spectral peak energy in units of [keV].
1956!WSL_GFORTRAN_BUG !> \param[in] alpha : The lower spectral exponent of the Band model.
1957!WSL_GFORTRAN_BUG !> \param[in] beta : The upper spectral exponent of the Band model.
1958!WSL_GFORTRAN_BUG !> \param[in] tolerance : The relative accuracy tolerance of the integration.
1959!WSL_GFORTRAN_BUG !> \param[out] energyFluence : The fluence in units of the input energy (keV) within the input energy range.
1960!WSL_GFORTRAN_BUG !> \param[out] Err : An object of class [err_type](@ref pm_err::err_type) containing error-handling information.
1961!WSL_GFORTRAN_BUG subroutine getFluenceEnergy(lowerLim,upperLim,epk,alpha,beta,tolerance,energyFluence,Err)
1962!WSL_GFORTRAN_BUG #if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1963!WSL_GFORTRAN_BUG !DEC$ ATTRIBUTES DLLEXPORT :: getFluenceEnergy
1964!WSL_GFORTRAN_BUG #endif
1965!WSL_GFORTRAN_BUG
1966!WSL_GFORTRAN_BUG !WSL_GFORTRAN_BUG use QuadPackRK1_pmod, only: qag
1967!WSL_GFORTRAN_BUG use pm_err, only: err_type
1968!WSL_GFORTRAN_BUG implicit none
1969!WSL_GFORTRAN_BUG real(RK) , intent(in) :: lowerLim, upperLim, epk, alpha, beta, tolerance
1970!WSL_GFORTRAN_BUG type(err_type) , intent(out) :: Err
1971!WSL_GFORTRAN_BUG real(RK) , intent(out) :: energyFluence
1972!WSL_GFORTRAN_BUG character(*, SK), parameter :: PROCEDURE_NAME = "@getFluenceEnergy()"
1973!WSL_GFORTRAN_BUG real(RK) :: ebrk, alphaPlusTwo
1974!WSL_GFORTRAN_BUG real(RK) :: thisUpperLim, alphaPlusTwoOverEpk, betaPlusTwo
1975!WSL_GFORTRAN_BUG real(RK) :: alphaMinusBand, coef, alphaPlusOne
1976!WSL_GFORTRAN_BUG real(RK) :: abserr
1977!WSL_GFORTRAN_BUG integer(IK) :: neval
1978!WSL_GFORTRAN_BUG integer(IK) :: ierr
1979!WSL_GFORTRAN_BUG
1980!WSL_GFORTRAN_BUG err%occurred = .false._LK
1981!WSL_GFORTRAN_BUG
1982!WSL_GFORTRAN_BUG if (lowerLim>=upperLim) then
1983!WSL_GFORTRAN_BUG energyFluence = 0._RK
1984!WSL_GFORTRAN_BUG return
1985!WSL_GFORTRAN_BUG end if
1986!WSL_GFORTRAN_BUG
1987!WSL_GFORTRAN_BUG ! check if the photon indices are consistent with the mathematical rules
1988!WSL_GFORTRAN_BUG if (alpha<beta .or. alpha<-2._RK) then
1989!WSL_GFORTRAN_BUG energyFluence = -HUGE_RK
1990!WSL_GFORTRAN_BUG err%occurred = .true._LK
1991!WSL_GFORTRAN_BUG err%msg = MODULE_NAME//PROCEDURE_NAME//SK_": Error occurred: alpha<beta .or. alpha<-2._RK"
1992!WSL_GFORTRAN_BUG return
1993!WSL_GFORTRAN_BUG end if
1994!WSL_GFORTRAN_BUG
1995!WSL_GFORTRAN_BUG ! integrate the spectrum
1996!WSL_GFORTRAN_BUG alphaPlusTwo = alpha + 2._RK
1997!WSL_GFORTRAN_BUG alphaMinusBand = alpha - beta
1998!WSL_GFORTRAN_BUG ebrk = epk*alphaMinusBand/alphaPlusTwo
1999!WSL_GFORTRAN_BUG
2000!WSL_GFORTRAN_BUG if (lowerLim>ebrk) then
2001!WSL_GFORTRAN_BUG
2002!WSL_GFORTRAN_BUG ! there is only the high energy component in the energyFluence
2003!WSL_GFORTRAN_BUG betaPlusTwo = beta + 2._RK
2004!WSL_GFORTRAN_BUG coef = ebrk**alphaMinusBand * exp(-alphaMinusBand);
2005!WSL_GFORTRAN_BUG energyFluence = coef * ( upperLim**betaPlusTwo - lowerLim**betaPlusTwo ) / betaPlusTwo
2006!WSL_GFORTRAN_BUG return
2007!WSL_GFORTRAN_BUG
2008!WSL_GFORTRAN_BUG !#if WSL_ENABLED && CODECOV_ENABLED
2009!WSL_GFORTRAN_BUG !! LCOV_EXCL_START
2010!WSL_GFORTRAN_BUG !#endif
2011!WSL_GFORTRAN_BUG
2012!WSL_GFORTRAN_BUG elseif (lowerLim<ebrk) then
2013!WSL_GFORTRAN_BUG
2014!WSL_GFORTRAN_BUG alphaPlusTwoOverEpk = alphaPlusTwo / epk
2015!WSL_GFORTRAN_BUG thisUpperLim = min(upperLim,ebrk)
2016!WSL_GFORTRAN_BUG alphaPlusOne = alpha + 1._RK
2017!WSL_GFORTRAN_BUG ! use brute-force integration
2018!WSL_GFORTRAN_BUG call qag( f = getBandCompLowEnergy &
2019!WSL_GFORTRAN_BUG , a = lowerLim &
2020!WSL_GFORTRAN_BUG , b = thisUpperLim &
2021!WSL_GFORTRAN_BUG , epsabs = 0._RK &
2022!WSL_GFORTRAN_BUG , epsrel = tolerance &
2023!WSL_GFORTRAN_BUG , key = 1_IK &
2024!WSL_GFORTRAN_BUG , result = energyFluence &
2025!WSL_GFORTRAN_BUG , abserr = abserr &
2026!WSL_GFORTRAN_BUG , neval = neval &
2027!WSL_GFORTRAN_BUG , ier = ierr &
2028!WSL_GFORTRAN_BUG )
2029!WSL_GFORTRAN_BUG
2030!WSL_GFORTRAN_BUG if (ierr/=0_IK) then
2031!WSL_GFORTRAN_BUG energyFluence = -HUGE_RK
2032!WSL_GFORTRAN_BUG err%occurred = .true._LK
2033!WSL_GFORTRAN_BUG err%stat = ierr
2034!WSL_GFORTRAN_BUG err%msg = MODULE_NAME//PROCEDURE_NAME//SK_": Error occurred at QuadPack routine. Check the error code to identify the root cause."
2035!WSL_GFORTRAN_BUG return
2036!WSL_GFORTRAN_BUG end if
2037!WSL_GFORTRAN_BUG
2038!WSL_GFORTRAN_BUG if (upperLim>ebrk) then ! add the remaining part of the energyFluence from the high-energy component
2039!WSL_GFORTRAN_BUG betaPlusTwo = beta + 2._RK
2040!WSL_GFORTRAN_BUG alphaMinusBand = alpha - beta
2041!WSL_GFORTRAN_BUG coef = ebrk**alphaMinusBand * exp(-alphaMinusBand);
2042!WSL_GFORTRAN_BUG energyFluence = energyFluence + coef * ( upperLim**betaPlusTwo - ebrk**betaPlusTwo ) / betaPlusTwo
2043!WSL_GFORTRAN_BUG return
2044!WSL_GFORTRAN_BUG end if
2045!WSL_GFORTRAN_BUG
2046!WSL_GFORTRAN_BUG end if
2047!WSL_GFORTRAN_BUG
2048!WSL_GFORTRAN_BUG contains
2049!WSL_GFORTRAN_BUG
2050!WSL_GFORTRAN_BUG pure function getBandCompLowEnergy(energy) result(bandCompLow)
2051!WSL_GFORTRAN_BUG implicit none
2052!WSL_GFORTRAN_BUG real(RK) , intent(in) :: energy
2053!WSL_GFORTRAN_BUG real(RK) :: bandCompLow
2054!WSL_GFORTRAN_BUG bandCompLow = energy**alphaPlusOne * exp(-alphaPlusTwoOverEpk*energy)
2055!WSL_GFORTRAN_BUG end function getBandCompLowEnergy
2056!WSL_GFORTRAN_BUG
2057!WSL_GFORTRAN_BUG !#if WSL_ENABLED && CODECOV_ENABLED
2058!WSL_GFORTRAN_BUG !! LCOV_EXCL_STOP
2059!WSL_GFORTRAN_BUG !#endif
2060!WSL_GFORTRAN_BUG
2061!WSL_GFORTRAN_BUG end subroutine getFluenceEnergy
2062
2063!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2064
2065end module pm_distBand ! LCOV_EXCL_LINE
Generate and return the spectral break energy parameter of the Band spectral model/distribution from ...
Generate and return the spectral peak energy parameter of the Band spectral model/distribution from t...
Generate and return the unnormalized density function (UDF) of the Band spectral model/distribution.
Generate and return the coefficient of continuity of the Band spectral model/distribution from the Ba...
Generate and return the energy integral (the energy fluence in units of the input break energy) of th...
Generate and return the mean of the Band distribution for an input set of parameters.
Generate and return the photon integral (the photon fluence in units of photon counts) of the Band mo...
Generate and return the unnormalized cumulative distribution function (UCDF) of the Band spectral mod...
This module contains procedures and generic interfaces for computing the Band photon distribution wid...
Definition: pm_distBand.F90:97
character(*, SK), parameter MODULE_NAME
real(RKB), parameter MEAN_ALPHA
The scalar constant of type real of kind RKB, containing the average reported value for the low-energ...
real(RKB), parameter MEAN_BETA
The scalar constant of type real of kind RKB, containing the average reported value for the high-ener...
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter RK5
Definition: pm_kind.F90:478
integer, parameter RK4
Definition: pm_kind.F90:489
integer, parameter RK2
Definition: pm_kind.F90:511
integer, parameter RK3
Definition: pm_kind.F90:500
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter RKB
The scalar integer constant of intrinsic default kind, representing the Best-precision real kind supp...
Definition: pm_kind.F90:1371
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RK1
Definition: pm_kind.F90:522
This is the derived type for signifying distributions that are of type Band as defined in the descrip...