ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_mathGamma.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
67
69
70 use pm_kind, only: SK, IK, LK
71
72 implicit none
73
74 character(*, SK), parameter :: MODULE_NAME = "@pm_mathGamma"
75
76!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
77
161
162 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
163
164#if RK5_ENABLED
165 impure elemental module function getGammaIncLow_RK5(x, kappa, tol) result(gammaIncLow)
166#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
167 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK5
168#endif
169 use pm_kind, only: RKG => RK5
170 real(RKG) , intent(in) :: x, kappa
171 real(RKG) , intent(in) , optional :: tol
172 real(RKG) :: gammaIncLow
173 end function
174#endif
175
176#if RK4_ENABLED
177 impure elemental module function getGammaIncLow_RK4(x, kappa, tol) result(gammaIncLow)
178#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
179 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK4
180#endif
181 use pm_kind, only: RKG => RK4
182 real(RKG) , intent(in) :: x, kappa
183 real(RKG) , intent(in) , optional :: tol
184 real(RKG) :: gammaIncLow
185 end function
186#endif
187
188#if RK3_ENABLED
189 impure elemental module function getGammaIncLow_RK3(x, kappa, tol) result(gammaIncLow)
190#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
191 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK3
192#endif
193 use pm_kind, only: RKG => RK3
194 real(RKG) , intent(in) :: x, kappa
195 real(RKG) , intent(in) , optional :: tol
196 real(RKG) :: gammaIncLow
197 end function
198#endif
199
200#if RK2_ENABLED
201 impure elemental module function getGammaIncLow_RK2(x, kappa, tol) result(gammaIncLow)
202#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
203 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK2
204#endif
205 use pm_kind, only: RKG => RK2
206 real(RKG) , intent(in) :: x, kappa
207 real(RKG) , intent(in) , optional :: tol
208 real(RKG) :: gammaIncLow
209 end function
210#endif
211
212#if RK1_ENABLED
213 impure elemental module function getGammaIncLow_RK1(x, kappa, tol) result(gammaIncLow)
214#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
215 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK1
216#endif
217 use pm_kind, only: RKG => RK1
218 real(RKG) , intent(in) :: x, kappa
219 real(RKG) , intent(in) , optional :: tol
220 real(RKG) :: gammaIncLow
221 end function
222#endif
223
224 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
225
226 end interface
227
228!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
229
312
313 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
314
315#if RK5_ENABLED
316 impure elemental module function getGammaIncUpp_RK5(x, kappa, tol) result(gammaIncUpp)
317#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
318 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK5
319#endif
320 use pm_kind, only: RKG => RK5
321 real(RKG) , intent(in) :: x, kappa
322 real(RKG) , intent(in) , optional :: tol
323 real(RKG) :: gammaIncUpp
324 end function
325#endif
326
327#if RK4_ENABLED
328 impure elemental module function getGammaIncUpp_RK4(x, kappa, tol) result(gammaIncUpp)
329#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
330 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK4
331#endif
332 use pm_kind, only: RKG => RK4
333 real(RKG) , intent(in) :: x, kappa
334 real(RKG) , intent(in) , optional :: tol
335 real(RKG) :: gammaIncUpp
336 end function
337#endif
338
339#if RK3_ENABLED
340 impure elemental module function getGammaIncUpp_RK3(x, kappa, tol) result(gammaIncUpp)
341#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
342 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK3
343#endif
344 use pm_kind, only: RKG => RK3
345 real(RKG) , intent(in) :: x, kappa
346 real(RKG) , intent(in) , optional :: tol
347 real(RKG) :: gammaIncUpp
348 end function
349#endif
350
351#if RK2_ENABLED
352 impure elemental module function getGammaIncUpp_RK2(x, kappa, tol) result(gammaIncUpp)
353#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
354 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK2
355#endif
356 use pm_kind, only: RKG => RK2
357 real(RKG) , intent(in) :: x, kappa
358 real(RKG) , intent(in) , optional :: tol
359 real(RKG) :: gammaIncUpp
360 end function
361#endif
362
363#if RK1_ENABLED
364 impure elemental module function getGammaIncUpp_RK1(x, kappa, tol) result(gammaIncUpp)
365#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
366 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK1
367#endif
368 use pm_kind, only: RKG => RK1
369 real(RKG) , intent(in) :: x, kappa
370 real(RKG) , intent(in) , optional :: tol
371 real(RKG) :: gammaIncUpp
372 end function
373#endif
374
375 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
376
377 end interface
378
379!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
380
471
472 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
473
474#if RK5_ENABLED
475 PURE elemental module subroutine setGammaIncLowDef_RK5(gammaIncLow, x, logGammaKappa, kappa, info, tol)
476#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
477 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowDef_RK5
478#endif
479 use pm_kind, only: RKG => RK5
480 real(RKG) , intent(out) :: gammaIncLow
481 real(RKG) , intent(in) :: x, logGammaKappa, kappa
482 integer(IK) , intent(out) :: info
483 real(RKG) , intent(in) , optional :: tol
484 end subroutine
485#endif
486
487#if RK4_ENABLED
488 PURE elemental module subroutine setGammaIncLowDef_RK4(gammaIncLow, x, logGammaKappa, kappa, info, tol)
489#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
490 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowDef_RK4
491#endif
492 use pm_kind, only: RKG => RK4
493 real(RKG) , intent(out) :: gammaIncLow
494 real(RKG) , intent(in) :: x, logGammaKappa, kappa
495 integer(IK) , intent(out) :: info
496 real(RKG) , intent(in) , optional :: tol
497 end subroutine
498#endif
499
500#if RK3_ENABLED
501 PURE elemental module subroutine setGammaIncLowDef_RK3(gammaIncLow, x, logGammaKappa, kappa, info, tol)
502#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
503 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowDef_RK3
504#endif
505 use pm_kind, only: RKG => RK3
506 real(RKG) , intent(out) :: gammaIncLow
507 real(RKG) , intent(in) :: x, logGammaKappa, kappa
508 integer(IK) , intent(out) :: info
509 real(RKG) , intent(in) , optional :: tol
510 end subroutine
511#endif
512
513#if RK2_ENABLED
514 PURE elemental module subroutine setGammaIncLowDef_RK2(gammaIncLow, x, logGammaKappa, kappa, info, tol)
515#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
516 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowDef_RK2
517#endif
518 use pm_kind, only: RKG => RK2
519 real(RKG) , intent(out) :: gammaIncLow
520 real(RKG) , intent(in) :: x, logGammaKappa, kappa
521 integer(IK) , intent(out) :: info
522 real(RKG) , intent(in) , optional :: tol
523 end subroutine
524#endif
525
526#if RK1_ENABLED
527 PURE elemental module subroutine setGammaIncLowDef_RK1(gammaIncLow, x, logGammaKappa, kappa, info, tol)
528#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
529 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowDef_RK1
530#endif
531 use pm_kind, only: RKG => RK1
532 real(RKG) , intent(out) :: gammaIncLow
533 real(RKG) , intent(in) :: x, logGammaKappa, kappa
534 integer(IK) , intent(out) :: info
535 real(RKG) , intent(in) , optional :: tol
536 end subroutine
537#endif
538
539 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
540
541 end interface
542
543!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
544
638
639 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
640
641#if RK5_ENABLED
642 PURE elemental module subroutine setGammaIncUppDef_RK5(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
643#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
644 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppDef_RK5
645#endif
646 use pm_kind, only: RKG => RK5
647 real(RKG) , intent(out) :: gammaIncUpp
648 real(RKG) , intent(in) :: x, logGammaKappa, kappa
649 integer(IK) , intent(out) :: info
650 real(RKG) , intent(in) , optional :: tol
651 end subroutine
652#endif
653
654#if RK4_ENABLED
655 PURE elemental module subroutine setGammaIncUppDef_RK4(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
656#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
657 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppDef_RK4
658#endif
659 use pm_kind, only: RKG => RK4
660 real(RKG) , intent(out) :: gammaIncUpp
661 real(RKG) , intent(in) :: x, logGammaKappa, kappa
662 integer(IK) , intent(out) :: info
663 real(RKG) , intent(in) , optional :: tol
664 end subroutine
665#endif
666
667#if RK3_ENABLED
668 PURE elemental module subroutine setGammaIncUppDef_RK3(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
669#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
670 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppDef_RK3
671#endif
672 use pm_kind, only: RKG => RK3
673 real(RKG) , intent(out) :: gammaIncUpp
674 real(RKG) , intent(in) :: x, logGammaKappa, kappa
675 integer(IK) , intent(out) :: info
676 real(RKG) , intent(in) , optional :: tol
677 end subroutine
678#endif
679
680#if RK2_ENABLED
681 PURE elemental module subroutine setGammaIncUppDef_RK2(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
682#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
683 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppDef_RK2
684#endif
685 use pm_kind, only: RKG => RK2
686 real(RKG) , intent(out) :: gammaIncUpp
687 real(RKG) , intent(in) :: x, logGammaKappa, kappa
688 integer(IK) , intent(out) :: info
689 real(RKG) , intent(in) , optional :: tol
690 end subroutine
691#endif
692
693#if RK1_ENABLED
694 PURE elemental module subroutine setGammaIncUppDef_RK1(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
695#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
696 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppDef_RK1
697#endif
698 use pm_kind, only: RKG => RK1
699 real(RKG) , intent(out) :: gammaIncUpp
700 real(RKG) , intent(in) :: x, logGammaKappa, kappa
701 integer(IK) , intent(out) :: info
702 real(RKG) , intent(in) , optional :: tol
703 end subroutine
704#endif
705
706 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
707
708 end interface
709
710!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
711
796
797#if RK5_ENABLED
798 PURE elemental module subroutine setGammaIncLowSeries_RK5(gammaIncLow, x, logGammaKappa, kappa, info, tol)
799#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
800 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowSeries_RK5
801#endif
802 use pm_kind, only: RKG => RK5
803 real(RKG) , intent(out) :: gammaIncLow
804 real(RKG) , intent(in) :: x, logGammaKappa, kappa
805 integer(IK) , intent(out) :: info
806 real(RKG) , intent(in) , optional :: tol
807 end subroutine
808#endif
809
810#if RK4_ENABLED
811 PURE elemental module subroutine setGammaIncLowSeries_RK4(gammaIncLow, x, logGammaKappa, kappa, info, tol)
812#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
813 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowSeries_RK4
814#endif
815 use pm_kind, only: RKG => RK4
816 real(RKG) , intent(out) :: gammaIncLow
817 real(RKG) , intent(in) :: x, logGammaKappa, kappa
818 integer(IK) , intent(out) :: info
819 real(RKG) , intent(in) , optional :: tol
820 end subroutine
821#endif
822
823#if RK3_ENABLED
824 PURE elemental module subroutine setGammaIncLowSeries_RK3(gammaIncLow, x, logGammaKappa, kappa, info, tol)
825#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
826 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowSeries_RK3
827#endif
828 use pm_kind, only: RKG => RK3
829 real(RKG) , intent(out) :: gammaIncLow
830 real(RKG) , intent(in) :: x, logGammaKappa, kappa
831 integer(IK) , intent(out) :: info
832 real(RKG) , intent(in) , optional :: tol
833 end subroutine
834#endif
835
836#if RK2_ENABLED
837 PURE elemental module subroutine setGammaIncLowSeries_RK2(gammaIncLow, x, logGammaKappa, kappa, info, tol)
838#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
839 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowSeries_RK2
840#endif
841 use pm_kind, only: RKG => RK2
842 real(RKG) , intent(out) :: gammaIncLow
843 real(RKG) , intent(in) :: x, logGammaKappa, kappa
844 integer(IK) , intent(out) :: info
845 real(RKG) , intent(in) , optional :: tol
846 end subroutine
847#endif
848
849#if RK1_ENABLED
850 PURE elemental module subroutine setGammaIncLowSeries_RK1(gammaIncLow, x, logGammaKappa, kappa, info, tol)
851#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
852 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncLowSeries_RK1
853#endif
854 use pm_kind, only: RKG => RK1
855 real(RKG) , intent(out) :: gammaIncLow
856 real(RKG) , intent(in) :: x, logGammaKappa, kappa
857 integer(IK) , intent(out) :: info
858 real(RKG) , intent(in) , optional :: tol
859 end subroutine
860#endif
861
862 end interface
863
864!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
865
952
953#if RK5_ENABLED
954 PURE elemental module subroutine setGammaIncUppContFrac_RK5(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
955#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
956 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppContFrac_RK5
957#endif
958 use pm_kind, only: RKG => RK5
959 real(RKG) , intent(out) :: gammaIncUpp
960 real(RKG) , intent(in) :: x, logGammaKappa, kappa
961 integer(IK) , intent(out) :: info
962 real(RKG) , intent(in) , optional :: tol
963 end subroutine
964#endif
965
966#if RK4_ENABLED
967 PURE elemental module subroutine setGammaIncUppContFrac_RK4(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
968#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
969 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppContFrac_RK4
970#endif
971 use pm_kind, only: RKG => RK4
972 real(RKG) , intent(out) :: gammaIncUpp
973 real(RKG) , intent(in) :: x, logGammaKappa, kappa
974 integer(IK) , intent(out) :: info
975 real(RKG) , intent(in) , optional :: tol
976 end subroutine
977#endif
978
979#if RK3_ENABLED
980 PURE elemental module subroutine setGammaIncUppContFrac_RK3(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
981#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
982 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppContFrac_RK3
983#endif
984 use pm_kind, only: RKG => RK3
985 real(RKG) , intent(out) :: gammaIncUpp
986 real(RKG) , intent(in) :: x, logGammaKappa, kappa
987 integer(IK) , intent(out) :: info
988 real(RKG) , intent(in) , optional :: tol
989 end subroutine
990#endif
991
992#if RK2_ENABLED
993 PURE elemental module subroutine setGammaIncUppContFrac_RK2(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
994#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
995 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppContFrac_RK2
996#endif
997 use pm_kind, only: RKG => RK2
998 real(RKG) , intent(out) :: gammaIncUpp
999 real(RKG) , intent(in) :: x, logGammaKappa, kappa
1000 integer(IK) , intent(out) :: info
1001 real(RKG) , intent(in) , optional :: tol
1002 end subroutine
1003#endif
1004
1005#if RK1_ENABLED
1006 PURE elemental module subroutine setGammaIncUppContFrac_RK1(gammaIncUpp, x, logGammaKappa, kappa, info, tol)
1007#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1008 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaIncUppContFrac_RK1
1009#endif
1010 use pm_kind, only: RKG => RK1
1011 real(RKG) , intent(out) :: gammaIncUpp
1012 real(RKG) , intent(in) :: x, logGammaKappa, kappa
1013 integer(IK) , intent(out) :: info
1014 real(RKG) , intent(in) , optional :: tol
1015 end subroutine
1016#endif
1017
1018 end interface
1019
1020!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1021
1022!contains
1023
1024!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1025
1026! !> \brief
1027! !> Return the Gamma function for a half-integer input as real of kind \RK.
1028! !>
1029! !> \param[in] positiveHalfInteger : The input half integer as a real number
1030! !>
1031! !> \return
1032! !> `gammaHalfInt` : The Gamma function for a half integer input.
1033! !>
1034! !> \remark
1035! !> The equation for half-integer Gamma-function is given as,
1036! !> \f{equation}{
1037! !> \Gamma \left( \frac{n}{2} \right) = \sqrt \pi \frac{ (n-2)!! }{ 2^\frac{n-1}{2} } ~,
1038! !> \f}
1039! pure function getGamHalfInt(positiveHalfInteger) result(gammaHalfInt)
1040!#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1041! !DEC$ ATTRIBUTES DLLEXPORT :: getGamHalfInt
1042!#endif
1043! use pm_kind, only: IK, RK, SQRT_PI
1044! implicit none
1045! real(RK), intent(in) :: positiveHalfInteger
1046! real(RKG) :: gammaHalfInt
1047! integer(IK) :: i,k
1048! gammaHalfInt = SQRT_PI
1049! k = nint(positiveHalfInteger-0.5_RK,kind=IK) ! positiveHalfInteger = k + 1/2
1050! do i = k+1, 2*k
1051! gammaHalfInt = gammaHalfInt * 0.25_RK * i
1052! end do
1053! end function getGamHalfInt
1054!
1055!!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1056!
1057! !> \brief
1058! !> Return the natural logarithm of the Gamma function for a half-integer input as real of kind \RK.
1059! !>
1060! !> \param[in] positiveHalfInteger : The input half integer as a real number
1061! !>
1062! !> \return
1063! !> `gammaHalfInt` : The Gamma function for a half integer input.
1064! !>
1065! !> \remark
1066! !> The equation for half-integer Gamma-function is given as,
1067! !> \f{equation}{
1068! !> \Gamma \left( \frac{n}{2} \right) = \sqrt \pi \frac{ (n-2)!! }{ 2^\frac{n-1}{2} } ~,
1069! !> \f}
1070! pure function getLogGammaHalfInt(positiveHalfInteger) result(logGammaHalfInt)
1071!#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1072! !DEC$ ATTRIBUTES DLLEXPORT :: getLogGammaHalfInt
1073!#endif
1074! use pm_kind, only: IK, RK, SQRT_PI
1075! implicit none
1076! real(RK), intent(in) :: positiveHalfInteger
1077! real(RK), parameter :: COEF = log(0.25_RK)
1078! real(RK), parameter :: LOG_SQRTPI = log(SQRT_PI)
1079! real(RKG) :: logGammaHalfInt
1080! integer(IK) :: i, k
1081! k = nint(positiveHalfInteger-0.5_RK,kind=IK) ! positiveHalfInteger = k + 1/2
1082! logGammaHalfInt = LOG_SQRTPI
1083! do i = k+1, 2*k
1084! logGammaHalfInt = logGammaHalfInt + COEF + log(real(i,kind=RK))
1085! end do
1086! end function getLogGammaHalfInt
1087
1088!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1089
1090end module pm_mathGamma ! LCOV_EXCL_LINE
Generate and return the regularized Lower Incomplete Gamma function for the specified shape parameter...
Generate and return the regularized Upper Incomplete Gamma function for the specified shape parameter...
Return the regularized Lower Incomplete Gamma function for the specified upper limit x and shape para...
Return the regularized Lower Incomplete Gamma function for the specified shape parameter ( ) and uppe...
Return the regularized Upper Incomplete Gamma function for the specified lower limit x and shape para...
Return the regularized Upper Incomplete Gamma function for the specified shape parameter ( ) and lowe...
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 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 module contains procedures and generic interfaces for the Lower and Upper Incomplete Gamma funct...
character(*, SK), parameter MODULE_NAME