ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
test_pm_quantile.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
19
21
23 use pm_err, only: err_type
24 use pm_test, only: test_type, LK
25 implicit none
26
27 private
28 public :: setTest
29 type(test_type) :: test
30
31 integer(IK) , parameter :: lenRnd = 50_IK
32
33 !real(RK) :: UnifRnd(lenRnd) = [ 0.0759666916908419_RK &
34 ! , 0.2399161535536580_RK &
35 ! , 0.1233189348351660_RK &
36 ! , 0.1839077882824170_RK &
37 ! , 0.2399525256649030_RK &
38 ! , 0.4172670690843700_RK &
39 ! , 0.0496544303257421_RK &
40 ! , 0.9027161099152810_RK &
41 ! , 0.9447871897216460_RK &
42 ! , 0.4908640924680800_RK &
43 ! , 0.4892526384000190_RK &
44 ! , 0.3377194098213770_RK &
45 ! , 0.9000538464176620_RK &
46 ! , 0.3692467811202150_RK &
47 ! , 0.1112027552937870_RK &
48 ! , 0.7802520683211380_RK &
49 ! , 0.3897388369612530_RK &
50 ! , 0.2416912859138330_RK &
51 ! , 0.4039121455881150_RK &
52 ! , 0.0964545251683886_RK &
53 ! , 0.1319732926063350_RK &
54 ! , 0.9420505907754850_RK &
55 ! , 0.9561345402298020_RK &
56 ! , 0.5752085950784660_RK &
57 ! , 0.0597795429471558_RK &
58 ! , 0.2347799133724060_RK &
59 ! , 0.3531585712220710_RK &
60 ! , 0.8211940401979590_RK &
61 ! , 0.0154034376515551_RK &
62 ! , 0.0430238016578078_RK &
63 ! , 0.1689900294627040_RK &
64 ! , 0.6491154749564520_RK &
65 ! , 0.7317223856586700_RK &
66 ! , 0.6477459631363070_RK &
67 ! , 0.4509237064309450_RK &
68 ! , 0.5470088922863450_RK &
69 ! , 0.2963208056077730_RK &
70 ! , 0.7446928070741560_RK &
71 ! , 0.1889550150325450_RK &
72 ! , 0.6867754333653150_RK &
73 ! , 0.1835111557372700_RK &
74 ! , 0.3684845964903370_RK &
75 ! , 0.6256185607296900_RK &
76 ! , 0.7802274351513770_RK &
77 ! , 0.0811257688657853_RK &
78 ! , 0.9293859709687300_RK &
79 ! , 0.7757126786084020_RK &
80 ! , 0.4867916324031720_RK &
81 ! , 0.4358585885809190_RK &
82 ! , 0.4467837494298060_RK ]
83 !
84 real(RK) :: StdNormRnd1(lenRnd) = [ 0.537667139546100_RK &
85 , 1.83388501459509_RK &
86 , -2.25884686100365_RK &
87 , 0.862173320368121_RK &
88 , 0.318765239858981_RK &
89 , -1.30768829630527_RK &
90 , -0.433592022305684_RK &
91 , 0.342624466538650_RK &
92 , 3.57839693972576_RK &
93 , 2.76943702988488_RK &
94 , -1.34988694015652_RK &
95 , 3.03492346633185_RK &
96 , 0.725404224946106_RK &
97 , -0.0630548731896562_RK &
98 , 0.714742903826096_RK &
99 , -0.204966058299775_RK &
100 , -0.124144348216312_RK &
101 , 1.48969760778547_RK &
102 , 1.40903448980048_RK &
103 , 1.41719241342961_RK &
104 , 0.671497133608081_RK &
105 , -1.20748692268504_RK &
106 , 0.717238651328839_RK &
107 , 1.63023528916473_RK &
108 , 0.488893770311789_RK &
109 , 1.03469300991786_RK &
110 , 0.726885133383238_RK &
111 , -0.303440924786016_RK &
112 , 0.293871467096658_RK &
113 , -0.787282803758638_RK &
114 , 0.888395631757642_RK &
115 , -1.14707010696915_RK &
116 , -1.06887045816803_RK &
117 , -0.809498694424876_RK &
118 , -2.94428416199490_RK &
119 , 1.43838029281510_RK &
120 , 0.325190539456198_RK &
121 , -0.754928319169703_RK &
122 , 1.37029854009523_RK &
123 , -1.71151641885370_RK &
124 , -0.102242446085491_RK &
125 , -0.241447041607358_RK &
126 , 0.319206739165502_RK &
127 , 0.312858596637428_RK &
128 , -0.864879917324457_RK &
129 , -0.0300512961962686_RK &
130 , -0.164879019209038_RK &
131 , 0.627707287528727_RK &
132 , 1.09326566903948_RK &
133 , 1.10927329761440_RK ]
134
135 !real(RK) :: StdNormRnd2(lenRnd) = [ -0.863652821988714_RK &
136 ! , 0.0773590911304249_RK &
137 ! , -1.21411704361541_RK &
138 ! , -1.11350074148676_RK &
139 ! , -0.00684932810334806_RK &
140 ! , 1.53263030828475_RK &
141 ! , -0.769665913753682_RK &
142 ! , 0.371378812760058_RK &
143 ! , -0.225584402271252_RK &
144 ! , 1.11735613881447_RK &
145 ! , -1.08906429505224_RK &
146 ! , 0.0325574641649735_RK &
147 ! , 0.552527021112224_RK &
148 ! , 1.10061021788087_RK &
149 ! , 1.54421189550395_RK &
150 ! , 0.0859311331754255_RK &
151 ! , -1.49159031063761_RK &
152 ! , -0.742301837259857_RK &
153 ! , -1.06158173331999_RK &
154 ! , 2.35045722400204_RK &
155 ! , -0.615601881466894_RK &
156 ! , 0.748076783703985_RK &
157 ! , -0.192418510588264_RK &
158 ! , 0.888610425420721_RK &
159 ! , -0.764849236567874_RK &
160 ! , -1.40226896933876_RK &
161 ! , -1.42237592509150_RK &
162 ! , 0.488193909859941_RK &
163 ! , -0.177375156618825_RK &
164 ! , -0.196053487807333_RK &
165 ! , 1.41931015064255_RK &
166 ! , 0.291584373984183_RK &
167 ! , 0.197811053464361_RK &
168 ! , 1.58769908997406_RK &
169 ! , -0.804465956349547_RK &
170 ! , 0.696624415849607_RK &
171 ! , 0.835088165072682_RK &
172 ! , -0.243715140377952_RK &
173 ! , 0.215670086403744_RK &
174 ! , -1.16584393148205_RK &
175 ! , -1.14795277889859_RK &
176 ! , 0.104874716016494_RK &
177 ! , 0.722254032225002_RK &
178 ! , 2.58549125261624_RK &
179 ! , -0.666890670701386_RK &
180 ! , 0.187331024578940_RK &
181 ! , -0.0824944253709554_RK &
182 ! , -1.93302291785099_RK &
183 ! , -0.438966153934773_RK &
184 ! , -1.79467884145512_RK ]
185
186 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
187
188 integer, parameter :: ndata = 50
189 !integer(IK), parameter :: DataUnsorted_IK(ndata)= &
190 ! [ 1201_IK &
191 ! , 1187_IK &
192 ! , 1188_IK &
193 ! , 1193_IK &
194 ! , 1177_IK &
195 ! , 1153_IK &
196 ! , 1134_IK &
197 ! , 1146_IK &
198 ! , 1172_IK &
199 ! , 1181_IK &
200 ! , 1197_IK &
201 ! , 1172_IK &
202 ! , 1141_IK &
203 ! , 1216_IK &
204 ! , 1158_IK &
205 ! , 1174_IK &
206 ! , 1189_IK &
207 ! , 1211_IK &
208 ! , 1157_IK &
209 ! , 1184_IK &
210 ! , 1177_IK &
211 ! , 1157_IK &
212 ! , 1191_IK &
213 ! , 1176_IK &
214 ! , 1196_IK &
215 ! , 1150_IK &
216 ! , 1185_IK &
217 ! , 1190_IK &
218 ! , 1172_IK &
219 ! , 1161_IK &
220 ! , 1179_IK &
221 ! , 1189_IK &
222 ! , 1136_IK &
223 ! , 1148_IK &
224 ! , 1176_IK &
225 ! , 1142_IK &
226 ! , 1146_IK &
227 ! , 1202_IK &
228 ! , 1156_IK &
229 ! , 1170_IK &
230 ! , 1146_IK &
231 ! , 1170_IK &
232 ! , 1169_IK &
233 ! , 1211_IK &
234 ! , 1168_IK &
235 ! , 1189_IK &
236 ! , 1170_IK &
237 ! , 1162_IK &
238 ! , 1167_IK &
239 ! , 1180_IK ]
240 real(RK), parameter :: DataUnsorted_RK(ndata) = &
241 [ 5.28935260000000_RK &
242 , 5.50145870000000_RK &
243 , 5.89022390000000_RK &
244 , 5.06549460000000_RK &
245 , 5.62128260000000_RK &
246 , 4.49246930000000_RK &
247 , 3.54559920000000_RK &
248 , 4.17171310000000_RK &
249 , 5.34432780000000_RK &
250 , 4.30855910000000_RK &
251 , 6.12466330000000_RK &
252 , 4.45103540000000_RK &
253 , 4.08259680000000_RK &
254 , 7.64761290000000_RK &
255 , 6.53095480000000_RK &
256 , 6.07550490000000_RK &
257 , 7.32100850000000_RK &
258 , 5.82501650000000_RK &
259 , 4.19347540000000_RK &
260 , 4.89687790000000_RK &
261 , 5.61290890000000_RK &
262 , 5.70994940000000_RK &
263 , 5.00047920000000_RK &
264 , 5.47741520000000_RK &
265 , 4.99151560000000_RK &
266 , 5.08172850000000_RK &
267 , 5.98773500000000_RK &
268 , 6.97849360000000_RK &
269 , 6.91612860000000_RK &
270 , 4.90595890000000_RK &
271 , 5.71852950000000_RK &
272 , 4.12146660000000_RK &
273 , 5.51241440000000_RK &
274 , 5.26293780000000_RK &
275 , 5.14932990000000_RK &
276 , 4.14738170000000_RK &
277 , 5.55786790000000_RK &
278 , 7.08800450000000_RK &
279 , 6.08987380000000_RK &
280 , 4.73697940000000_RK &
281 , 3.80934450000000_RK &
282 , 6.03942270000000_RK &
283 , 5.96600840000000_RK &
284 , 6.06674510000000_RK &
285 , 5.84361600000000_RK &
286 , 6.19013970000000_RK &
287 , 4.43891700000000_RK &
288 , 4.45833300000000_RK &
289 , 5.47659170000000_RK &
290 , 4.65761920000000_RK ]
291
292 !integer(IK) , parameter :: DataUnsorted2_IK(ndata) = DataUnsorted_IK
293 !real(RK) , parameter :: DataUnsorted2_RK(ndata) = DataUnsorted_RK
294
295!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
296
297contains
298
299!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
300
301 subroutine setTest()
302
304 call test%run(test_getQuantile_1, SK_"test_getQuantile_1")
305 call test%run(test_getQuantile_2, SK_"test_getQuantile_2")
306 call test%run(test_getMedian_RK_1, SK_"test_getMedian_RK_1")
307 call test%run(test_getMedian_RK_2, SK_"test_getMedian_RK_2")
308 call test%summarize()
309
310 end subroutine setTest
311
312!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
313
314 function test_getQuantile_1() result(assertion)
315
316 use pm_kind, only: RK, IK
317 use pm_arrayVerbose, only: getVerbose
318
319 implicit none
320
321 logical(LK) :: assertion
322 integer(IK) , parameter :: nq = 9_IK
323 integer(IK) , parameter :: Weight(lenRnd) = [ 2_IK, 2_IK, 1_IK, 0_IK, 0_IK, 1_IK, 2_IK, 2_IK, 1_IK, 0_IK, 0_IK, 0_IK, 1_IK, 2_IK &
324 , 2_IK, 1_IK, 0_IK, 0_IK, 1_IK, 2_IK, 2_IK, 1_IK, 0_IK, 0_IK, 1_IK, 2_IK, 2_IK, 1_IK &
325 , 0_IK, 0_IK, 1_IK, 2_IK, 2_IK, 2_IK, 1_IK, 0_IK, 0_IK, 1_IK, 2_IK, 2_IK, 1_IK, 0_IK &
326 , 0_IK, 1_IK, 2_IK, 2_IK, 1_IK, 0_IK, 0_IK, 1_IK ]
327 real(RK) , parameter :: SortedQuantileProbability(nq) = [0._RK, 0.05_RK, 0.1_RK, 0.25_RK, 0.5_RK, 0.75_RK, 0.9_RK, 0.95_RK, 1._RK]
328 real(RK) :: Quantile_ref(nq)
329 real(RK) :: Quantile(nq)
330 real(RK) :: Difference(nq)
331
332 ! First, getVerbose the input `Point` vector to compute the reference quantiles.
333
334 call getQuantile( SortedQuantileProbability = SortedQuantileProbability & ! LCOV_EXCL_LINE
335 , Point = getVerbose(StdNormRnd1, Weight) & ! LCOV_EXCL_LINE
336 , Quantile = Quantile_ref & ! LCOV_EXCL_LINE
337 )
338 if (assertion) return
339 assertion = .not. assertion
340
341 ! Now computed the quantiles using the input weighted `Point`.
342
343 call getQuantile( SortedQuantileProbability = SortedQuantileProbability & ! LCOV_EXCL_LINE
344 , Point = StdNormRnd1 & ! LCOV_EXCL_LINE
345 , Weight = Weight & ! LCOV_EXCL_LINE
346 , sumWeight = sum(Weight) & ! LCOV_EXCL_LINE
347 , Quantile = Quantile & ! LCOV_EXCL_LINE
348 )
349 if (assertion) return
350 assertion = .not. assertion
351
352 Difference = (Quantile - Quantile_ref)
353 assertion = all(Difference == 0._RK)
354
355 if (test%traceable .and. .not. assertion) then
356 ! LCOV_EXCL_START
357 write(test%disp%unit,"(*(g0,:,', '))")
358 write(test%disp%unit,"(*(g0,:,', '))") "Quantile_ref ", Quantile_ref
359 write(test%disp%unit,"(*(g0,:,', '))") "Quantile ", Quantile
360 write(test%disp%unit,"(*(g0,:,', '))") "difference ", Difference
361 write(test%disp%unit,"(*(g0,:,', '))")
362 ! LCOV_EXCL_STOP
363 end if
364
365 end function test_getQuantile_1
366
367!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
368
369 function test_getQuantile_2() result(assertion)
370
371 use pm_kind, only: RK, IK
372
373 implicit none
374
375 logical(LK) :: assertion
376 integer(IK) , parameter :: nq = 9_IK
377 real(RK) , parameter :: SortedQuantileProbability(nq) = [0._RK, 0.05_RK, 0.1_RK, 0.25_RK, 0.5_RK, 0.75_RK, 0.9_RK, 0.95_RK, 1._RK]
378 real(RK) , parameter :: Quantile_ref(nq) = [ -2.944284161994900_RK &
379 , -1.711516418853700_RK &
380 , -1.307688296305270_RK &
381 , -.4335920223056840_RK &
382 , +.3192067391655020_RK &
383 , +1.034693009917860_RK &
384 , +1.489697607785470_RK &
385 , +2.769437029884880_RK &
386 , +3.578396939725760_RK ]
387 real(RK) :: Quantile(nq)
388 real(RK) :: Difference(nq)
389
390 call getQuantile( SortedQuantileProbability = SortedQuantileProbability & ! LCOV_EXCL_LINE
391 , Point = StdNormRnd1 & ! LCOV_EXCL_LINE
392 , Quantile = Quantile & ! LCOV_EXCL_LINE
393 )
394
395 Difference = (Quantile - Quantile_ref)
396 assertion = all(Difference == 0._RK)
397
398 if (test%traceable .and. .not. assertion) then
399 ! LCOV_EXCL_START
400 write(test%disp%unit,"(*(g0,:,', '))")
401 write(test%disp%unit,"(*(g0,:,', '))") "Quantile_ref ", Quantile_ref
402 write(test%disp%unit,"(*(g0,:,', '))") "Quantile ", Quantile
403 write(test%disp%unit,"(*(g0,:,', '))") "difference ", Difference
404 write(test%disp%unit,"(*(g0,:,', '))")
405 end if
406 ! LCOV_EXCL_STOP
407
408 end function test_getQuantile_2
409
410!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
411
412 function test_getMedian_RK_1() result(assertion)
413
414 use pm_err, only: err_type
415 implicit none
416
417 logical(LK) :: assertion
418 real(RK) , parameter :: median_ref = 5.477003450000000_RK
419 real(RK) :: median
420
421 call getMedian(lenArray = ndata, Array = DataUnsorted_RK, median = median)
422
423 assertion = median == median_ref
424
425 if (test%traceable .and. .not. assertion) then
426 ! LCOV_EXCL_START
427 write(test%disp%unit,"(*(g0,:,' '))")
428 write(test%disp%unit,"(*(g0,:,' '))") "median_ref =", median_ref
429 write(test%disp%unit,"(*(g0,:,' '))") "median =", median
430 write(test%disp%unit,"(*(g0,:,' '))")
431 ! LCOV_EXCL_STOP
432 end if
433
434 end function test_getMedian_RK_1
435
436!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
437
438 function test_getMedian_RK_2() result(assertion)
439
440 use pm_err, only: err_type
441 implicit none
442
443 logical(LK) :: assertion
444 real(RK) , parameter :: median_ref = 5.477415200000000_RK
445 real(RK) :: median
446
447 call getMedian_RK(lenArray = ndata-1, Array = DataUnsorted_RK(1:ndata-1), median = median)
448
449 assertion = median == median_ref
450
451 if (test%traceable .and. .not. assertion) then
452 ! LCOV_EXCL_START
453 write(test%disp%unit,"(*(g0,:,' '))")
454 write(test%disp%unit,"(*(g0,:,' '))") "median_ref =", median_ref
455 write(test%disp%unit,"(*(g0,:,' '))") "median =", median
456 write(test%disp%unit,"(*(g0,:,' '))")
457 ! LCOV_EXCL_STOP
458 end if
459
460 end function test_getMedian_RK_2
461
462!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
463
464end module test_pm_sampleQuan
Generate an equally-weighted (verbose or flattened) array of the input weighted array of rank 1 or 2.
This module contains procedures and generic interfaces for flattening (duplicating the elements of) a...
This module contains classes and procedures for reporting and handling errors.
Definition: pm_err.F90:52
character(*, SK), parameter MODULE_NAME
Definition: pm_err.F90:58
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter RK
The default real kind in the ParaMonte library: real64 in Fortran, c_double in C-Fortran Interoperati...
Definition: pm_kind.F90:543
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
This module contains procedures and data types for computing sample quantile.
This module contains a simple unit-testing framework for the Fortran libraries, including the ParaMon...
Definition: pm_test.F90:42
This module contains tests of the module pm_sampleQuan.
type(test_type) test
logical(LK) function test_getMedian_RK_1()
logical(LK) function test_getMedian_RK_2()
integer, parameter ndata
real(RK), dimension(ndata), parameter DataUnsorted_RK
real(RK), dimension(lenRnd) StdNormRnd1
logical(LK) function test_getQuantile_1()
integer(IK), parameter lenRnd
logical(LK) function test_getQuantile_2()
This is the derived type for generating objects to gracefully and verbosely handle runtime unexpected...
Definition: pm_err.F90:157
This is the derived type test_type for generating objects that facilitate testing of a series of proc...
Definition: pm_test.F90:209