Line data    Source code

       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             : 
      17             : !>  \brief
      18             : !>  This include file contains the implementation of procedures in [pm_distGenExpGamma](@ref pm_distGenExpGamma).
      19             : !>
      20             : !>  \finmain
      21             : !>
      22             : !>  \author
      23             : !>  \AmirShahmoradi, Oct 16, 2009, 12:20 PM, Michigan
      24             : 
      25             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      26             : 
      27             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      28             : #if     getGenExpGammaLogPDFNF_ENABLED && KD_ENABLED
      29             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      30             : 
      31       89361 :         CHECK_ASSERTION(__LINE__, kappa > 0._RKC, SK_"@getGenExpGammaLogPDFNF(): The condition `kappa > 0.` must hold. kappa = "//getStr(kappa)) ! fpp
      32       89361 :         logPDFNF = -log_gamma(kappa)
      33             : 
      34             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      35             : #elif   getGenExpGammaLogPDFNF_ENABLED && KO_ENABLED
      36             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      37             :         
      38       83638 :         CHECK_ASSERTION(__LINE__, invOmega > 0._RKC, SK_"@getGenExpGammaLogPDFNF(): The condition `invOmega > 0.` must hold. invOmega = "//getStr(invOmega)) ! fpp
      39       83638 :         logPDFNF = getGenExpGammaLogPDFNF(kappa)
      40       83638 :         if (invOmega /= 1._RKC) logPDFNF = logPDFNF + log(invOmega)
      41             : 
      42             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%
      43             : #elif   getGenExpGammaLogPDF_ENABLED
      44             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%
      45             : 
      46             :         real(RKC) :: kappa_def, invOmega_def, logSigma_def
      47        6406 :         kappa_def = 1._RKC; if (present(kappa)) kappa_def = kappa
      48        6406 :         invOmega_def = 1._RKC; if (present(invOmega)) invOmega_def = invOmega
      49        6406 :         logSigma_def = 0._RKC; if (present(logSigma)) logSigma_def = logSigma
      50        6406 :         call setGenExpGammaLogPDF(logPDF, x, getGenExpGammaLogPDFNF(kappa_def, invOmega_def), kappa_def, invOmega_def, logSigma_def)
      51             : 
      52             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      53             : #elif   setGenExpGammaLogPDF_ENABLED && DDDD_ENABLED
      54             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      55             : 
      56        1401 :         logPDF = x - exp(x)
      57             : 
      58             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      59             : #elif   setGenExpGammaLogPDF_ENABLED && NKDD_ENABLED
      60             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      61             : 
      62        1579 :         CHECK_ASSERTION(__LINE__, kappa > 0._RKC, SK_"@getGenExpGammaLogPDFNF(): The condition `kappa > 0.` must hold. kappa = "//getStr(kappa)) ! fpp
      63        4737 :         CHECK_ASSERTION(__LINE__, abs(getGenExpGammaLogPDFNF(kappa) - logPDFNF) <= 100 * epsilon(0._RKC), \
      64             :         SK_"@setGenExpGammaLogPDF(): The condition `abs(getGenExpGammaLogPDFNF(kappa) - logPDFNF) <= 100 * epsilon(0._RKC)` must hold. getGenExpGammaLogPDFNF(kappa), logPDFNF = " \
      65             :         //getStr([getGenExpGammaLogPDFNF(kappa), logPDFNF])) ! fpp
      66        1579 :         logPDF = logPDFNF + kappa * x - exp(x)
      67             : 
      68             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      69             : #elif   setGenExpGammaLogPDF_ENABLED && NKOD_ENABLED
      70             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      71             : 
      72             :         real(RKC) :: xscaled
      73             :         real(RKC), parameter :: LOG_SQRT_HUGE = log(sqrt(huge(0._RKC)))
      74       22292 :         xscaled = x * invOmega
      75       22292 :         CHECK_ASSERTION(__LINE__, kappa > 0._RKC, SK_"@getGenExpGammaLogPDFNF(): The condition `kappa > 0.` must hold. kappa = "//getStr(kappa)) ! fpp
      76       22292 :         CHECK_ASSERTION(__LINE__, invOmega > 0._RKC, SK_"@setGenExpGammaLogPDF(): The condition `invOmega > 0.` must hold. invOmega = "//getStr(invOmega)) ! fpp
      77       66876 :         CHECK_ASSERTION(__LINE__, abs(getGenExpGammaLogPDFNF(kappa, invOmega) - logPDFNF) <= 100 * epsilon(0._RKC), \
      78             :         SK_"@setGenExpGammaLogPDF(): The condition `abs(getGenExpGammaLogPDFNF(kappa, invOmega) - logPDFNF) <= 100 * epsilon(0._RKC)` must hold. getGenExpGammaLogPDFNF(kappa, invOmega), logPDFNF = " \
      79             :         //getStr([getGenExpGammaLogPDFNF(kappa, invOmega), logPDFNF])) ! fpp
      80       22292 :         if (xscaled < LOG_SQRT_HUGE) then
      81       22284 :             logPDF = logPDFNF + kappa * xscaled - exp(xscaled)
      82             :         else
      83           8 :             logPDF = -LOG_SQRT_HUGE
      84             :         end if
      85             : 
      86             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      87             : #elif   setGenExpGammaLogPDF_ENABLED && NKOS_ENABLED
      88             :         !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      89             : 
      90       43623 :         CHECK_ASSERTION(__LINE__, abs(getGenExpGammaLogPDFNF(kappa, invOmega) - logPDFNF) <= 100 * epsilon(0._RKC), \
      91             :         SK_"@setGenExpGammaLogPDF(): The condition `abs(getGenExpGammaLogPDFNF(kappa, invOmega) - logPDFNF) <= 100 * epsilon(0._RKC)` must hold. getGenExpGammaLogPDFNF(kappa, invOmega), logPDFNF = " \
      92             :         //getStr([getGenExpGammaLogPDFNF(kappa, invOmega), logPDFNF])) ! fpp
      93       14541 :         call setGenExpGammaLogPDF(logPDF, x - logSigma, logPDFNF, kappa, invOmega)
      94             : 
      95             :         !%%%%%%%%%%%%%%%%%%%%%%%%
      96             : #elif   getGenExpGammaCDF_ENABLED
      97             :         !%%%%%%%%%%%%%%%%%%%%%%%%
      98             : 
      99             :         integer(IK) :: info
     100             :         real(RKC) :: xnormed
     101        6046 :         if (present(logSigma)) then
     102           9 :             xnormed = x - logSigma
     103             :         else
     104        6037 :             xnormed = x
     105             :         end if
     106        6046 :         if (present(invOmega)) xnormed = xnormed * invOmega
     107        6046 :         if (present(kappa)) then
     108        6037 :             call setGenExpGammaCDF(cdf, xnormed, log_gamma(kappa), kappa, info)
     109             :         else
     110           9 :             call setGenExpGammaCDF(cdf, xnormed, info)
     111             :         end if
     112        6046 :         if (info < 0_IK) error stop MODULE_NAME//SK_"@getGenExpGammaCDF(): The computation of the regularized Lower Incomplete Gamma function failed. This can happen if `kappa` is too large."
     113             : 
     114             :         !%%%%%%%%%%%%%%%%%%%%%%%%
     115             : #elif   setGenExpGammaCDF_ENABLED
     116             :         !%%%%%%%%%%%%%%%%%%%%%%%%
     117             : 
     118             : #if     DDD_ENABLED
     119             :         real(RKC), parameter :: kappa = 1._RKC, logGammaKappa = log_gamma(kappa)
     120          10 :         call setGammaIncLow(cdf, exp(x), logGammaKappa, kappa, info)
     121             : #else
     122       12050 :         CHECK_ASSERTION(__LINE__, kappa > 0._RKC, SK_"@setGenExpGammaCDF(): The condition `kappa > 0.` must hold. kappa = "//getStr(kappa)) ! fpp
     123       36150 :         CHECK_ASSERTION(__LINE__, abs(log_gamma(kappa) - logGammaKappa) < 100 * epsilon(0._RKC), SK_"@setGenExpGammaCDF(): The condition `abs(log_gamma(kappa) - logGammaKappa) < 100 * epsilon(0._RKC)` must hold. log_gamma(kappa), logGammaKappa = "//getStr([log_gamma(kappa), logGammaKappa])) ! fpp
     124             : #if     KDD_ENABLED
     125        6048 :         call setGammaIncLow(cdf, exp(x), logGammaKappa, kappa, info)
     126             : #else
     127        6002 :         CHECK_ASSERTION(__LINE__, invOmega > 0._RKC, SK_"@setGenExpGammaCDF(): The condition `invOmega > 0.` must hold. invOmega = "//getStr(invOmega)) ! fpp
     128             : #if     KOD_ENABLED
     129        6001 :         call setGammaIncLow(cdf, exp(x * invOmega), logGammaKappa, kappa, info)
     130             : #elif   KOS_ENABLED
     131           1 :         call setGammaIncLow(cdf, exp((x - logSigma) * invOmega), logGammaKappa, kappa, info)
     132             : #else
     133             : #error  "Unrecognized interface."
     134             : #endif
     135             : #endif
     136             : #endif
     137             : 
     138             : #else
     139             :         !%%%%%%%%%%%%%%%%%%%%%%%%
     140             : #error  "Unrecognized interface."
     141             :         !%%%%%%%%%%%%%%%%%%%%%%%%
     142             : #endif