Line data Source code
1 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3 : !!!!
4 : !!!! MIT License
5 : !!!!
6 : !!!! ParaMonte: plain powerful parallel Monte Carlo library.
7 : !!!!
8 : !!!! Copyright (C) 2012-present, The Computational Data Science Lab
9 : !!!!
10 : !!!! This file is part of the ParaMonte library.
11 : !!!!
12 : !!!! Permission is hereby granted, free of charge, to any person obtaining a
13 : !!!! copy of this software and associated documentation files (the "Software"),
14 : !!!! to deal in the Software without restriction, including without limitation
15 : !!!! the rights to use, copy, modify, merge, publish, distribute, sublicense,
16 : !!!! and/or sell copies of the Software, and to permit persons to whom the
17 : !!!! Software is furnished to do so, subject to the following conditions:
18 : !!!!
19 : !!!! The above copyright notice and this permission notice shall be
20 : !!!! included in all copies or substantial portions of the Software.
21 : !!!!
22 : !!!! THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
23 : !!!! EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
24 : !!!! MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
25 : !!!! IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,
26 : !!!! DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR
27 : !!!! OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE
28 : !!!! OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
29 : !!!!
30 : !!!! ACKNOWLEDGMENT
31 : !!!!
32 : !!!! ParaMonte is an honor-ware and its currency is acknowledgment and citations.
33 : !!!! As per the ParaMonte library license agreement terms, if you use any parts of
34 : !!!! this library for any purposes, kindly acknowledge the use of ParaMonte in your
35 : !!!! work (education/research/industry/development/...) by citing the ParaMonte
36 : !!!! library as described on this page:
37 : !!!!
38 : !!!! https://github.com/cdslaborg/paramonte/blob/main/ACKNOWLEDGMENT.md
39 : !!!!
40 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
41 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
42 :
43 : !> \brief This module contains tests of the module [Math_mod](@ref math_mod).
44 : !> \author Amir Shahmoradi
45 :
46 : module Test_Math_mod
47 :
48 : use Math_mod
49 : use Err_mod, only: Err_type
50 : use Test_mod, only: Test_type
51 : implicit none
52 :
53 : private
54 : public :: test_Math
55 :
56 : type(Test_type) :: Test
57 :
58 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
59 :
60 : contains
61 :
62 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
63 :
64 1 : subroutine test_Math()
65 : implicit none
66 1 : Test = Test_type(moduleName=MODULE_NAME)
67 1 : call Test%run(test_getCumSum_IK_1, "test_getCumSum_IK_1")
68 1 : call Test%run(test_getCumSum_RK_1, "test_getCumSum_RK_1")
69 1 : call Test%run(test_getFactorial_1, "test_getFactorial_1")
70 1 : call Test%run(test_getDistanceSq_1, "test_getDistanceSq_1")
71 1 : call Test%run(test_getEllVolCoef_1, "test_getEllVolCoef_1")
72 1 : call Test%run(test_getEllVolCoef_2, "test_getEllVolCoef_2")
73 1 : call Test%run(test_getLowerGamma_1, "test_getLowerGamma_1")
74 1 : call Test%run(test_getLowerGamma_2, "test_getLowerGamma_2")
75 1 : call Test%run(test_getLowerGamma_3, "test_getLowerGamma_3")
76 1 : call Test%run(test_getLowerGamma_4, "test_getLowerGamma_4")
77 1 : call Test%run(test_getLowerGamma_5, "test_getLowerGamma_5")
78 1 : call Test%run(test_getUpperGamma_1, "test_getUpperGamma_1")
79 1 : call Test%run(test_getUpperGamma_2, "test_getUpperGamma_2")
80 1 : call Test%run(test_getUpperGamma_3, "test_getUpperGamma_3")
81 1 : call Test%run(test_getUpperGamma_4, "test_getUpperGamma_4")
82 1 : call Test%run(test_getUpperGamma_5, "test_getUpperGamma_5")
83 1 : call Test%run(test_getGammaSeries_1, "test_getGammaSeries_1")
84 1 : call Test%run(test_getLogSubExp_RK_1, "test_getLogSubExp_RK_1")
85 1 : call Test%run(test_getLogSumExp_RK_1, "test_getLogSumExp_RK_1")
86 1 : call Test%run(test_getLogSumExp_RK_2, "test_getLogSumExp_RK_2")
87 1 : call Test%run(test_getLogSumExp_CK_1, "test_getLogSumExp_CK_1")
88 1 : call Test%run(test_getLogSumExp_CK_2, "test_getLogSumExp_CK_2")
89 1 : call Test%run(test_getLogFactorial_1, "test_getLogFactorial_1")
90 1 : call Test%run(test_getGammaHalfInt_1, "test_getGammaHalfInt_1")
91 1 : call Test%run(test_getGammaContFrac_1, "test_getGammaContFrac_1")
92 1 : call Test%run(test_getLogEllVolCoef_1, "test_getLogEllVolCoef_1")
93 1 : call Test%run(test_getLogEllVolCoef_2, "test_getLogEllVolCoef_2")
94 1 : call Test%run(test_getLogEggBoxSD_RK_1, "test_getLogEggBoxSD_RK_1")
95 1 : call Test%run(test_getLogEggBoxSD_CK_1, "test_getLogEggBoxSD_CK_1")
96 1 : call Test%run(test_getLogEggBoxMD_RK_1, "test_getLogEggBoxMD_RK_1")
97 1 : call Test%run(test_getLogEggBoxMD_CK_1, "test_getLogEggBoxMD_CK_1")
98 1 : call Test%run(test_getLogVolUnitBall_1, "test_getLogVolUnitBall_1")
99 1 : call Test%run(test_getLogVolUnitBall_2, "test_getLogVolUnitBall_2")
100 1 : call Test%run(test_getLogGammaHalfInt_1, "test_getLogGammaHalfInt_1")
101 1 : call Test%run(test_getLogVolEllipsoid_1, "test_getLogVolEllipsoid_1")
102 1 : call Test%run(test_getLogVolEllipsoids_1, "test_getLogVolEllipsoids_1")
103 1 : call Test%run(test_getCumSumReverse_IK_1, "test_getCumSumReverse_IK_1")
104 1 : call Test%run(test_getCumSumReverse_RK_1, "test_getCumSumReverse_RK_1")
105 1 : call Test%run(test_getCorCeofFromFisherTrans_1, "test_getCorCeofFromFisherTrans_1")
106 1 : call Test%run(test_getFisherTransFromcorCoef_1, "test_getFisherTransFromcorCoef_1")
107 1 : call Test%finalize()
108 1 : end subroutine test_Math
109 :
110 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
111 :
112 1 : function test_getDistanceSq_1() result(assertion)
113 1 : use Constants_mod, only: RK, IK
114 : implicit none
115 : logical :: assertion
116 : real(RK), parameter :: Point1(*) = [0._RK, 1._RK, 2._RK, 3._RK, 4._RK]
117 : real(RK), parameter :: Point2(*) = Point1 + 1._RK
118 : real(RK), parameter :: tolerance = 1.e-10_RK
119 : real(RK), parameter :: distanceSq_ref = norm2(Point2-Point1)**2
120 1 : real(RK) :: distanceSq
121 1 : real(RK) :: difference
122 1 : distanceSq = getDistanceSq(nd = size(Point1), Point1 = Point1, Point2 = Point2)
123 1 : difference = abs(distanceSq - distanceSq_ref) / distanceSq_ref
124 1 : assertion = difference < tolerance
125 1 : if (Test%isDebugMode .and. .not. assertion) then
126 : ! LCOV_EXCL_START
127 : write(Test%outputUnit,"(*(g0,:,' '))")
128 : write(Test%outputUnit,"(*(g0,:,' '))") "distanceSq_ref = ", distanceSq_ref
129 : write(Test%outputUnit,"(*(g0,:,' '))") "distanceSq = ", distanceSq
130 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
131 : write(Test%outputUnit,"(*(g0,:,' '))")
132 : end if
133 : ! LCOV_EXCL_STOP
134 1 : end function test_getDistanceSq_1
135 :
136 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
137 :
138 1 : function test_getCorCeofFromFisherTrans_1() result(assertion)
139 1 : use Constants_mod, only: RK, IK
140 : implicit none
141 : logical :: assertion
142 : real(RK), parameter :: fisherTrans = 1.5_RK
143 : real(RK), parameter :: corCoef_ref = 0.905148253644866_RK ! tanh(fisherTrans)
144 : real(RK), parameter :: tolerance = 1.e-10_RK
145 1 : real(RK) :: corCoef
146 1 : real(RK) :: difference
147 1 : corCoef = getCorCeofFromFisherTrans(fisherTrans = fisherTrans)
148 1 : difference = abs(corCoef - corCoef_ref) / corCoef_ref
149 1 : assertion = difference < tolerance
150 1 : if (Test%isDebugMode .and. .not. assertion) then
151 : ! LCOV_EXCL_START
152 : write(Test%outputUnit,"(*(g0,:,' '))")
153 : write(Test%outputUnit,"(*(g0,:,' '))") "corCoef_ref = ", corCoef_ref
154 : write(Test%outputUnit,"(*(g0,:,' '))") "corCoef = ", corCoef
155 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
156 : write(Test%outputUnit,"(*(g0,:,' '))")
157 : end if
158 : ! LCOV_EXCL_STOP
159 1 : end function test_getCorCeofFromFisherTrans_1
160 :
161 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
162 :
163 1 : function test_getFisherTransFromcorCoef_1() result(assertion)
164 1 : use Constants_mod, only: RK, IK
165 : implicit none
166 : logical :: assertion
167 : real(RK), parameter :: corCoef = 0.905148253644866_RK ! tanh(fisherTrans)
168 : real(RK), parameter :: fisherTrans_ref = 1.5_RK
169 : real(RK), parameter :: tolerance = 1.e-10_RK
170 1 : real(RK) :: fisherTrans
171 1 : real(RK) :: difference
172 1 : fisherTrans = getFisherTransFromcorCoef(corCoef = corCoef)
173 1 : difference = abs(fisherTrans - fisherTrans_ref) / fisherTrans_ref
174 1 : assertion = difference < tolerance
175 1 : if (Test%isDebugMode .and. .not. assertion) then
176 : ! LCOV_EXCL_START
177 : write(Test%outputUnit,"(*(g0,:,' '))")
178 : write(Test%outputUnit,"(*(g0,:,' '))") "fisherTrans_ref = ", fisherTrans_ref
179 : write(Test%outputUnit,"(*(g0,:,' '))") "fisherTrans = ", fisherTrans
180 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
181 : write(Test%outputUnit,"(*(g0,:,' '))")
182 : end if
183 : ! LCOV_EXCL_STOP
184 1 : end function test_getFisherTransFromcorCoef_1
185 :
186 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
187 :
188 1 : function test_getCumSum_IK_1() result(assertion)
189 1 : use Constants_mod, only: RK, IK
190 : implicit none
191 : logical :: assertion
192 : integer(IK) :: i
193 : integer(IK), parameter :: CumSum_ref(*) = [1_IK, 3_IK, 6_IK, 10_IK, 15_IK, 21_IK, 28_IK, 36_IK, 45_IK, 55_IK]
194 : integer(IK), parameter :: Vector(*) = [(i,i=1,size(CumSum_ref))]
195 : integer(IK), allocatable :: CumSum(:)
196 : integer(IK), allocatable :: Difference(:)
197 1 : CumSum = getCumSum(vecLen = int(size(Vector),kind=IK), Vec = Vector)
198 : ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
199 1 : if (allocated(Difference)) deallocate(Difference); allocate(Difference, mold = CumSum)
200 12 : Difference = abs(CumSum - CumSum_ref)
201 11 : assertion = all(Difference == 0_IK)
202 1 : if (Test%isDebugMode .and. .not. assertion) then
203 : ! LCOV_EXCL_START
204 : write(Test%outputUnit,"(*(g0,:,' '))")
205 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSum_ref = ", CumSum_ref
206 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSum = ", CumSum
207 : write(Test%outputUnit,"(*(g0,:,' '))") "Difference = ", Difference
208 : write(Test%outputUnit,"(*(g0,:,' '))")
209 : end if
210 : ! LCOV_EXCL_STOP
211 1 : end function test_getCumSum_IK_1
212 :
213 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
214 :
215 1 : function test_getCumSum_RK_1() result(assertion)
216 1 : use Constants_mod, only: RK, IK
217 : implicit none
218 : logical :: assertion
219 : integer(IK) :: i
220 : real(RK), parameter :: CumSum_ref(*) = [1._RK, 3._RK, 6._RK, 10._RK, 15._RK, 21._RK, 28._RK, 36._RK, 45._RK, 55._RK]
221 : real(RK), parameter :: Vector(*) = [(real(i,kind=RK),i=1,size(CumSum_ref))]
222 : real(RK), parameter :: tolerance = 1.e-14_RK
223 : real(RK), allocatable :: CumSum(:)
224 : real(RK), allocatable :: Difference(:)
225 1 : CumSum = getCumSum(vecLen = int(size(Vector),kind=IK), Vec = Vector)
226 : ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
227 1 : if (allocated(Difference)) deallocate(Difference); allocate(Difference, mold = CumSum)
228 12 : Difference = abs(CumSum - CumSum_ref)
229 11 : assertion = all(Difference < tolerance)
230 1 : if (Test%isDebugMode .and. .not. assertion) then
231 : ! LCOV_EXCL_START
232 : write(Test%outputUnit,"(*(g0,:,' '))")
233 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSum_ref = ", CumSum_ref
234 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSum = ", CumSum
235 : write(Test%outputUnit,"(*(g0,:,' '))") "Difference = ", Difference
236 : write(Test%outputUnit,"(*(g0,:,' '))")
237 : end if
238 : ! LCOV_EXCL_STOP
239 1 : end function test_getCumSum_RK_1
240 :
241 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
242 :
243 1 : function test_getCumSumReverse_IK_1() result(assertion)
244 1 : use Constants_mod, only: RK, IK
245 : implicit none
246 : logical :: assertion
247 : integer(IK) :: i
248 : integer(IK), parameter :: CumSumReverse_ref(*) = [10_IK, 19_IK, 27_IK, 34_IK, 40_IK, 45_IK, 49_IK, 52_IK, 54_IK, 55_IK]
249 : integer(IK), parameter :: Vector(*) = [(i,i=1,size(CumSumReverse_ref),1)]
250 : integer(IK), allocatable :: CumSumReverse(:)
251 : integer(IK), allocatable :: Difference(:)
252 1 : CumSumReverse = getCumSumReverse_IK(vecLen = int(size(Vector),kind=IK), Vec = Vector)
253 : ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
254 1 : if (allocated(Difference)) deallocate(Difference); allocate(Difference, mold = CumSumReverse)
255 12 : Difference = abs(CumSumReverse - CumSumReverse_ref)
256 11 : assertion = all(Difference == 0_IK)
257 1 : if (Test%isDebugMode .and. .not. assertion) then
258 : ! LCOV_EXCL_START
259 : write(Test%outputUnit,"(*(g0,:,' '))")
260 : write(Test%outputUnit,"(*(g0,:,' '))") "Vector = ", Vector
261 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSumReverse_ref = ", CumSumReverse_ref
262 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSumReverse = ", CumSumReverse
263 : write(Test%outputUnit,"(*(g0,:,' '))") "Difference = ", Difference
264 : write(Test%outputUnit,"(*(g0,:,' '))")
265 : end if
266 : ! LCOV_EXCL_STOP
267 1 : end function test_getCumSumReverse_IK_1
268 :
269 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
270 :
271 1 : function test_getCumSumReverse_RK_1() result(assertion)
272 1 : use Constants_mod, only: RK, IK
273 : implicit none
274 : logical :: assertion
275 : integer(IK) :: i
276 : real(RK), parameter :: CumSumReverse_ref(*) = [10._RK, 19._RK, 27._RK, 34._RK, 40._RK, 45._RK, 49._RK, 52._RK, 54._RK, 55._RK]
277 : real(RK), parameter :: Vector(*) = [(real(i,kind=RK),i=1,size(CumSumReverse_ref))]
278 : real(RK), parameter :: tolerance = 1.e-14_RK
279 : real(RK), allocatable :: CumSumReverse(:)
280 : real(RK), allocatable :: Difference(:)
281 1 : CumSumReverse = getCumSumReverse_RK(vecLen = int(size(Vector),kind=IK), Vec = Vector)
282 : ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
283 1 : if (allocated(Difference)) deallocate(Difference); allocate(Difference, mold = CumSumReverse)
284 12 : Difference = abs(CumSumReverse - CumSumReverse_ref)
285 11 : assertion = all(Difference < tolerance)
286 1 : if (Test%isDebugMode .and. .not. assertion) then
287 : ! LCOV_EXCL_START
288 : write(Test%outputUnit,"(*(g0,:,' '))")
289 : write(Test%outputUnit,"(*(g0,:,' '))") "Vector = ", Vector
290 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSumReverse_ref = ", CumSumReverse_ref
291 : write(Test%outputUnit,"(*(g0,:,' '))") "CumSumReverse = ", CumSumReverse
292 : write(Test%outputUnit,"(*(g0,:,' '))") "Difference = ", Difference
293 : write(Test%outputUnit,"(*(g0,:,' '))")
294 : end if
295 : ! LCOV_EXCL_STOP
296 1 : end function test_getCumSumReverse_RK_1
297 :
298 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
299 :
300 1 : function test_getLogSubExp_RK_1() result(assertion)
301 1 : use Constants_mod, only: RK, IK
302 : implicit none
303 : logical :: assertion
304 : real(RK), parameter :: logTiny2 = log(2*tiny(1._RK))
305 : real(RK), parameter :: logTiny1 = log(2._RK) + logTiny2
306 : real(RK), parameter :: tolerance = 1.e-10_RK
307 : real(RK), parameter :: logSubExp_ref = -707.7032713517043_RK
308 1 : real(RK) :: logSubExp
309 1 : real(RK) :: difference
310 1 : logSubExp = getLogSubExp_RK(logValueLarger = logTiny1, logValueSamller = logTiny2)
311 1 : difference = abs(logSubExp - logSubExp_ref)
312 1 : assertion = difference < tolerance
313 1 : if (Test%isDebugMode .and. .not. assertion) then
314 : ! LCOV_EXCL_START
315 : write(Test%outputUnit,"(*(g0,:,' '))")
316 : write(Test%outputUnit,"(*(g0,:,' '))") "logSubExp_ref = ", logSubExp_ref
317 : write(Test%outputUnit,"(*(g0,:,' '))") "logSubExp = ", logSubExp
318 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
319 : write(Test%outputUnit,"(*(g0,:,' '))")
320 : end if
321 : ! LCOV_EXCL_STOP
322 1 : end function test_getLogSubExp_RK_1
323 :
324 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
325 :
326 1 : function test_getLogSumExp_RK_1() result(assertion)
327 1 : use Constants_mod, only: RK, IK
328 : implicit none
329 : logical :: assertion
330 : real(RK), parameter :: LogValue(*) = [ log(0.5*huge(1._RK)), log(0.9*huge(1._RK)), log(0.1*huge(1._RK)) ]
331 : real(RK), parameter :: tolerance = 1.e-10_RK
332 : real(RK), parameter :: logSumExp_ref = 710.1881779865910_RK
333 1 : real(RK) :: logSumExp
334 1 : real(RK) :: difference
335 1 : logSumExp = getLogSumExp_RK(lenLogValue = size(LogValue), LogValue = LogValue, maxLogValue = maxval(LogValue))
336 1 : difference = abs(logSumExp - logSumExp_ref)
337 1 : assertion = difference < tolerance
338 1 : if (Test%isDebugMode .and. .not. assertion) then
339 : ! LCOV_EXCL_START
340 : write(Test%outputUnit,"(*(g0,:,' '))")
341 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp_ref = ", logSumExp_ref
342 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp = ", logSumExp
343 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
344 : write(Test%outputUnit,"(*(g0,:,' '))")
345 : end if
346 : ! LCOV_EXCL_STOP
347 1 : end function test_getLogSumExp_RK_1
348 :
349 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
350 :
351 1 : function test_getLogSumExp_RK_2() result(assertion)
352 1 : use Constants_mod, only: RK, IK
353 : implicit none
354 : logical :: assertion
355 : real(RK), parameter :: LogValue(*) = [ log(0.5*huge(1._RK)), log(0.9*huge(1._RK)), log(0.1*huge(1._RK)) ]
356 : real(RK), parameter :: tolerance = 1.e-10_RK
357 : real(RK), parameter :: logSumExp_ref = 710.1881779865910_RK
358 1 : real(RK) :: logSumExp
359 1 : real(RK) :: difference
360 1 : logSumExp = getLogSumExp_RK(lenLogValue = size(LogValue), LogValue = LogValue)
361 1 : difference = abs(logSumExp - logSumExp_ref)
362 1 : assertion = difference < tolerance
363 1 : if (Test%isDebugMode .and. .not. assertion) then
364 : ! LCOV_EXCL_START
365 : write(Test%outputUnit,"(*(g0,:,' '))")
366 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp_ref = ", logSumExp_ref
367 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp = ", logSumExp
368 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
369 : write(Test%outputUnit,"(*(g0,:,' '))")
370 : end if
371 : ! LCOV_EXCL_STOP
372 1 : end function test_getLogSumExp_RK_2
373 :
374 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
375 :
376 1 : function test_getLogSumExp_CK_1() result(assertion)
377 1 : use Constants_mod, only: RK, IK
378 : implicit none
379 : logical :: assertion
380 : complex(RK), parameter :: LogValue(*) = [ cmplx( log(0.5*huge(1._RK)), 0._RK, kind = RK ) &
381 : , cmplx( log(0.9*huge(1._RK)), 0._RK, kind = RK ) &
382 : , cmplx( log(0.1*huge(1._RK)), 0._RK, kind = RK ) &
383 : ]
384 : real(RK), parameter :: logSumExp_ref = cmplx(710.1881779865910_RK, 0._RK, RK)
385 : complex(RK), parameter :: tolerance = cmplx(1.e-10_RK, 0._RK, RK)
386 : complex(RK) :: logSumExp
387 : complex(RK) :: difference
388 : logSumExp = getLogSumExp_CK ( lenLogValue = int(size(LogValue),kind=IK) &
389 : , LogValue = LogValue &
390 : , maxLogValue = cmplx( maxval(real(LogValue,kind=RK)), kind = RK ) &
391 5 : )
392 1 : difference = abs(logSumExp - logSumExp_ref)
393 1 : assertion = real(difference,RK) < real(tolerance,RK)
394 1 : if (Test%isDebugMode .and. .not. assertion) then
395 : ! LCOV_EXCL_START
396 : write(Test%outputUnit,"(*(g0,:,' '))")
397 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp_ref = ", logSumExp_ref
398 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp = ", logSumExp
399 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
400 : write(Test%outputUnit,"(*(g0,:,' '))")
401 : end if
402 : ! LCOV_EXCL_STOP
403 1 : end function test_getLogSumExp_CK_1
404 :
405 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
406 :
407 1 : function test_getLogSumExp_CK_2() result(assertion)
408 1 : use Constants_mod, only: RK, IK
409 : implicit none
410 : logical :: assertion
411 : complex(RK), parameter :: LogValue(*) = [ cmplx( log(0.5*huge(1._RK)), 0._RK, kind = RK ) &
412 : , cmplx( log(0.9*huge(1._RK)), 0._RK, kind = RK ) &
413 : , cmplx( log(0.1*huge(1._RK)), 0._RK, kind = RK ) &
414 : ]
415 : real(RK), parameter :: logSumExp_ref = cmplx(710.1881779865910_RK, 0._RK, RK)
416 : complex(RK), parameter :: tolerance = cmplx(1.e-10_RK, 0._RK, RK)
417 : complex(RK) :: logSumExp
418 : complex(RK) :: difference
419 1 : logSumExp = getLogSumExp_CK(lenLogValue = int(size(LogValue),IK), LogValue = LogValue)
420 1 : difference = abs(logSumExp - logSumExp_ref)
421 1 : assertion = real(difference,RK) < real(tolerance,RK)
422 1 : if (Test%isDebugMode .and. .not. assertion) then
423 : ! LCOV_EXCL_START
424 : write(Test%outputUnit,"(*(g0,:,' '))")
425 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp_ref = ", logSumExp_ref
426 : write(Test%outputUnit,"(*(g0,:,' '))") "logSumExp = ", logSumExp
427 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
428 : write(Test%outputUnit,"(*(g0,:,' '))")
429 : end if
430 : ! LCOV_EXCL_STOP
431 1 : end function test_getLogSumExp_CK_2
432 :
433 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
434 :
435 1 : function test_getLogEggBoxSD_RK_1() result(assertion)
436 1 : use Constants_mod, only: RK, IK
437 : implicit none
438 : logical :: assertion
439 : real(RK), parameter :: coef = 1._RK
440 : real(RK), parameter :: point = 1.e1_RK
441 : real(RK), parameter :: constant = 1._RK
442 : real(RK), parameter :: exponent = 1._RK
443 : real(RK), parameter :: tolerance = 1.e-10_RK
444 : real(RK), parameter :: logEggBox_ref = 0.16092847092354756_RK
445 1 : real(RK) :: logEggBox
446 1 : real(RK) :: difference
447 1 : logEggBox = getLogEggBoxSD_RK(constant = constant, exponent = exponent, coef = coef, point = point)
448 1 : difference = abs(logEggBox - logEggBox_ref)
449 1 : assertion = difference < tolerance
450 1 : if (Test%isDebugMode .and. .not. assertion) then
451 : ! LCOV_EXCL_START
452 : write(Test%outputUnit,"(*(g0,:,' '))")
453 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox_ref = ", logEggBox_ref
454 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox = ", logEggBox
455 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
456 : write(Test%outputUnit,"(*(g0,:,' '))")
457 : end if
458 : ! LCOV_EXCL_STOP
459 1 : end function test_getLogEggBoxSD_RK_1
460 :
461 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
462 :
463 1 : function test_getLogEggBoxSD_CK_1() result(assertion)
464 1 : use Constants_mod, only: RK, IK
465 : implicit none
466 : logical :: assertion
467 : complex(RK), parameter :: coef = cmplx(1._RK, 0._RK, RK)
468 : complex(RK), parameter :: point = cmplx(1.e1_RK, 0._RK, RK)
469 : complex(RK), parameter :: constant = cmplx(1._RK, 0._RK, RK)
470 : complex(RK), parameter :: exponent = cmplx(1._RK, 0._RK, RK)
471 : complex(RK), parameter :: tolerance = cmplx(1.e-10_RK, 0._RK, RK)
472 : complex(RK), parameter :: logEggBox_ref = cmplx(0.16092847092354756_RK, 0._RK, RK)
473 : complex(RK) :: logEggBox
474 : complex(RK) :: difference
475 1 : logEggBox = getLogEggBoxSD_CK(constant = constant, exponent = exponent, coef = coef, point = point)
476 1 : difference = abs(logEggBox - logEggBox_ref)
477 1 : assertion = real(difference,RK) < real(tolerance,RK)
478 1 : if (Test%isDebugMode .and. .not. assertion) then
479 : ! LCOV_EXCL_START
480 : write(Test%outputUnit,"(*(g0,:,' '))")
481 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox_ref = ", logEggBox_ref
482 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox = ", logEggBox
483 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
484 : write(Test%outputUnit,"(*(g0,:,' '))")
485 : end if
486 : ! LCOV_EXCL_STOP
487 1 : end function test_getLogEggBoxSD_CK_1
488 :
489 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
490 :
491 1 : function test_getLogEggBoxMD_RK_1() result(assertion)
492 1 : use Constants_mod, only: RK, IK
493 : implicit none
494 : logical :: assertion
495 : integer(IK) , parameter :: nd = 3_IK
496 : real(RK) , parameter :: coef(nd) = [1._RK, 2._RK, 3._RK]
497 : real(RK) , parameter :: point(nd) = [1.e1_RK, 1.e2_RK, 1.e3_RK]
498 : real(RK) , parameter :: constant = 1._RK
499 : real(RK) , parameter :: exponent = 1._RK
500 : real(RK) , parameter :: tolerance = 1.e-10_RK
501 : real(RK) , parameter :: logEggBox_ref = 1.3988445480199623_RK
502 1 : real(RK) :: logEggBox
503 1 : real(RK) :: difference
504 1 : logEggBox = getLogEggBoxMD_RK(nd = nd, constant = constant, exponent = exponent, coef = coef, point = point)
505 1 : difference = abs(logEggBox - logEggBox_ref)
506 1 : assertion = difference < tolerance
507 1 : if (Test%isDebugMode .and. .not. assertion) then
508 : ! LCOV_EXCL_START
509 : write(Test%outputUnit,"(*(g0,:,' '))")
510 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox_ref = ", logEggBox_ref
511 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox = ", logEggBox
512 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
513 : write(Test%outputUnit,"(*(g0,:,' '))")
514 : end if
515 : ! LCOV_EXCL_STOP
516 1 : end function test_getLogEggBoxMD_RK_1
517 :
518 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
519 :
520 1 : function test_getLogEggBoxMD_CK_1() result(assertion)
521 1 : use Constants_mod, only: RK, IK
522 : implicit none
523 : logical :: assertion
524 : integer(IK), parameter :: nd = 3_IK
525 : complex(RK), parameter :: coef(nd) = [cmplx(1._RK, 0._RK, RK), cmplx(2._RK, 0._RK, RK), cmplx(3._RK, 0._RK, RK)]
526 : complex(RK), parameter :: point(nd) = [cmplx(1.e1_RK, 0._RK, RK), cmplx(1.e2_RK, 0._RK, RK), cmplx(1.e3_RK, 0._RK, RK)]
527 : complex(RK), parameter :: constant = cmplx(1._RK, 0._RK, RK)
528 : complex(RK), parameter :: exponent = cmplx(1._RK, 0._RK, RK)
529 : complex(RK), parameter :: tolerance = cmplx(1.e-10_RK, 0._RK, RK)
530 : complex(RK), parameter :: logEggBox_ref = cmplx(1.3988445480199623_RK, 0._RK, RK)
531 : complex(RK) :: logEggBox
532 : complex(RK) :: difference
533 1 : logEggBox = getLogEggBoxMD_CK(nd = nd, constant = constant, exponent = exponent, coef = coef, point = point)
534 1 : difference = abs(logEggBox - logEggBox_ref)
535 1 : assertion = real(difference,RK) < real(tolerance,RK)
536 1 : if (Test%isDebugMode .and. .not. assertion) then
537 : ! LCOV_EXCL_START
538 : write(Test%outputUnit,"(*(g0,:,' '))")
539 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox_ref = ", logEggBox_ref
540 : write(Test%outputUnit,"(*(g0,:,' '))") "logEggBox = ", logEggBox
541 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
542 : write(Test%outputUnit,"(*(g0,:,' '))")
543 : end if
544 : ! LCOV_EXCL_STOP
545 1 : end function test_getLogEggBoxMD_CK_1
546 :
547 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
548 :
549 1 : function test_getFactorial_1() result(assertion)
550 : use Constants_mod, only: RK, IK ! LCOV_EXCL_LINE
551 : implicit none
552 : logical :: assertion
553 : integer(IK), parameter :: factorial_ref = 3628800_IK
554 : integer(IK), parameter :: positiveInteger = 10_IK
555 1 : real(RK) :: difference
556 1 : real(RK) :: factorial
557 1 : factorial = getFactorial(positiveInteger = positiveInteger)
558 1 : difference = abs(factorial - factorial_ref)
559 1 : assertion = difference == 0_IK
560 1 : if (Test%isDebugMode .and. .not. assertion) then
561 : ! LCOV_EXCL_START
562 : write(Test%outputUnit,"(*(g0,:,' '))")
563 : write(Test%outputUnit,"(*(g0,:,' '))") "factorial_ref = ", factorial_ref
564 : write(Test%outputUnit,"(*(g0,:,' '))") "factorial = ", factorial
565 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
566 : write(Test%outputUnit,"(*(g0,:,' '))")
567 : end if
568 : ! LCOV_EXCL_STOP
569 1 : end function test_getFactorial_1
570 :
571 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
572 :
573 1 : function test_getLogFactorial_1() result(assertion)
574 1 : use Constants_mod, only: RK, IK
575 : implicit none
576 : logical :: assertion
577 : integer(IK) , parameter :: positiveInteger = 10_IK
578 : real(RK) , parameter :: logFactorial_ref = 15.10441257307552_RK
579 : real(RK) , parameter :: tolerance = 1.e-10_RK
580 1 : real(RK) :: logFactorial
581 1 : real(RK) :: difference
582 1 : logFactorial = getLogFactorial(positiveInteger = positiveInteger)
583 1 : difference = abs(logFactorial - logFactorial_ref)
584 1 : assertion = difference < tolerance
585 1 : if (Test%isDebugMode .and. .not. assertion) then
586 : ! LCOV_EXCL_START
587 : write(Test%outputUnit,"(*(g0,:,' '))")
588 : write(Test%outputUnit,"(*(g0,:,' '))") "logFactorial_ref = ", logFactorial_ref
589 : write(Test%outputUnit,"(*(g0,:,' '))") "logFactorial = ", logFactorial
590 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
591 : write(Test%outputUnit,"(*(g0,:,' '))")
592 : end if
593 : ! LCOV_EXCL_STOP
594 1 : end function test_getLogFactorial_1
595 :
596 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
597 :
598 1 : function test_getEllVolCoef_1() result(assertion)
599 1 : use Constants_mod, only: RK, IK, PI
600 : implicit none
601 : logical :: assertion
602 : integer(IK) , parameter :: nd = 10_IK
603 : real(RK) , parameter :: tolerance = 1.e-10_RK
604 : real(RK) , parameter :: ellVolCoef_ref = PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ! 2.550164039877345_RK
605 1 : real(RK) :: ellVolCoef
606 1 : real(RK) :: difference
607 1 : ellVolCoef = getEllVolCoef(nd = nd)
608 1 : difference = abs(ellVolCoef - ellVolCoef_ref)
609 1 : assertion = difference < tolerance
610 1 : if (Test%isDebugMode .and. .not. assertion) then
611 : ! LCOV_EXCL_START
612 : write(Test%outputUnit,"(*(g0,:,' '))")
613 : write(Test%outputUnit,"(*(g0,:,' '))") "ellVolCoef_ref = ", ellVolCoef_ref
614 : write(Test%outputUnit,"(*(g0,:,' '))") "ellVolCoef = ", ellVolCoef
615 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
616 : write(Test%outputUnit,"(*(g0,:,' '))")
617 : end if
618 : ! LCOV_EXCL_STOP
619 1 : end function test_getEllVolCoef_1
620 :
621 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
622 :
623 1 : function test_getEllVolCoef_2() result(assertion)
624 1 : use Constants_mod, only: RK, IK, PI
625 : implicit none
626 : logical :: assertion
627 : integer(IK) , parameter :: nd = 11_IK
628 : real(RK) , parameter :: tolerance = 1.e-10_RK
629 : real(RK) , parameter :: ellVolCoef_ref = PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ! 1.884103879389900_RK
630 1 : real(RK) :: ellVolCoef
631 1 : real(RK) :: difference
632 : !integer(IK) :: i
633 : !do i = 1, 10000000
634 1 : ellVolCoef = getEllVolCoef(nd = nd)
635 : !end do
636 1 : difference = abs(ellVolCoef - ellVolCoef_ref)
637 1 : assertion = difference < tolerance
638 1 : if (Test%isDebugMode .and. .not. assertion) then
639 : ! LCOV_EXCL_START
640 : write(Test%outputUnit,"(*(g0,:,' '))")
641 : write(Test%outputUnit,"(*(g0,:,' '))") "ellVolCoef_ref = ", ellVolCoef_ref
642 : write(Test%outputUnit,"(*(g0,:,' '))") "ellVolCoef = ", ellVolCoef
643 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
644 : write(Test%outputUnit,"(*(g0,:,' '))")
645 : end if
646 : ! LCOV_EXCL_STOP
647 1 : end function test_getEllVolCoef_2
648 :
649 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
650 :
651 1 : function test_getLogEllVolCoef_1() result(assertion)
652 1 : use Constants_mod, only: RK, IK, PI
653 : implicit none
654 : logical :: assertion
655 : integer(IK) , parameter :: nd = 10_IK
656 : real(RK) , parameter :: tolerance = 1.e-10_RK
657 : real(RK) , parameter :: logEllVolCoef_ref = log( PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ) ! 0.936157686464955_RK
658 1 : real(RK) :: logEllVolCoef
659 1 : real(RK) :: difference
660 1 : logEllVolCoef = getLogEllVolCoef(nd = nd)
661 1 : difference = abs(logEllVolCoef - logEllVolCoef_ref)
662 1 : assertion = difference < tolerance
663 1 : if (Test%isDebugMode .and. .not. assertion) then
664 : ! LCOV_EXCL_START
665 : write(Test%outputUnit,"(*(g0,:,' '))")
666 : write(Test%outputUnit,"(*(g0,:,' '))") "logEllVolCoef_ref = ", logEllVolCoef_ref
667 : write(Test%outputUnit,"(*(g0,:,' '))") "logEllVolCoef = ", logEllVolCoef
668 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
669 : write(Test%outputUnit,"(*(g0,:,' '))")
670 : end if
671 : ! LCOV_EXCL_STOP
672 1 : end function test_getLogEllVolCoef_1
673 :
674 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
675 :
676 1 : function test_getLogEllVolCoef_2() result(assertion)
677 1 : use Constants_mod, only: RK, IK, PI
678 : implicit none
679 : logical :: assertion
680 : integer(IK) , parameter :: nd = 11_IK
681 : real(RK) , parameter :: tolerance = 1.e-10_RK
682 : real(RK) , parameter :: logEllVolCoef_ref = log( PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ) ! 0.633452312314559_RK
683 1 : real(RK) :: logEllVolCoef
684 1 : real(RK) :: difference
685 : !integer(IK) :: i
686 : !do i = 1, 10000000
687 1 : logEllVolCoef = getLogEllVolCoef(nd = nd)
688 : !end do
689 1 : difference = abs(logEllVolCoef - logEllVolCoef_ref)
690 1 : assertion = difference < tolerance
691 1 : if (Test%isDebugMode .and. .not. assertion) then
692 : ! LCOV_EXCL_START
693 : write(Test%outputUnit,"(*(g0,:,' '))")
694 : write(Test%outputUnit,"(*(g0,:,' '))") "logEllVolCoef_ref = ", logEllVolCoef_ref
695 : write(Test%outputUnit,"(*(g0,:,' '))") "logEllVolCoef = ", logEllVolCoef
696 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
697 : write(Test%outputUnit,"(*(g0,:,' '))")
698 : end if
699 : ! LCOV_EXCL_STOP
700 1 : end function test_getLogEllVolCoef_2
701 :
702 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
703 :
704 1 : function test_getLogVolUnitBall_1() result(assertion)
705 1 : use Constants_mod, only: RK, IK, PI
706 : implicit none
707 : logical :: assertion
708 : integer(IK) , parameter :: nd = 10_IK
709 : real(RK) , parameter :: tolerance = 1.e-10_RK
710 : real(RK) , parameter :: logVolUnitBall_ref = log( PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ) ! .9361576864649548_RK
711 1 : real(RK) :: logVolUnitBall
712 1 : real(RK) :: difference
713 1 : logVolUnitBall = getLogVolUnitBall(nd = nd)
714 1 : difference = abs(logVolUnitBall - logVolUnitBall_ref)
715 1 : assertion = difference < tolerance
716 1 : if (Test%isDebugMode .and. .not. assertion) then
717 : ! LCOV_EXCL_START
718 : write(Test%outputUnit,"(*(g0,:,' '))")
719 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolUnitBall_ref = ", logVolUnitBall_ref
720 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolUnitBall = ", logVolUnitBall
721 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
722 : write(Test%outputUnit,"(*(g0,:,' '))")
723 : end if
724 : ! LCOV_EXCL_STOP
725 1 : end function test_getLogVolUnitBall_1
726 :
727 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
728 :
729 1 : function test_getLogVolUnitBall_2() result(assertion)
730 1 : use Constants_mod, only: RK, IK, PI
731 : implicit none
732 : logical :: assertion
733 : integer(IK) , parameter :: nd = 11_IK
734 : real(RK) , parameter :: tolerance = 1.e-10_RK
735 : real(RK) , parameter :: logVolUnitBall_ref = log( PI**(0.5_RK*nd) / gamma(0.5_RK*nd+1._RK) ) ! .6334523123145592_RK
736 1 : real(RK) :: logVolUnitBall
737 1 : real(RK) :: difference
738 : !integer(IK) :: i
739 : !do i = 1, 10000000
740 1 : logVolUnitBall = getLogVolUnitBall(nd = nd)
741 : !end do
742 1 : difference = abs(logVolUnitBall - logVolUnitBall_ref)
743 1 : assertion = difference < tolerance
744 1 : if (Test%isDebugMode .and. .not. assertion) then
745 : ! LCOV_EXCL_START
746 : write(Test%outputUnit,"(*(g0,:,' '))")
747 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolUnitBall_ref = ", logVolUnitBall_ref
748 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolUnitBall = ", logVolUnitBall
749 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
750 : write(Test%outputUnit,"(*(g0,:,' '))")
751 : end if
752 : ! LCOV_EXCL_STOP
753 1 : end function test_getLogVolUnitBall_2
754 :
755 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
756 :
757 1 : function test_getLogVolEllipsoid_1() result(assertion)
758 1 : use Constants_mod, only: RK, IK
759 : implicit none
760 : logical :: assertion
761 : integer(IK) , parameter :: nd = 10_IK
762 : real(RK) , parameter :: tolerance = 1.e-10_RK
763 : real(RK) , parameter :: logSqrtDetCovMat = 2._RK
764 : real(RK) , parameter :: logVolEllipsoid_ref = 2.936157686464955_RK
765 1 : real(RK) :: logVolEllipsoid
766 1 : real(RK) :: difference
767 1 : logVolEllipsoid = getLogVolEllipsoid(nd = nd, logSqrtDetCovMat = logSqrtDetCovMat)
768 1 : difference = abs(logVolEllipsoid - logVolEllipsoid_ref)
769 1 : assertion = difference < tolerance
770 1 : if (Test%isDebugMode .and. .not. assertion) then
771 : ! LCOV_EXCL_START
772 : write(Test%outputUnit,"(*(g0,:,' '))")
773 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolEllipsoid_ref = ", logVolEllipsoid_ref
774 : write(Test%outputUnit,"(*(g0,:,' '))") "logVolEllipsoid = ", logVolEllipsoid
775 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", difference
776 : write(Test%outputUnit,"(*(g0,:,' '))")
777 : end if
778 : ! LCOV_EXCL_STOP
779 1 : end function test_getLogVolEllipsoid_1
780 :
781 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
782 :
783 1 : function test_getLogVolEllipsoids_1() result(assertion)
784 1 : use Constants_mod, only: RK, IK
785 : implicit none
786 : integer(IK) :: i
787 : logical :: assertion
788 : integer(IK) , parameter :: nEllipsoid = 2_IK
789 : real(RK) , parameter :: tolerance = 1.e-10_RK
790 : integer(IK) , parameter :: nd = 2_IK
791 : real(RK) , parameter :: LogSqrtDetCovMat(nEllipsoid) = [ (log(real(i,RK)), i = 1, nEllipsoid) ]
792 : real(RK) , parameter :: EllipsoidVolume_ref(*) = [ 1.144729885849400_RK, 1.837877066409345_RK ]
793 : real(RK), allocatable :: EllipsoidVolume(:)
794 : real(RK), allocatable :: Difference(:)
795 1 : EllipsoidVolume = getLogVolEllipsoids(nd = nd, nEllipsoid = nEllipsoid, LogSqrtDetCovMat = LogSqrtDetCovMat)
796 : ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
797 1 : if (allocated(Difference)) deallocate(Difference); allocate(Difference, mold = EllipsoidVolume)
798 4 : Difference = abs(EllipsoidVolume - EllipsoidVolume_ref)
799 3 : assertion = all(Difference < tolerance)
800 1 : if (Test%isDebugMode .and. .not. assertion) then
801 : ! LCOV_EXCL_START
802 : write(Test%outputUnit,"(*(g0,:,' '))")
803 : write(Test%outputUnit,"(*(g0,:,' '))") "EllipsoidVolume_ref = ", EllipsoidVolume_ref
804 : write(Test%outputUnit,"(*(g0,:,' '))") "EllipsoidVolume = ", EllipsoidVolume
805 : write(Test%outputUnit,"(*(g0,:,' '))") "difference = ", Difference
806 : write(Test%outputUnit,"(*(g0,:,' '))")
807 : end if
808 : ! LCOV_EXCL_STOP
809 1 : end function test_getLogVolEllipsoids_1
810 :
811 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
812 :
813 : !> \brief
814 : !> Test [getLowerGamma](@ref math_mod::getlowergamma) with a small `tolerance` input optional argument.
815 1 : function test_getLowerGamma_1() result(assertion)
816 1 : use Constants_mod, only: RK, IK
817 : implicit none
818 : integer(IK) :: i
819 : logical :: assertion
820 : integer(IK) , parameter :: ntest = 3
821 : real(RK) , parameter :: tolerance = 1.e-10_RK ! 1.e-7_RK
822 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
823 : real(RK) , parameter :: UpperLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
824 : real(RK) , parameter :: LowerGamma_ref(ntest) = [ 0.632120558828558_RK, 0.184736755476228_RK, 0.999817189367018_RK ]
825 : real(RK) :: LowerGamma(ntest)
826 1 : real(RK) :: difference
827 4 : do i = 1, ntest
828 3 : LowerGamma(i) = getLowerGamma ( exponent = Exponent(i) &
829 3 : , logGammaExponent = log_gamma(Exponent(i)) &
830 3 : , upperLim = UpperLim(i) &
831 : , tolerance = tolerance &
832 3 : )
833 3 : difference = 2 * abs(LowerGamma(i) - LowerGamma_ref(i)) / LowerGamma_ref(i)
834 3 : assertion = difference < tolerance
835 4 : if (Test%isDebugMode .and. .not. assertion) then
836 : ! LCOV_EXCL_START
837 : write(Test%outputUnit,"(*(g0,:,', '))")
838 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference LowerGamma, Computed LowerGamma, difference:"
839 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), UpperLim(i), LowerGamma(i), LowerGamma_ref(i), difference
840 : write(Test%outputUnit,"(*(g0,:,', '))")
841 : end if
842 : ! LCOV_EXCL_STOP
843 : end do
844 1 : end function test_getLowerGamma_1
845 :
846 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
847 :
848 : !> \brief
849 : !> Test [getLowerGamma](@ref math_mod::getlowergamma) with a medium `tolerance` input optional argument.
850 1 : function test_getLowerGamma_2() result(assertion)
851 1 : use Constants_mod, only: RK, IK
852 : implicit none
853 : integer(IK) :: i
854 : logical :: assertion
855 : integer(IK) , parameter :: ntest = 3
856 : real(RK) , parameter :: tolerance = 1.e-7_RK
857 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
858 : real(RK) , parameter :: UpperLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
859 : real(RK) , parameter :: LowerGamma_ref(ntest) = [ 0.632120558828558_RK, 0.184736755476228_RK, 0.999817189367018_RK ]
860 : real(RK) :: LowerGamma(ntest)
861 1 : real(RK) :: difference
862 4 : do i = 1, ntest
863 3 : LowerGamma(i) = getLowerGamma ( exponent = Exponent(i) &
864 3 : , logGammaExponent = log_gamma(Exponent(i)) &
865 3 : , upperLim = UpperLim(i) &
866 : , tolerance = tolerance &
867 3 : )
868 3 : difference = 2 * abs(LowerGamma(i) - LowerGamma_ref(i)) / LowerGamma_ref(i)
869 3 : assertion = difference < tolerance
870 4 : if (Test%isDebugMode .and. .not. assertion) then
871 : ! LCOV_EXCL_START
872 : write(Test%outputUnit,"(*(g0,:,', '))")
873 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference LowerGamma, Computed LowerGamma, difference:"
874 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), UpperLim(i), LowerGamma(i), LowerGamma_ref(i), difference
875 : write(Test%outputUnit,"(*(g0,:,', '))")
876 : end if
877 : ! LCOV_EXCL_STOP
878 : end do
879 1 : end function test_getLowerGamma_2
880 :
881 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
882 :
883 : !> \brief
884 : !> Test [getLowerGamma](@ref math_mod::getlowergamma) with a large `tolerance` input optional argument.
885 1 : function test_getLowerGamma_3() result(assertion)
886 1 : use Constants_mod, only: RK, IK
887 : implicit none
888 : integer(IK) :: i
889 : logical :: assertion
890 : integer(IK) , parameter :: ntest = 3
891 : real(RK) , parameter :: tolerance = 1.e-3_RK
892 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
893 : real(RK) , parameter :: UpperLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
894 : real(RK) , parameter :: LowerGamma_ref(ntest) = [ 0.632120558828558_RK, 0.184736755476228_RK, 0.999817189367018_RK ]
895 : real(RK) :: LowerGamma(ntest)
896 1 : real(RK) :: difference
897 4 : do i = 1, ntest
898 3 : LowerGamma(i) = getLowerGamma ( exponent = Exponent(i) &
899 3 : , logGammaExponent = log_gamma(Exponent(i)) &
900 3 : , upperLim = UpperLim(i) &
901 : , tolerance = tolerance &
902 3 : )
903 3 : difference = 2 * abs(LowerGamma(i) - LowerGamma_ref(i)) / LowerGamma_ref(i)
904 3 : assertion = difference < tolerance
905 4 : if (Test%isDebugMode .and. .not. assertion) then
906 : ! LCOV_EXCL_START
907 : write(Test%outputUnit,"(*(g0,:,', '))")
908 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference LowerGamma, Computed LowerGamma, difference:"
909 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), UpperLim(i), LowerGamma(i), LowerGamma_ref(i), difference
910 : write(Test%outputUnit,"(*(g0,:,', '))")
911 : end if
912 : ! LCOV_EXCL_STOP
913 : end do
914 1 : end function test_getLowerGamma_3
915 :
916 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
917 :
918 : !> \brief
919 : !> Test [getLowerGamma](@ref math_mod::getlowergamma) without the `tolerance` input optional argument,
920 : !> in which case, the procedure should default to `epsilon` for the value of tolerance.
921 1 : function test_getLowerGamma_4() result(assertion)
922 1 : use Constants_mod, only: RK, IK
923 : implicit none
924 : integer(IK) :: i
925 : logical :: assertion
926 : integer(IK) , parameter :: ntest = 3
927 : real(RK) , parameter :: tolerance = epsilon(1._RK)
928 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
929 : real(RK) , parameter :: UpperLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
930 : real(RK) , parameter :: LowerGamma_ref(ntest) = [ 0.632120558828558_RK, 0.184736755476228_RK, 0.999817189367018_RK ]
931 : real(RK) :: LowerGamma(ntest)
932 1 : real(RK) :: difference
933 4 : do i = 1, ntest
934 3 : LowerGamma(i) = getLowerGamma ( exponent = Exponent(i) &
935 3 : , logGammaExponent = log_gamma(Exponent(i)) &
936 3 : , upperLim = UpperLim(i) &
937 3 : )
938 3 : difference = 2 * abs(LowerGamma(i) - LowerGamma_ref(i)) / LowerGamma_ref(i)
939 3 : assertion = difference < 1000 * tolerance
940 4 : if (Test%isDebugMode .and. .not. assertion) then
941 : ! LCOV_EXCL_START
942 : write(Test%outputUnit,"(*(g0,:,', '))")
943 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference LowerGamma, Computed LowerGamma, difference:"
944 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), UpperLim(i), LowerGamma(i), LowerGamma_ref(i), difference
945 : write(Test%outputUnit,"(*(g0,:,', '))")
946 : end if
947 : ! LCOV_EXCL_STOP
948 : end do
949 1 : end function test_getLowerGamma_4
950 :
951 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
952 :
953 : !> \brief
954 : !> Test [getLowerGamma](@ref math_mod::getlowergamma) with a wrong positive value for the input argument `upperLim`.
955 1 : function test_getLowerGamma_5() result(assertion)
956 1 : use Constants_mod, only: RK, IK, HUGE_RK
957 : implicit none
958 : logical :: assertion
959 : real(RK) , parameter :: exponent = +1.0_RK
960 : real(RK) , parameter :: upperlim = -1.0_RK
961 : real(RK) , parameter :: lowerGamma_ref = -HUGE_RK
962 1 : real(RK) :: lowerGamma
963 : lowerGamma = getLowerGamma ( exponent = exponent &
964 : , logGammaExponent = log_gamma(exponent) &
965 : , upperLim = upperLim &
966 1 : )
967 1 : assertion = lowerGamma == lowerGamma_ref
968 1 : if (Test%isDebugMode .and. .not. assertion) then
969 : ! LCOV_EXCL_START
970 : write(Test%outputUnit,"(*(g0,:,', '))")
971 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference LowerGamma, Computed LowerGamma:"
972 : write(Test%outputUnit,"(*(g0,:,', '))") exponent, upperlim, lowerGamma, lowerGamma_ref
973 : write(Test%outputUnit,"(*(g0,:,', '))")
974 : end if
975 : ! LCOV_EXCL_STOP
976 1 : end function test_getLowerGamma_5
977 :
978 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
979 :
980 : !> \brief
981 : !> Test [getUpperGamma](@ref math_mod::getuppergamma) with a small `tolerance` input optional argument.
982 1 : function test_getUpperGamma_1() result(assertion)
983 1 : use Constants_mod, only: RK, IK
984 : implicit none
985 : integer(IK) :: i
986 : logical :: assertion
987 : integer(IK) , parameter :: ntest = 3
988 : real(RK) , parameter :: tolerance = 1.e-10_RK ! 1.e-7_RK
989 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
990 : real(RK) , parameter :: LowerLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
991 : real(RK) , parameter :: UpperGamma_ref(ntest) = [0.367879441171442_RK, 0.815263244523772_RK, 1.828106329818355e-04_RK]
992 : real(RK) :: UpperGamma(ntest)
993 1 : real(RK) :: difference
994 4 : do i = 1, ntest
995 3 : UpperGamma(i) = getUpperGamma ( exponent = Exponent(i) &
996 3 : , logGammaExponent = log_gamma(Exponent(i)) &
997 3 : , lowerLim = LowerLim(i) &
998 : , tolerance = tolerance &
999 3 : )
1000 3 : difference = 2 * abs(UpperGamma(i) - UpperGamma_ref(i)) / UpperGamma_ref(i)
1001 3 : assertion = difference < tolerance
1002 4 : if (Test%isDebugMode .and. .not. assertion) then
1003 : ! LCOV_EXCL_START
1004 : write(Test%outputUnit,"(*(g0,:,', '))")
1005 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference UpperGamma, Computed UpperGamma, difference:"
1006 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), LowerLim(i), UpperGamma(i), UpperGamma_ref(i), difference
1007 : write(Test%outputUnit,"(*(g0,:,', '))")
1008 : end if
1009 : ! LCOV_EXCL_STOP
1010 : end do
1011 1 : end function test_getUpperGamma_1
1012 :
1013 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1014 :
1015 : !> \brief
1016 : !> Test [getUpperGamma](@ref math_mod::getuppergamma) with a medium `tolerance` input optional argument.
1017 1 : function test_getUpperGamma_2() result(assertion)
1018 1 : use Constants_mod, only: RK, IK
1019 : implicit none
1020 : integer(IK) :: i
1021 : logical :: assertion
1022 : integer(IK) , parameter :: ntest = 3
1023 : real(RK) , parameter :: tolerance = 1.e-7_RK
1024 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
1025 : real(RK) , parameter :: LowerLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
1026 : real(RK) , parameter :: UpperGamma_ref(ntest) = [0.367879441171442_RK, 0.815263244523772_RK, 1.828106329818355e-04_RK]
1027 : real(RK) :: UpperGamma(ntest)
1028 1 : real(RK) :: difference
1029 4 : do i = 1, ntest
1030 3 : UpperGamma(i) = getUpperGamma ( exponent = Exponent(i) &
1031 3 : , logGammaExponent = log_gamma(Exponent(i)) &
1032 3 : , lowerLim = LowerLim(i) &
1033 : , tolerance = tolerance &
1034 3 : )
1035 3 : difference = 2 * abs(UpperGamma(i) - UpperGamma_ref(i)) / UpperGamma_ref(i)
1036 3 : assertion = difference < tolerance
1037 4 : if (Test%isDebugMode .and. .not. assertion) then
1038 : ! LCOV_EXCL_START
1039 : write(Test%outputUnit,"(*(g0,:,', '))")
1040 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference UpperGamma, Computed UpperGamma, difference:"
1041 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), LowerLim(i), UpperGamma(i), UpperGamma_ref(i), difference
1042 : write(Test%outputUnit,"(*(g0,:,', '))")
1043 : end if
1044 : ! LCOV_EXCL_STOP
1045 : end do
1046 1 : end function test_getUpperGamma_2
1047 :
1048 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1049 :
1050 : !> \brief
1051 : !> Test [getUpperGamma](@ref math_mod::getuppergamma) with a coarse `tolerance` input optional argument.
1052 1 : function test_getUpperGamma_3() result(assertion)
1053 1 : use Constants_mod, only: RK, IK
1054 : implicit none
1055 : integer(IK) :: i
1056 : logical :: assertion
1057 : integer(IK) , parameter :: ntest = 3
1058 : real(RK) , parameter :: tolerance = 1.e-3_RK
1059 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
1060 : real(RK) , parameter :: LowerLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
1061 : real(RK) , parameter :: UpperGamma_ref(ntest) = [0.367879441171442_RK, 0.815263244523772_RK, 1.828106329818355e-04_RK]
1062 : real(RK) :: UpperGamma(ntest)
1063 1 : real(RK) :: difference
1064 4 : do i = 1, ntest
1065 3 : UpperGamma(i) = getUpperGamma ( exponent = Exponent(i) &
1066 3 : , logGammaExponent = log_gamma(Exponent(i)) &
1067 3 : , lowerLim = LowerLim(i) &
1068 : , tolerance = tolerance &
1069 3 : )
1070 3 : difference = 2 * abs(UpperGamma(i) - UpperGamma_ref(i)) / UpperGamma_ref(i)
1071 3 : assertion = difference < tolerance
1072 4 : if (Test%isDebugMode .and. .not. assertion) then
1073 : ! LCOV_EXCL_START
1074 : write(Test%outputUnit,"(*(g0,:,', '))")
1075 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference UpperGamma, Computed UpperGamma, difference:"
1076 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), LowerLim(i), UpperGamma(i), UpperGamma_ref(i), difference
1077 : write(Test%outputUnit,"(*(g0,:,', '))")
1078 : end if
1079 : ! LCOV_EXCL_STOP
1080 : end do
1081 1 : end function test_getUpperGamma_3
1082 :
1083 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1084 :
1085 : !> \brief
1086 : !> Test [getUpperGamma](@ref math_mod::getuppergamma) without the `tolerance` input optional argument, in which case
1087 : !> the procedure should default to `epsilon` for the value of tolerance.
1088 1 : function test_getUpperGamma_4() result(assertion)
1089 1 : use Constants_mod, only: RK, IK
1090 : implicit none
1091 : integer(IK) :: i
1092 : logical :: assertion
1093 : integer(IK) , parameter :: ntest = 3
1094 : real(RK) , parameter :: tolerance = epsilon(1._RK)
1095 : real(RK) , parameter :: Exponent(ntest) = [1.0_RK, 5.0_RK, 0.5_RK]
1096 : real(RK) , parameter :: LowerLim(ntest) = [1.0_RK, 3.0_RK, 7.0_RK]
1097 : real(RK) , parameter :: UpperGamma_ref(ntest) = [0.367879441171442_RK, 0.815263244523772_RK, 1.828106329818355e-04_RK]
1098 : real(RK) :: UpperGamma(ntest)
1099 1 : real(RK) :: difference
1100 4 : do i = 1, ntest
1101 3 : UpperGamma(i) = getUpperGamma ( exponent = Exponent(i) &
1102 3 : , logGammaExponent = log_gamma(Exponent(i)) &
1103 3 : , lowerLim = LowerLim(i) &
1104 3 : )
1105 3 : difference = 2 * abs(UpperGamma(i) - UpperGamma_ref(i)) / UpperGamma_ref(i)
1106 3 : assertion = difference < 1000 * tolerance
1107 4 : if (Test%isDebugMode .and. .not. assertion) then
1108 : ! LCOV_EXCL_START
1109 : write(Test%outputUnit,"(*(g0,:,', '))")
1110 : write(Test%outputUnit,"(*(g0,:,', '))") "Exponent, UpperLim, Reference UpperGamma, Computed UpperGamma, difference:"
1111 : write(Test%outputUnit,"(*(g0,:,', '))") Exponent(i), LowerLim(i), UpperGamma(i), UpperGamma_ref(i), difference
1112 : write(Test%outputUnit,"(*(g0,:,', '))")
1113 : end if
1114 : ! LCOV_EXCL_STOP
1115 : end do
1116 1 : end function test_getUpperGamma_4
1117 :
1118 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1119 :
1120 : !> \brief
1121 : !> Test [getUpperGamma](@ref math_mod::getuppergamma) with a wrong positive value for the input argument `lowerLim`.
1122 1 : function test_getUpperGamma_5() result(assertion)
1123 1 : use Constants_mod, only: RK, IK, HUGE_RK
1124 : implicit none
1125 : logical :: assertion
1126 : real(RK) , parameter :: exponent = +1.0_RK
1127 : real(RK) , parameter :: lowerLim = -1.0_RK
1128 : real(RK) , parameter :: upperGamma_ref = -HUGE_RK
1129 1 : real(RK) :: upperGamma
1130 : upperGamma = getUpperGamma ( exponent = exponent &
1131 : , logGammaExponent = log_gamma(exponent) &
1132 : , lowerLim = lowerLim &
1133 1 : )
1134 1 : assertion = upperGamma == upperGamma_ref
1135 1 : if (Test%isDebugMode .and. .not. assertion) then
1136 : ! LCOV_EXCL_START
1137 : write(Test%outputUnit,"(*(g0,:,', '))")
1138 : write(Test%outputUnit,"(*(g0,:,', '))") "exponent ", exponent
1139 : write(Test%outputUnit,"(*(g0,:,', '))") "lowerLim ", lowerLim
1140 : write(Test%outputUnit,"(*(g0,:,', '))") "upperGamma ", upperGamma
1141 : write(Test%outputUnit,"(*(g0,:,', '))") "upperGamma_ref ", upperGamma_ref
1142 : write(Test%outputUnit,"(*(g0,:,', '))")
1143 : end if
1144 : ! LCOV_EXCL_STOP
1145 1 : end function test_getUpperGamma_5
1146 :
1147 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1148 :
1149 : !> \brief
1150 : !> Test [getGammaSeries](@ref math_mod::getgammaseries) with a zero value for the input argument `upperLim`.
1151 1 : function test_getGammaSeries_1() result(assertion)
1152 1 : use Constants_mod, only: RK, IK
1153 : implicit none
1154 : logical :: assertion
1155 : real(RK) , parameter :: exponent = +1.0_RK
1156 : real(RK) , parameter :: upperLim = 0._RK
1157 : real(RK) , parameter :: tolerance = 1.e-3_RK
1158 : real(RK) , parameter :: gammaSeries_ref = 0._RK
1159 1 : real(RK) :: gammaSeries
1160 : gammaSeries = getGammaSeries( exponent = exponent &
1161 : , logGammaExponent = log_gamma(exponent) &
1162 : , upperLim = upperLim &
1163 : , tolerance = tolerance &
1164 1 : )
1165 1 : assertion = gammaSeries == gammaSeries_ref
1166 1 : if (Test%isDebugMode .and. .not. assertion) then
1167 : ! LCOV_EXCL_START
1168 : write(Test%outputUnit,"(*(g0,:,', '))")
1169 : write(Test%outputUnit,"(*(g0,:,', '))") "exponent ", exponent
1170 : write(Test%outputUnit,"(*(g0,:,', '))") "upperLim ", upperLim
1171 : write(Test%outputUnit,"(*(g0,:,', '))") "tolerance ", tolerance
1172 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaSeries ", gammaSeries
1173 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaSeries_ref", gammaSeries_ref
1174 : write(Test%outputUnit,"(*(g0,:,', '))")
1175 : end if
1176 : ! LCOV_EXCL_STOP
1177 1 : end function test_getGammaSeries_1
1178 :
1179 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1180 :
1181 : !> \brief
1182 : !> Test [getGammaContFrac](@ref math_mod::getgammacontfrac) with a zero value for the input argument `lowerLim`.
1183 1 : function test_getGammaContFrac_1() result(assertion)
1184 1 : use Constants_mod, only: RK, IK
1185 : implicit none
1186 : logical :: assertion
1187 : real(RK) , parameter :: exponent = +1.0_RK
1188 : real(RK) , parameter :: lowerLim = 0._RK
1189 : real(RK) , parameter :: tolerance = 1.e-3_RK
1190 : real(RK) , parameter :: gammaContFrac_ref = 1._RK
1191 1 : real(RK) :: gammaContFrac
1192 : gammaContFrac = getGammaContFrac( exponent = exponent &
1193 : , logGammaExponent = log_gamma(exponent) &
1194 : , lowerLim = lowerLim &
1195 : , tolerance = tolerance &
1196 1 : )
1197 1 : assertion = gammaContFrac == gammaContFrac_ref
1198 1 : if (Test%isDebugMode .and. .not. assertion) then
1199 : ! LCOV_EXCL_START
1200 : write(Test%outputUnit,"(*(g0,:,', '))")
1201 : write(Test%outputUnit,"(*(g0,:,', '))") "exponent ", exponent
1202 : write(Test%outputUnit,"(*(g0,:,', '))") "lowerLim ", lowerLim
1203 : write(Test%outputUnit,"(*(g0,:,', '))") "tolerance ", tolerance
1204 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaContFrac ", gammaContFrac
1205 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaContFrac_ref ", gammaContFrac_ref
1206 : write(Test%outputUnit,"(*(g0,:,', '))")
1207 : end if
1208 : ! LCOV_EXCL_STOP
1209 1 : end function test_getGammaContFrac_1
1210 :
1211 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1212 :
1213 : !> \brief
1214 : !> Test the accuracy of [getLogGammaHalfInt](@ref math_mod::getloggammahalfint).
1215 1 : function test_getGammaHalfInt_1() result(assertion)
1216 1 : use Constants_mod, only: RK, IK
1217 : implicit none
1218 : logical :: assertion
1219 : real(RK) , parameter :: tolerance = 1.e-10_RK
1220 : real(RK) , parameter :: positiveHalfInteger = 0.5_RK
1221 : real(RK) , parameter :: gammaHalfInt_ref = gamma(positiveHalfInteger)
1222 1 : real(RK) :: gammaHalfInt
1223 1 : real(RK) :: difference
1224 1 : gammaHalfInt = getGammaHalfInt(positiveHalfInteger)
1225 1 : difference = abs(gammaHalfInt - gammaHalfInt_ref) / abs(gammaHalfInt_ref)
1226 1 : assertion = difference < tolerance
1227 1 : if (Test%isDebugMode .and. .not. assertion) then
1228 : ! LCOV_EXCL_START
1229 : write(Test%outputUnit,"(*(g0,:,', '))")
1230 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaHalfInt_ref ", gammaHalfInt_ref
1231 : write(Test%outputUnit,"(*(g0,:,', '))") "gammaHalfInt ", gammaHalfInt
1232 : write(Test%outputUnit,"(*(g0,:,', '))") "difference ", difference
1233 : write(Test%outputUnit,"(*(g0,:,', '))") "tolerance ", tolerance
1234 : write(Test%outputUnit,"(*(g0,:,', '))")
1235 : end if
1236 : ! LCOV_EXCL_STOP
1237 1 : end function test_getGammaHalfInt_1
1238 :
1239 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1240 :
1241 : !> \brief
1242 : !> Test the accuracy of [getLogGammaHalfInt](@ref math_mod::getloggammahalfint).
1243 1 : function test_getLogGammaHalfInt_1() result(assertion)
1244 1 : use Constants_mod, only: RK, IK
1245 : implicit none
1246 : logical :: assertion
1247 : real(RK) , parameter :: tolerance = 1.e-10_RK
1248 : real(RK) , parameter :: positiveHalfInteger = 0.5_RK
1249 : real(RK) , parameter :: logGammaHalfInt_ref = log_gamma(positiveHalfInteger)
1250 1 : real(RK) :: logGammaHalfInt
1251 1 : real(RK) :: difference
1252 1 : logGammaHalfInt = getLogGammaHalfInt(positiveHalfInteger)
1253 1 : difference = abs(logGammaHalfInt - logGammaHalfInt_ref) / abs(logGammaHalfInt_ref)
1254 1 : assertion = difference < tolerance
1255 1 : if (Test%isDebugMode .and. .not. assertion) then
1256 : ! LCOV_EXCL_START
1257 : write(Test%outputUnit,"(*(g0,:,', '))")
1258 : write(Test%outputUnit,"(*(g0,:,', '))") "logGammaHalfInt_ref", logGammaHalfInt_ref
1259 : write(Test%outputUnit,"(*(g0,:,', '))") "logGammaHalfInt ", logGammaHalfInt
1260 : write(Test%outputUnit,"(*(g0,:,', '))") "difference ", difference
1261 : write(Test%outputUnit,"(*(g0,:,', '))") "tolerance ", tolerance
1262 : write(Test%outputUnit,"(*(g0,:,', '))")
1263 : end if
1264 : ! LCOV_EXCL_STOP
1265 1 : end function test_getLogGammaHalfInt_1
1266 :
1267 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1268 :
1269 : end module Test_Math_mod ! LCOV_EXCL_LINE
|
