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
36
38
39 use pm_kind, only: SK, IK, LK
40
41 implicit none
42
43 character(*, SK), parameter :: MODULE_NAME = "@pm_mathGamma"
44
45!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
46
125
126 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
127
128#if RK5_ENABLED
129 impure elemental module function getGammaIncLow_RK5(x, kappa) result(gammaIncLow)
130#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
131 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK5
132#endif
133 use pm_kind, only: RKG => RK5
134 real(RKG) , intent(in) :: x, kappa
135 real(RKG) :: gammaIncLow
136 end function
137#endif
138
139#if RK4_ENABLED
140 impure elemental module function getGammaIncLow_RK4(x, kappa) result(gammaIncLow)
141#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
142 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK4
143#endif
144 use pm_kind, only: RKG => RK4
145 real(RKG) , intent(in) :: x, kappa
146 real(RKG) :: gammaIncLow
147 end function
148#endif
149
150#if RK3_ENABLED
151 impure elemental module function getGammaIncLow_RK3(x, kappa) result(gammaIncLow)
152#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
153 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK3
154#endif
155 use pm_kind, only: RKG => RK3
156 real(RKG) , intent(in) :: x, kappa
157 real(RKG) :: gammaIncLow
158 end function
159#endif
160
161#if RK2_ENABLED
162 impure elemental module function getGammaIncLow_RK2(x, kappa) result(gammaIncLow)
163#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
164 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK2
165#endif
166 use pm_kind, only: RKG => RK2
167 real(RKG) , intent(in) :: x, kappa
168 real(RKG) :: gammaIncLow
169 end function
170#endif
171
172#if RK1_ENABLED
173 impure elemental module function getGammaIncLow_RK1(x, kappa) result(gammaIncLow)
174#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
175 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncLow_RK1
176#endif
177 use pm_kind, only: RKG => RK1
178 real(RKG) , intent(in) :: x, kappa
179 real(RKG) :: gammaIncLow
180 end function
181#endif
182
183 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
184
185 end interface
186
187!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
188
265
266 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
267
268#if RK5_ENABLED
269 impure elemental module function getGammaIncUpp_RK5(x, kappa) result(gammaIncUpp)
270#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
271 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK5
272#endif
273 use pm_kind, only: RKG => RK5
274 real(RKG) , intent(in) :: x, kappa
275 real(RKG) :: gammaIncUpp
276 end function
277#endif
278
279#if RK4_ENABLED
280 impure elemental module function getGammaIncUpp_RK4(x, kappa) result(gammaIncUpp)
281#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
282 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK4
283#endif
284 use pm_kind, only: RKG => RK4
285 real(RKG) , intent(in) :: x, kappa
286 real(RKG) :: gammaIncUpp
287 end function
288#endif
289
290#if RK3_ENABLED
291 impure elemental module function getGammaIncUpp_RK3(x, kappa) result(gammaIncUpp)
292#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
293 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK3
294#endif
295 use pm_kind, only: RKG => RK3
296 real(RKG) , intent(in) :: x, kappa
297 real(RKG) :: gammaIncUpp
298 end function
299#endif
300
301#if RK2_ENABLED
302 impure elemental module function getGammaIncUpp_RK2(x, kappa) result(gammaIncUpp)
303#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
304 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK2
305#endif
306 use pm_kind, only: RKG => RK2
307 real(RKG) , intent(in) :: x, kappa
308 real(RKG) :: gammaIncUpp
309 end function
310#endif
311
312#if RK1_ENABLED
313 impure elemental module function getGammaIncUpp_RK1(x, kappa) result(gammaIncUpp)
314#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
315 !DEC$ ATTRIBUTES DLLEXPORT :: getGammaIncUpp_RK1
316#endif
317 use pm_kind, only: RKG => RK1
318 real(RKG) , intent(in) :: x, kappa
319 real(RKG) :: gammaIncUpp
320 end function
321#endif
322
323 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
324
325 end interface
326
327!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
328
416 interface setGammaInc
417
418 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
419
420#if RK5_ENABLED
421 PURE elemental module subroutine setGammaInc_RK5(gammaIncLow, gammaIncUpp, x, kappa, info)
422#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
423 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaInc_RK5
424#endif
425 use pm_kind, only: RKG => RK5
426 real(RKG) , intent(out) :: gammaIncLow, gammaIncUpp
427 real(RKG) , intent(in) :: x, kappa
428 integer(IK) , intent(out) :: info
429 end subroutine
430#endif
431
432#if RK4_ENABLED
433 PURE elemental module subroutine setGammaInc_RK4(gammaIncLow, gammaIncUpp, x, kappa, info)
434#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
435 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaInc_RK4
436#endif
437 use pm_kind, only: RKG => RK4
438 real(RKG) , intent(out) :: gammaIncLow, gammaIncUpp
439 real(RKG) , intent(in) :: x, kappa
440 integer(IK) , intent(out) :: info
441 end subroutine
442#endif
443
444#if RK3_ENABLED
445 PURE elemental module subroutine setGammaInc_RK3(gammaIncLow, gammaIncUpp, x, kappa, info)
446#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
447 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaInc_RK3
448#endif
449 use pm_kind, only: RKG => RK3
450 real(RKG) , intent(out) :: gammaIncLow, gammaIncUpp
451 real(RKG) , intent(in) :: x, kappa
452 integer(IK) , intent(out) :: info
453 end subroutine
454#endif
455
456#if RK2_ENABLED
457 PURE elemental module subroutine setGammaInc_RK2(gammaIncLow, gammaIncUpp, x, kappa, info)
458#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
459 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaInc_RK2
460#endif
461 use pm_kind, only: RKG => RK2
462 real(RKG) , intent(out) :: gammaIncLow, gammaIncUpp
463 real(RKG) , intent(in) :: x, kappa
464 integer(IK) , intent(out) :: info
465 end subroutine
466#endif
467
468#if RK1_ENABLED
469 PURE elemental module subroutine setGammaInc_RK1(gammaIncLow, gammaIncUpp, x, kappa, info)
470#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
471 !DEC$ ATTRIBUTES DLLEXPORT :: setGammaInc_RK1
472#endif
473 use pm_kind, only: RKG => RK1
474 real(RKG) , intent(out) :: gammaIncLow, gammaIncUpp
475 real(RKG) , intent(in) :: x, kappa
476 integer(IK) , intent(out) :: info
477 end subroutine
478#endif
479
480 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
481
482 end interface
483
484!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
485
486!contains
487
488!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
489
490! !> \brief
491! !> Return the Gamma function for a half-integer input as real of kind \RK.
492! !>
493! !> \param[in] positiveHalfInteger : The input half integer as a real number
494! !>
495! !> \return
496! !> `gammaHalfInt` : The Gamma function for a half integer input.
497! !>
498! !> \remark
499! !> The equation for half-integer Gamma-function is given as,
500! !> \f{equation}{
501! !> \Gamma \left( \frac{n}{2} \right) = \sqrt \pi \frac{ (n-2)!! }{ 2^\frac{n-1}{2} } ~,
502! !> \f}
503! pure function getGamHalfInt(positiveHalfInteger) result(gammaHalfInt)
504!#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
505! !DEC$ ATTRIBUTES DLLEXPORT :: getGamHalfInt
506!#endif
507! use pm_kind, only: IK, RK, SQRT_PI
508! implicit none
509! real(RK), intent(in) :: positiveHalfInteger
510! real(RKG) :: gammaHalfInt
511! integer(IK) :: i,k
512! gammaHalfInt = SQRT_PI
513! k = nint(positiveHalfInteger-0.5_RK,kind=IK) ! positiveHalfInteger = k + 1/2
514! do i = k+1, 2*k
515! gammaHalfInt = gammaHalfInt * 0.25_RK * i
516! end do
517! end function getGamHalfInt
518!
519!!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
520!
521! !> \brief
522! !> Return the natural logarithm of the Gamma function for a half-integer input as real of kind \RK.
523! !>
524! !> \param[in] positiveHalfInteger : The input half integer as a real number
525! !>
526! !> \return
527! !> `gammaHalfInt` : The Gamma function for a half integer input.
528! !>
529! !> \remark
530! !> The equation for half-integer Gamma-function is given as,
531! !> \f{equation}{
532! !> \Gamma \left( \frac{n}{2} \right) = \sqrt \pi \frac{ (n-2)!! }{ 2^\frac{n-1}{2} } ~,
533! !> \f}
534! pure function getLogGammaHalfInt(positiveHalfInteger) result(logGammaHalfInt)
535!#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
536! !DEC$ ATTRIBUTES DLLEXPORT :: getLogGammaHalfInt
537!#endif
538! use pm_kind, only: IK, RK, SQRT_PI
539! implicit none
540! real(RK), intent(in) :: positiveHalfInteger
541! real(RK), parameter :: COEF = log(0.25_RK)
542! real(RK), parameter :: LOG_SQRTPI = log(SQRT_PI)
543! real(RKG) :: logGammaHalfInt
544! integer(IK) :: i, k
545! k = nint(positiveHalfInteger-0.5_RK,kind=IK) ! positiveHalfInteger = k + 1/2
546! logGammaHalfInt = LOG_SQRTPI
547! do i = k+1, 2*k
548! logGammaHalfInt = logGammaHalfInt + COEF + log(real(i,kind=RK))
549! end do
550! end function getLogGammaHalfInt
551
552!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
553
554end 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 shape parameter ( ) and uppe...
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