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 [Integration_mod](@ref integration_mod).
44 : !> \author Amir Shahmoradi
45 :
46 : module Test_Integration_mod
47 :
48 : use Integration_mod
49 : !use QuadPackDouble_mod, only: dqagi
50 : use Constants_mod, only: IK, RK
51 : use Err_mod, only: Err_type
52 : use Test_mod, only: Test_type
53 : implicit none
54 :
55 : private
56 : public :: test_Integration
57 :
58 : type(Test_type) :: Test
59 :
60 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
61 :
62 : contains
63 :
64 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
65 :
66 1 : subroutine test_Integration()
67 : implicit none
68 1 : Test = Test_type(moduleName=MODULE_NAME)
69 1 : call Test%run(test_doQuadRombOpen_1, "test_doQuadRombOpen_1")
70 1 : call Test%run(test_doQuadRombOpen_2, "test_doQuadRombOpen_2")
71 1 : call Test%run(test_doQuadRombOpen_3, "test_doQuadRombOpen_3")
72 1 : call Test%run(test_doQuadRombOpen_4, "test_doQuadRombOpen_4")
73 1 : call Test%run(test_doQuadRombClosed_1, "test_doQuadRombClosed_1")
74 1 : call Test%finalize()
75 1 : end subroutine test_Integration
76 :
77 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
78 :
79 : ! integrate the Gaussian PDF via the midexp routine on the open interval \f$[0, +\infty)\f$.
80 1 : function test_doQuadRombOpen_1() result(assertion)
81 :
82 1 : use Constants_mod, only: RK, IK, POSINF_RK
83 : use Integration_mod
84 :
85 : implicit none
86 :
87 : logical :: assertion
88 1 : real(RK) :: integral, relativeError, difference
89 : real(RK), parameter :: integral_ref = 0.5_RK
90 : real(RK), parameter :: tolerance = 1.e-10_RK
91 : integer(IK) :: numFuncEval, ierr
92 :
93 : call doQuadRombOpen ( getFunc = getTestFuncOpenInterval_1 &
94 : , integrate = midexp &
95 : , lowerLim = 0._RK &
96 : , upperLim = POSINF_RK & ! do not set this to huge() as GNU Fortran in debug mode crashes
97 : , maxRelativeError = 0.1*tolerance &
98 : , nRefinement = 10_IK &
99 : , integral = integral &
100 : , relativeError = relativeError &
101 : , numFuncEval = numFuncEval &
102 : , ierr = ierr &
103 1 : )
104 :
105 1 : assertion = ierr == 0_IK
106 :
107 1 : if (Test%isDebugMode .and. .not. assertion) then
108 : ! LCOV_EXCL_START
109 : if (Test%isDebugMode .and. .not. assertion) then
110 : write(Test%outputUnit,"(*(g0))") "ierr = ", ierr, " /= 0"
111 : end if
112 : return
113 : ! LCOV_EXCL_STOP
114 : end if
115 :
116 1 : difference = abs(integral - integral_ref) / abs(integral_ref)
117 1 : assertion = difference < tolerance
118 1 : assertion = assertion .and. relativeError <= tolerance
119 :
120 : ! LCOV_EXCL_START
121 : if (Test%isDebugMode .and. .not. assertion) then
122 : write(Test%outputUnit,"(*(g0))")
123 : write(Test%outputUnit,"(*(g0))") "integral_ref = ", integral_ref
124 : write(Test%outputUnit,"(*(g0))") "integral = ", integral
125 : write(Test%outputUnit,"(*(g0))") "difference = ", difference
126 : write(Test%outputUnit,"(*(g0))") "relativeError = ", relativeError
127 : write(Test%outputUnit,"(*(g0))")
128 : end if
129 : ! LCOV_EXCL_STOP
130 :
131 1 : end function test_doQuadRombOpen_1
132 :
133 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
134 :
135 1 : function test_doQuadRombOpen_2() result(assertion)
136 :
137 1 : use Constants_mod, only: RK, IK, POSINF_RK
138 : use Integration_mod
139 :
140 : implicit none
141 :
142 : logical :: assertion
143 1 : real(RK) :: relativeError, difference
144 1 : real(RK) :: integral
145 : real(RK), parameter :: integral_ref = 0.078649603525143_RK
146 : real(RK), parameter :: tolerance = 1.e-10_RK
147 : integer(IK) :: numFuncEval, ierr
148 :
149 : call doQuadRombOpen ( getFunc = getTestFuncOpenInterval_1 &
150 : , integrate = midinf &
151 : , lowerLim = 1._RK &
152 : , upperLim = POSINF_RK & ! do not set this to huge() as GNU Fortran in debug mode crashes
153 : , maxRelativeError = 0.1*tolerance &
154 : , nRefinement = 10_IK &
155 : , integral = integral &
156 : , relativeError = relativeError &
157 : , numFuncEval = numFuncEval &
158 : , ierr = ierr &
159 1 : )
160 :
161 1 : assertion = ierr == 0_IK
162 :
163 1 : if (Test%isDebugMode .and. .not. assertion) then
164 : ! LCOV_EXCL_START
165 : if (Test%isDebugMode .and. .not. assertion) then
166 : write(Test%outputUnit,"(*(g0))") "ierr = ", ierr, " /= 0"
167 : end if
168 : return
169 : ! LCOV_EXCL_STOP
170 : end if
171 :
172 1 : difference = abs(integral - integral_ref) / abs(integral_ref)
173 1 : assertion = difference < tolerance
174 1 : assertion = assertion .and. relativeError <= tolerance
175 :
176 1 : if (Test%isDebugMode .and. .not. assertion) then
177 : ! LCOV_EXCL_START
178 : write(Test%outputUnit,"(*(g0))")
179 : write(Test%outputUnit,"(*(g0))") "integral_ref = ", integral_ref
180 : write(Test%outputUnit,"(*(g0))") "integral = ", integral
181 : write(Test%outputUnit,"(*(g0))") "difference = ", difference
182 : write(Test%outputUnit,"(*(g0))") "relativeError = ", relativeError
183 : write(Test%outputUnit,"(*(g0))")
184 : end if
185 : ! LCOV_EXCL_STOP
186 :
187 1 : end function test_doQuadRombOpen_2
188 :
189 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
190 :
191 1 : function test_doQuadRombOpen_3() result(assertion)
192 :
193 1 : use Constants_mod, only: RK, IK, POSINF_RK
194 : use Integration_mod
195 :
196 : implicit none
197 :
198 : logical :: assertion
199 1 : real(RK) :: relativeError, difference
200 1 : real(RK) :: integral
201 : real(RK), parameter :: integral_ref = 0.078649603525143_RK
202 : real(RK), parameter :: tolerance = 1.e-4_RK
203 : integer(IK) :: numFuncEval, ierr
204 :
205 : call doQuadRombOpen ( getFunc = getTestFuncOpenInterval_1 &
206 : , integrate = midinf &
207 : , lowerLim = POSINF_RK & ! do not set this to huge() as GNU Fortran in debug mode crashes
208 : , upperLim = -1._RK &
209 : , maxRelativeError = 1.e-10_RK &
210 : , nRefinement = 10_IK &
211 : , integral = integral &
212 : , relativeError = relativeError &
213 : , numFuncEval = numFuncEval &
214 : , ierr = ierr &
215 1 : )
216 :
217 1 : assertion = ierr == 0_IK
218 1 : if (Test%isDebugMode .and. .not. assertion) then
219 : ! LCOV_EXCL_START
220 : if (Test%isDebugMode .and. .not. assertion) then
221 : write(Test%outputUnit,"(*(g0))") "ierr = ", ierr, " /= 0"
222 : end if
223 : return
224 : ! LCOV_EXCL_STOP
225 : end if
226 :
227 1 : difference = abs(integral - integral_ref) / abs(integral_ref)
228 1 : assertion = difference < tolerance
229 1 : assertion = assertion .and. relativeError <= tolerance
230 :
231 1 : if (Test%isDebugMode .and. .not. assertion) then
232 : ! LCOV_EXCL_START
233 : write(Test%outputUnit,"(*(g0))")
234 : write(Test%outputUnit,"(*(g0))") "integral_ref = ", integral_ref
235 : write(Test%outputUnit,"(*(g0))") "integral = ", integral
236 : write(Test%outputUnit,"(*(g0))") "difference = ", difference
237 : write(Test%outputUnit,"(*(g0))") "relativeError = ", relativeError
238 : write(Test%outputUnit,"(*(g0))")
239 : end if
240 : ! LCOV_EXCL_STOP
241 :
242 1 : end function test_doQuadRombOpen_3
243 :
244 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
245 :
246 1 : function test_doQuadRombOpen_4() result(assertion)
247 :
248 1 : use Constants_mod, only: RK, IK
249 : use Integration_mod
250 :
251 : implicit none
252 :
253 : logical :: assertion
254 1 : real(RK) :: relativeError, difference
255 1 : real(RK) :: integral
256 : real(RK), parameter :: integral_ref = 1._RK
257 : real(RK), parameter :: tolerance = 1.e-10_RK
258 : integer(IK) :: numFuncEval, ierr
259 :
260 : call doQuadRombOpen ( getFunc = getTestFuncOpenInterval_1 &
261 : , integrate = midpnt &
262 : , lowerLim = -5._RK &
263 : , upperLim = +5._RK &
264 : , maxRelativeError = 1.e-10_RK &
265 : , nRefinement = 10_IK &
266 : , integral = integral &
267 : , relativeError = relativeError &
268 : , numFuncEval = numFuncEval &
269 : , ierr = ierr &
270 1 : )
271 :
272 1 : assertion = ierr == 0_IK
273 1 : if (Test%isDebugMode .and. .not. assertion) then
274 : ! LCOV_EXCL_START
275 : if (Test%isDebugMode .and. .not. assertion) then
276 : write(Test%outputUnit,"(*(g0))") "ierr = ", ierr, " /= 0"
277 : end if
278 : return
279 : ! LCOV_EXCL_STOP
280 : end if
281 :
282 1 : difference = abs(integral - integral_ref) / abs(integral_ref)
283 1 : assertion = difference < tolerance
284 1 : assertion = assertion .and. relativeError <= tolerance
285 :
286 1 : if (Test%isDebugMode .and. .not. assertion) then
287 : ! LCOV_EXCL_START
288 : write(Test%outputUnit,"(*(g0))")
289 : write(Test%outputUnit,"(*(g0))") "integral_ref = ", integral_ref
290 : write(Test%outputUnit,"(*(g0))") "integral = ", integral
291 : write(Test%outputUnit,"(*(g0))") "difference = ", difference
292 : write(Test%outputUnit,"(*(g0))") "relativeError = ", relativeError
293 : write(Test%outputUnit,"(*(g0))")
294 : end if
295 : ! LCOV_EXCL_STOP
296 :
297 1 : end function test_doQuadRombOpen_4
298 :
299 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
300 :
301 1 : function test_doQuadRombClosed_1() result(assertion)
302 :
303 1 : use Constants_mod, only: RK, IK
304 : use Integration_mod
305 :
306 : implicit none
307 :
308 : logical :: assertion
309 1 : real(RK) :: relativeError, difference
310 1 : real(RK) :: integral
311 : real(RK), parameter :: integral_ref = 1._RK
312 : real(RK), parameter :: tolerance = 1.e-10_RK
313 : integer(IK) :: numFuncEval, ierr
314 :
315 : call doQuadRombClosed ( getFunc = getTestFuncOpenInterval_1 &
316 : , lowerLim = -5._RK &
317 : , upperLim = +5._RK &
318 : , maxRelativeError = 0.1*tolerance &
319 : , nRefinement = 10_IK &
320 : , integral = integral &
321 : , relativeError = relativeError &
322 : , numFuncEval = numFuncEval &
323 : , ierr = ierr &
324 1 : )
325 1 : assertion = ierr == 0_IK
326 1 : if (Test%isDebugMode .and. .not. assertion) then
327 : ! LCOV_EXCL_START
328 : if (Test%isDebugMode .and. .not. assertion) then
329 : write(Test%outputUnit,"(*(g0))") "ierr = ", ierr, " /= 0"
330 : end if
331 : return
332 : ! LCOV_EXCL_STOP
333 : end if
334 :
335 1 : difference = abs(integral - integral_ref) / abs(integral_ref)
336 1 : assertion = difference < tolerance
337 1 : assertion = assertion .and. relativeError <= tolerance
338 :
339 1 : if (Test%isDebugMode .and. .not. assertion) then
340 : ! LCOV_EXCL_START
341 : write(Test%outputUnit,"(*(g0))")
342 : write(Test%outputUnit,"(*(g0))") "integral_ref = ", integral_ref
343 : write(Test%outputUnit,"(*(g0))") "integral = ", integral
344 : write(Test%outputUnit,"(*(g0))") "difference = ", difference
345 : write(Test%outputUnit,"(*(g0))") "relativeError = ", relativeError
346 : write(Test%outputUnit,"(*(g0))")
347 : end if
348 : ! LCOV_EXCL_STOP
349 :
350 1 : end function test_doQuadRombClosed_1
351 :
352 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
353 :
354 : ! The Gaussian PDF with sigma = 1/sqrt(2)
355 79245 : function getTestFuncOpenInterval_1(x) result(funcVal)
356 1 : use Constants_mod, only: RK, INVSQRTPI, NEGLOGINF_RK
357 : implicit none
358 : real(RK), intent(in) :: x
359 : real(RK) :: funcVal
360 : real(RK), parameter :: LOGINVSQRTPI = log(INVSQRTPI)
361 79245 : funcVal = LOGINVSQRTPI - x**2
362 79245 : if (funcVal<NEGLOGINF_RK) then ! This takes care of the GNU Fortran 9.1 test crash in debug mode.
363 1480 : funcVal = 0._RK
364 : else
365 77765 : funcVal = exp(funcVal)
366 : end if
367 79245 : end function getTestFuncOpenInterval_1
368 :
369 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
370 :
371 : ! requires `midsql`.
372 : ! The MATLAB test function: https://www.mathworks.com/help/matlab/ref/integral.html
373 : ! LCOV_EXCL_START
374 : function getTestFuncOpenIntervalMATLAB(x) result(funcVal)
375 : use Constants_mod, only: RK
376 : implicit none
377 : real(RK), intent(in) :: x
378 : real(RK) :: funcVal
379 : funcVal = exp(-x**2) * log(x)**2
380 : end function getTestFuncOpenIntervalMATLAB
381 : ! LCOV_EXCL_STOP
382 :
383 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
384 :
385 : end module Test_Integration_mod ! LCOV_EXCL_LINE
|
