Line data Source code
1 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3 : !!!! !!!!
4 : !!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5 : !!!! !!!!
6 : !!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7 : !!!! !!!!
8 : !!!! This file is part of the ParaMonte library. !!!!
9 : !!!! !!!!
10 : !!!! LICENSE !!!!
11 : !!!! !!!!
12 : !!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13 : !!!! !!!!
14 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16 :
17 : !> \brief
18 : !> This file contains procedure implementations of tests of [pm_sampleVar](@ref pm_sampleVar).
19 : !>
20 : !> \fintest
21 : !>
22 : !> \author
23 : !> \FatemehBagheri, Wednesday 5:03 PM, August 11, 2021, Dallas, TX
24 :
25 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
26 :
27 : real(TKC), parameter :: rtol = epsilon(1._TKC) * 100
28 : ! Define the sample type.
29 : #if CK_ENABLED
30 : #define TYPE_OF_SAMPLE complex(TKC)
31 : complex(TKC), parameter :: ONE = (1._TKC, 1._TKC), TWO = 2 * (1._TKC, 1._TKC), ctol = (rtol, rtol)
32 : #elif RK_ENABLED
33 : #define TYPE_OF_SAMPLE real(TKC)
34 : real(TKC), parameter :: ONE = 1._TKC, TWO = 2._TKC, ctol = rtol
35 : #else
36 : #error "Unrecognized interface."
37 : #endif
38 : !%%%%%%%%%%%%%%%%%%%%%%%
39 : #if getVarCorrection_ENABLED
40 : !%%%%%%%%%%%%%%%%%%%%%%%
41 :
42 : integer(IK) :: itry
43 : real(TKC) :: weisum, weisqs, correction, correction_ref
44 44 : do itry = 1, 10
45 :
46 : ! Test reliability weight.
47 :
48 40 : weisum = getUnifRand(2._TKC, 3._TKC)
49 40 : weisqs = getUnifRand(3._TKC, 4._TKC)
50 40 : correction = getVarCorrection(weisum, weisqs)
51 40 : correction_ref = weisum**2 / (weisum**2 - weisqs)
52 40 : assertion = logical(correction == correction_ref, LK)
53 40 : call report(__LINE__, SK_"The Bessel Reliability correction must be correctly computed.")
54 :
55 40 : weisum = getUnifRand(2._TKC, 3._TKC)
56 40 : correction = getVarCorrection(weisum)
57 40 : correction_ref = weisum / (weisum - 1._TKC)
58 40 : assertion = logical(correction == correction_ref, LK)
59 44 : call report(__LINE__, SK_"The Bessel Frequency correction must be correctly computed.")
60 :
61 : end do
62 :
63 : contains
64 :
65 80 : subroutine report(line, msg)
66 : integer, intent(in) :: line
67 : character(*, SK), intent(in) :: msg
68 80 : if (test%traceable .and. .not. assertion) then
69 : ! LCOV_EXCL_START
70 : call test%disp%skip()
71 : call test%disp%show("[weisum, weisqs, correction, correction_ref]")
72 : call test%disp%show( [weisum, weisqs, correction, correction_ref] )
73 : call test%disp%skip()
74 : ! LCOV_EXCL_STOP
75 : end if
76 80 : call test%assert(assertion, msg, int(line, IK))
77 80 : end subroutine
78 :
79 : !%%%%%%%%%%%%%
80 : #elif getVar_ENABLED
81 : !%%%%%%%%%%%%%
82 :
83 : real(TKC) :: rweisum
84 : integer(IK) :: iweisum
85 : integer(IK) :: itry, nsam, ndim, dim
86 : real(TKC), allocatable :: rweight(:)
87 : integer(IK), allocatable :: iweight(:)
88 : TYPE_OF_SAMPLE, allocatable :: sample(:,:), mean(:)
89 : real(TKC), allocatable :: var(:), var_ref(:), vdiff(:)
90 8 : assertion = .true._LK
91 :
92 408 : do itry = 1, 50
93 :
94 400 : dim = 0
95 400 : nsam = getUnifRand(2_IK, 7_IK)
96 :
97 : ! test sample D2 DIM interface.
98 :
99 : block
100 :
101 400 : dim = getChoice([1, 2])
102 400 : ndim = getUnifRand(1_IK, 5_IK)
103 400 : call setResized(var_ref, ndim)
104 400 : call setResized(var, ndim)
105 400 : if (dim == 2) then
106 4197 : sample = getUnifRand(ONE, TWO, ndim, nsam)
107 : else
108 3482 : sample = getUnifRand(ONE, TWO, nsam, ndim)
109 : end if
110 2583 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
111 2583 : iweight = getUnifRand(1_IK, 9_IK, nsam)
112 2183 : iweisum = sum(iweight)
113 2183 : rweisum = sum(rweight)
114 :
115 : ! integer weighted
116 :
117 1991 : var = getVar(sample, dim, iweight)
118 1991 : mean = getMean(sample, dim, iweight)
119 400 : call setVar(var_ref, mean, sample, dim, iweight, iweisum)
120 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
121 :
122 3182 : var = getVar(sample, dim, iweight, fweight_type())
123 1591 : var_ref = var_ref * getVarCorrection(real(iweisum, TKC))
124 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
125 :
126 : ! real weighted
127 :
128 1991 : var = getVar(sample, dim, rweight)
129 1991 : mean = getMean(sample, dim, rweight)
130 400 : call setVar(var_ref, mean, sample, dim, rweight, rweisum)
131 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
132 :
133 3182 : var = getVar(sample, dim, rweight, rweight_type())
134 3374 : var_ref = var_ref * getVarCorrection(rweisum, sum(rweight**2))
135 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
136 :
137 : ! unweighted
138 :
139 1991 : var = getVar(sample, dim)
140 1991 : mean = getMean(sample, dim)
141 400 : call setVar(var_ref, mean, sample, dim)
142 400 : call setAssertedVar(__LINE__, SK_"unweighted")
143 :
144 3182 : var = getVar(sample, dim, fweight_type())
145 1591 : var_ref = var_ref * getVarCorrection(real(size(sample, dim), TKC))
146 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
147 :
148 : end block
149 :
150 : ! test sample D2 ALL interface.
151 :
152 : block
153 :
154 400 : dim = 0
155 400 : ndim = getUnifRand(1_IK, 5_IK)
156 400 : call setResized(var_ref, 1_IK)
157 400 : call setResized(mean, 1_IK)
158 400 : call setResized(var, 1_IK)
159 400 : if (getUnifRand()) then
160 3820 : sample = getUnifRand(ONE, TWO, ndim, nsam)
161 : else
162 3593 : sample = getUnifRand(ONE, TWO, nsam, ndim)
163 : end if
164 5935 : rweight = getUnifRand(1._TKC, 9._TKC, nsam * ndim)
165 5935 : iweight = getUnifRand(1_IK, 9_IK, nsam * ndim)
166 5535 : iweisum = sum(iweight)
167 5535 : rweisum = sum(rweight)
168 :
169 : ! integer weighted
170 :
171 400 : var(1) = getVar(sample, iweight)
172 400 : mean(1) = getMean(sample, iweight)
173 400 : call setVar(var_ref(1), mean(1), sample, iweight, iweisum)
174 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
175 :
176 400 : var(1) = getVar(sample, iweight, fweight_type())
177 800 : var_ref = var_ref * getVarCorrection(real(iweisum, TKC))
178 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
179 :
180 : ! real weighted
181 :
182 400 : var(1) = getVar(sample, rweight)
183 400 : mean(1) = getMean(sample, rweight)
184 400 : call setVar(var_ref(1), mean(1), sample, rweight, rweisum)
185 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
186 :
187 400 : var(1) = getVar(sample, rweight, rweight_type())
188 5935 : var_ref = var_ref * getVarCorrection(rweisum, sum(rweight**2))
189 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
190 :
191 : ! unweighted
192 :
193 400 : var(1) = getVar(sample)
194 400 : mean(1) = getMean(sample)
195 400 : call setVar(var_ref(1), mean(1), sample)
196 400 : call setAssertedVar(__LINE__, SK_"unweighted")
197 :
198 400 : var(1) = getVar(sample)
199 800 : var_ref = var_ref * getVarCorrection(real(size(sample), TKC))
200 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
201 :
202 : end block
203 :
204 : ! test sample D1 DIM interface.
205 :
206 : block
207 :
208 400 : dim = 1_IK
209 400 : ndim = 1_IK
210 400 : call setResized(var, 1_IK)
211 400 : call setResized(var_ref, 1_IK)
212 2983 : sample = getUnifRand(ONE, TWO, nsam, ndim)
213 2583 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
214 2583 : iweight = getUnifRand(1_IK, 9_IK, nsam)
215 2183 : iweisum = sum(iweight)
216 2183 : rweisum = sum(rweight)
217 :
218 : ! integer weighted
219 :
220 400 : var(1) = getVar(sample(:,1), dim, iweight)
221 400 : mean(1) = getMean(sample(:,1), dim, iweight)
222 400 : call setVar(var_ref(1), mean(1), sample(:,1), dim, iweight, iweisum)
223 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
224 :
225 400 : var(1) = getVar(sample(:,1), dim, iweight, fweight_type())
226 800 : var_ref = var_ref * getVarCorrection(real(iweisum, TKC))
227 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
228 :
229 : ! real weighted
230 :
231 400 : var(1) = getVar(sample(:,1), dim, rweight)
232 400 : mean(1) = getMean(sample(:,1), dim, rweight)
233 400 : call setVar(var_ref(1), mean(1), sample(:,1), dim, rweight, rweisum)
234 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
235 :
236 400 : var(1) = getVar(sample(:,1), dim, rweight, rweight_type())
237 2583 : var_ref = var_ref * getVarCorrection(rweisum, sum(rweight**2))
238 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
239 :
240 : ! unweighted
241 :
242 400 : var(1) = getVar(sample(:,1), dim)
243 400 : mean(1) = getMean(sample(:,1), dim)
244 400 : call setVar(var_ref(1), mean(1), sample(:,1), dim)
245 400 : call setAssertedVar(__LINE__, SK_"unweighted")
246 :
247 400 : var(1) = getVar(sample(:,1), dim)
248 800 : var_ref = var_ref * getVarCorrection(real(size(sample), TKC))
249 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
250 :
251 : end block
252 :
253 : ! test sample D1 ALL interface.
254 :
255 8 : block
256 :
257 400 : dim = 1_IK
258 400 : ndim = 1_IK
259 400 : call setResized(var, 1_IK)
260 400 : call setResized(var_ref, 1_IK)
261 2983 : sample = getUnifRand(ONE, TWO, nsam, ndim)
262 2583 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
263 2583 : iweight = getUnifRand(1_IK, 9_IK, nsam)
264 2183 : iweisum = sum(iweight)
265 2183 : rweisum = sum(rweight)
266 :
267 : ! integer weighted
268 :
269 400 : var(1) = getVar(sample(:,1), iweight)
270 400 : mean(1) = getMean(sample(:,1), iweight)
271 400 : call setVar(var_ref(1), mean(1), sample(:,1), iweight, iweisum)
272 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
273 :
274 400 : var(1) = getVar(sample(:,1), iweight, fweight_type())
275 800 : var_ref = var_ref * getVarCorrection(real(iweisum, TKC))
276 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
277 :
278 : ! real weighted
279 :
280 400 : var(1) = getVar(sample(:,1), rweight)
281 400 : mean(1) = getMean(sample(:,1), rweight)
282 400 : call setVar(var_ref(1), mean(1), sample(:,1), rweight, rweisum)
283 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
284 :
285 400 : var(1) = getVar(sample(:,1), rweight, rweight_type())
286 2583 : var_ref = var_ref * getVarCorrection(rweisum, sum(rweight**2))
287 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
288 :
289 : ! unweighted
290 :
291 400 : var(1) = getVar(sample(:,1))
292 400 : mean(1) = getMean(sample(:,1))
293 400 : call setVar(var_ref(1), mean(1), sample(:,1))
294 400 : call setAssertedVar(__LINE__, SK_"unweighted")
295 :
296 400 : var(1) = getVar(sample(:,1))
297 800 : var_ref = var_ref * getVarCorrection(real(size(sample), TKC))
298 400 : call report(__LINE__, SK_"The unbiased `var` must be computed correctly.")
299 :
300 : end block
301 :
302 : end do
303 :
304 : contains
305 :
306 4800 : subroutine setAssertedVar(line, this)
307 : character(*, SK), intent(in) :: this
308 : integer, intent(in) :: line
309 16773 : vdiff = abs(var - var_ref)
310 11973 : assertion = assertion .and. all(vdiff < rtol)
311 4800 : call report(line, SK_"The `var` must be computed correctly for "//this//SK_" sample.")
312 4800 : end subroutine
313 :
314 9600 : subroutine report(line, msg)
315 : integer, intent(in) :: line
316 : character(*, SK), intent(in) :: msg
317 9600 : if (test%traceable .and. .not. assertion) then
318 : ! LCOV_EXCL_START
319 : call test%disp%skip()
320 : call test%disp%show("[ndim, nsam, dim, shape(sample, IK)]")
321 : call test%disp%show( [ndim, nsam, dim, shape(sample, IK)] )
322 : call test%disp%show("sample")
323 : call test%disp%show( sample )
324 : call test%disp%show("iweight")
325 : call test%disp%show( iweight )
326 : call test%disp%show("rweight")
327 : call test%disp%show( rweight )
328 : call test%disp%show("var_ref")
329 : call test%disp%show( var_ref )
330 : call test%disp%show("var")
331 : call test%disp%show( var )
332 : call test%disp%show("vdiff")
333 : call test%disp%show( vdiff )
334 : call test%disp%skip()
335 : ! LCOV_EXCL_STOP
336 : end if
337 9600 : call test%assert(assertion, msg, int(line, IK))
338 9600 : end subroutine
339 :
340 : !%%%%%%%%%%%%%
341 : #elif setVar_ENABLED
342 : !%%%%%%%%%%%%%
343 :
344 : real(TKC) :: rweisum
345 : integer(IK) :: iweisum
346 : real(TKC), allocatable :: rweight(:)
347 : integer(IK), allocatable :: iweight(:)
348 : integer(IK) :: itry, nsam, ndim, dim
349 : TYPE_OF_SAMPLE, allocatable :: sample(:,:), mean(:)
350 : real(TKC), allocatable :: var(:), var_ref(:), vdiff(:)
351 8 : assertion = .true._LK
352 :
353 408 : do itry = 1, 50
354 :
355 400 : dim = 0
356 400 : nsam = getUnifRand(2_IK, 7_IK)
357 :
358 : ! test sample D2 DIM interface.
359 :
360 : block
361 :
362 400 : dim = getChoice([1, 2])
363 400 : ndim = getUnifRand(1_IK, 5_IK)
364 400 : call setResized(var, ndim)
365 400 : if (dim == 2) then
366 4161 : sample = getUnifRand(ONE, TWO, ndim, nsam)
367 : else
368 3458 : sample = getUnifRand(ONE, TWO, nsam, ndim)
369 : end if
370 2585 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
371 2585 : iweight = getUnifRand(1_IK, 9_IK, nsam)
372 2185 : iweisum = sum(iweight)
373 2185 : rweisum = sum(rweight)
374 :
375 : ! integer weighted
376 :
377 1979 : mean = getMean(sample, dim, iweight)
378 3764 : var_ref = getVarD2(mean, sample, dim, real(iweight, TKC))
379 :
380 400 : call setVar(var, mean, sample, dim, iweight, iweisum)
381 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
382 :
383 1579 : call setVar(var, getShifted(sample, dim, -mean), dim, iweight, iweisum)
384 400 : call setAssertedVar(__LINE__, SK_"integer-weighted shifted")
385 :
386 : ! real weighted
387 :
388 1979 : mean = getMean(sample, dim, rweight)
389 1979 : var_ref = getVarD2(mean, sample, dim, rweight)
390 :
391 400 : call setVar(var, mean, sample, dim, rweight, rweisum)
392 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
393 :
394 1579 : call setVar(var, getShifted(sample, dim, -mean), dim, rweight, rweisum)
395 400 : call setAssertedVar(__LINE__, SK_"real-weighted shifted")
396 :
397 : ! unweighted
398 :
399 1979 : mean = getMean(sample, dim)
400 1979 : var_ref = getVarD2(mean, sample, dim)
401 :
402 400 : call setVar(var, mean, sample, dim)
403 400 : call setAssertedVar(__LINE__, SK_"unweighted")
404 :
405 1579 : call setVar(var, getShifted(sample, dim, -mean), dim)
406 400 : call setAssertedVar(__LINE__, SK_"unweighted shifted")
407 :
408 : end block
409 :
410 : ! test sample D2 ALL interface.
411 :
412 : block
413 :
414 400 : dim = 0
415 400 : ndim = getUnifRand(1_IK, 5_IK)
416 400 : call setResized(var_ref, 1_IK)
417 400 : call setResized(mean, 1_IK)
418 400 : call setResized(var, 1_IK)
419 400 : if (getUnifRand()) then
420 3741 : sample = getUnifRand(ONE, TWO, ndim, nsam)
421 : else
422 4233 : sample = getUnifRand(ONE, TWO, nsam, ndim)
423 : end if
424 6416 : rweight = getUnifRand(1._TKC, 9._TKC, nsam * ndim)
425 6416 : iweight = getUnifRand(1_IK, 9_IK, nsam * ndim)
426 6016 : iweisum = sum(iweight)
427 6016 : rweisum = sum(rweight)
428 :
429 : ! integer weighted
430 :
431 800 : mean = getMean(sample, iweight)
432 6816 : var_ref = getVarD2(mean, sample, weight = real(iweight, TKC))
433 :
434 400 : call setVar(var(1), mean(1), sample, iweight, iweisum)
435 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
436 :
437 400 : call setVar(var(1), getShifted(sample, -mean(1)), iweight, iweisum)
438 400 : call setAssertedVar(__LINE__, SK_"integer-weighted shifted")
439 :
440 : ! real weighted
441 :
442 400 : mean(1) = getMean(sample, rweight)
443 1200 : var_ref = getVarD2(mean, sample, weight = rweight)
444 :
445 400 : call setVar(var(1), mean(1), sample, rweight, rweisum)
446 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
447 :
448 400 : call setVar(var(1), getShifted(sample, -mean(1)), rweight, rweisum)
449 400 : call setAssertedVar(__LINE__, SK_"real-weighted shifted")
450 :
451 : ! unweighted
452 :
453 400 : mean(1) = getMean(sample)
454 1200 : var_ref = getVarD2(mean, sample)
455 :
456 400 : call setVar(var(1), mean(1), sample)
457 400 : call setAssertedVar(__LINE__, SK_"unweighted")
458 :
459 400 : call setVar(var(1), getShifted(sample, -mean(1)))
460 400 : call setAssertedVar(__LINE__, SK_"unweighted shifted")
461 :
462 : end block
463 :
464 : ! test sample D1 DIM interface.
465 :
466 : block
467 :
468 400 : dim = 1_IK
469 400 : ndim = 1_IK
470 400 : call setResized(var, 1_IK)
471 2985 : sample = getUnifRand(ONE, TWO, nsam, ndim)
472 2585 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
473 2585 : iweight = getUnifRand(1_IK, 9_IK, nsam)
474 2185 : iweisum = sum(iweight)
475 2185 : rweisum = sum(rweight)
476 :
477 : ! integer weighted
478 :
479 800 : mean = getMean(sample(:,1), dim, iweight)
480 2585 : var_ref = getVarD1(mean(1), sample(:,1), real(iweight, TKC))
481 :
482 400 : call setVar(var(1), mean(1), sample(:,1), dim, iweight, iweisum)
483 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
484 :
485 400 : call setVar(var(1), getShifted(sample(:,1), dim, -mean(1)), dim, iweight, iweisum)
486 400 : call setAssertedVar(__LINE__, SK_"integer-weighted shifted")
487 :
488 : ! real weighted
489 :
490 800 : mean = getMean(sample(:,1), rweight)
491 800 : var_ref = getVarD1(mean(1), sample(:,1), weight = rweight)
492 :
493 400 : call setVar(var(1), mean(1), sample(:,1), dim, rweight, rweisum)
494 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
495 :
496 400 : call setVar(var(1), getShifted(sample(:,1), dim, -mean(1)), dim, rweight, rweisum)
497 400 : call setAssertedVar(__LINE__, SK_"real-weighted shifted")
498 :
499 : ! unweighted
500 :
501 800 : mean = getMean(sample(:,1))
502 800 : var_ref = getVarD1(mean(1), sample(:,1))
503 :
504 400 : call setVar(var(1), mean(1), sample(:,1), dim)
505 400 : call setAssertedVar(__LINE__, SK_"unweighted")
506 :
507 400 : call setVar(var(1), getShifted(sample(:,1), dim, -mean(1)), dim)
508 400 : call setAssertedVar(__LINE__, SK_"unweighted shifted")
509 :
510 : end block
511 :
512 : ! test sample D1 ALL interface.
513 :
514 8 : block
515 :
516 400 : dim = 1_IK
517 400 : ndim = 1_IK
518 400 : call setResized(var, 1_IK)
519 2985 : sample = getUnifRand(ONE, TWO, nsam, ndim)
520 2585 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
521 2585 : iweight = getUnifRand(1_IK, 9_IK, nsam)
522 2185 : iweisum = sum(iweight)
523 2185 : rweisum = sum(rweight)
524 :
525 : ! integer weighted
526 :
527 800 : mean = getMean(sample(:,1), iweight)
528 2585 : var_ref = getVarD1(mean(1), sample(:,1), real(iweight, TKC))
529 :
530 400 : call setVar(var(1), mean(1), sample(:,1), iweight, iweisum)
531 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
532 :
533 400 : call setVar(var(1), getShifted(sample(:,1), -mean(1)), iweight, iweisum)
534 400 : call setAssertedVar(__LINE__, SK_"integer-weighted shifted")
535 :
536 : ! real weighted
537 :
538 800 : mean = getMean(sample(:,1), rweight)
539 800 : var_ref = getVarD1(mean(1), sample(:,1), weight = rweight)
540 :
541 400 : call setVar(var(1), mean(1), sample(:,1), rweight, rweisum)
542 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
543 :
544 400 : call setVar(var(1), getShifted(sample(:,1), -mean(1)), rweight, rweisum)
545 400 : call setAssertedVar(__LINE__, SK_"real-weighted shifted")
546 :
547 : ! unweighted
548 :
549 800 : mean = getMean(sample(:,1))
550 800 : var_ref = getVarD1(mean(1), sample(:,1))
551 400 : call setVar(var(1), mean(1), sample(:,1))
552 400 : call setAssertedVar(__LINE__, SK_"unweighted")
553 :
554 400 : call setVar(var(1), getShifted(sample(:,1), -mean(1)))
555 400 : call setAssertedVar(__LINE__, SK_"unweighted shifted")
556 :
557 : end block
558 :
559 : end do
560 :
561 : contains
562 :
563 9600 : subroutine setAssertedVar(line, this)
564 : character(*, SK), intent(in) :: this
565 : integer, intent(in) :: line
566 33474 : vdiff = abs(var - var_ref)
567 23874 : assertion = assertion .and. all(vdiff < rtol)
568 9600 : call report(line, SK_"The `var` must be computed correctly for "//this//SK_" sample.")
569 9600 : end subroutine
570 :
571 9600 : subroutine report(line, msg)
572 : integer, intent(in) :: line
573 : character(*, SK), intent(in) :: msg
574 9600 : if (test%traceable .and. .not. assertion) then
575 : ! LCOV_EXCL_START
576 : call test%disp%skip()
577 : call test%disp%show("[ndim, nsam, dim, shape(sample, IK)]")
578 : call test%disp%show( [ndim, nsam, dim, shape(sample, IK)] )
579 : call test%disp%show("sample")
580 : call test%disp%show( sample )
581 : call test%disp%show("iweight")
582 : call test%disp%show( iweight )
583 : call test%disp%show("rweight")
584 : call test%disp%show( rweight )
585 : call test%disp%show("var_ref")
586 : call test%disp%show( var_ref )
587 : call test%disp%show("var")
588 : call test%disp%show( var )
589 : call test%disp%show("vdiff")
590 : call test%disp%show( vdiff )
591 : call test%disp%skip()
592 : ! LCOV_EXCL_STOP
593 : end if
594 9600 : call test%assert(assertion, msg, int(line, IK))
595 9600 : end subroutine
596 :
597 5937 : PURE function getVarD1(mean, sample, weight) result(var)
598 : real(TKC), intent(in), optional :: weight(:)
599 : TYPE_OF_SAMPLE, intent(in) :: sample(:), mean
600 11874 : TYPE_OF_SAMPLE :: sampleShifted(size(sample, 1, IK))
601 : real(TKC) :: var
602 32613 : sampleShifted = getShifted(sample, -mean)
603 5937 : if (present(weight)) then
604 43484 : var = real(dot_product(sampleShifted, sampleShifted * weight), TKC) / sum(weight)
605 : else
606 10871 : var = real(dot_product(sampleShifted, sampleShifted), TKC) / size(sample, 1, IK)
607 : end if
608 5937 : end function
609 :
610 2400 : PURE function getVarD2(mean, sample, dim, weight) result(var)
611 : real(TKC), intent(in), optional :: weight(:)
612 : TYPE_OF_SAMPLE, intent(in) :: sample(:,:), mean(:)
613 : integer(IK), intent(in), optional :: dim
614 : real(TKC), allocatable :: var(:)
615 : integer(IK) :: idim
616 2400 : if (present(dim)) then
617 1200 : allocate(var(size(sample, 3 - dim, IK)))
618 4737 : do idim = 1, size(sample, 3 - dim, IK)
619 4737 : if (dim == 2) then
620 2460 : var(idim) = getVarD1(mean(idim), sample(idim,:), weight)
621 : else
622 2256 : var(idim) = getVarD1(mean(idim), sample(:,idim), weight)
623 : end if
624 : end do
625 : else
626 : block
627 1200 : TYPE_OF_SAMPLE :: sampleShifted(size(sample, kind = IK))
628 19248 : sampleShifted = reshape(sample, shape(sampleShifted)) - mean(1)
629 1200 : if (present(weight)) then
630 25664 : var = [real(dot_product(sampleShifted, sampleShifted * weight), TKC) / sum(weight)]
631 : else
632 6816 : var = [real(dot_product(sampleShifted, sampleShifted), TKC) / size(sampleShifted, 1, IK)]
633 : end if
634 : end block
635 : end if
636 2400 : end function
637 :
638 : !%%%%%%%%%%%%%%%%%
639 : #elif setVarMean_ENABLED
640 : !%%%%%%%%%%%%%%%%%
641 :
642 : real(TKC) :: rweisum, rweisum_ref
643 : integer(IK) :: iweisum, iweisum_ref
644 : real(TKC), allocatable :: rweight(:)
645 : integer(IK), allocatable :: iweight(:)
646 : integer(IK) :: itry, nsam, ndim, dim
647 : TYPE_OF_SAMPLE, allocatable :: sample(:,:), mean(:), mean_ref(:), meang(:), mdiff(:)
648 : real(TKC), allocatable :: var(:), var_ref(:), vdiff(:)
649 8 : assertion = .true._LK
650 :
651 408 : do itry = 1, 50
652 :
653 400 : dim = 0
654 400 : nsam = getUnifRand(2_IK, 7_IK)
655 :
656 : ! test sample D2 DIM interface.
657 :
658 : block
659 :
660 400 : dim = getChoice([1, 2])
661 400 : ndim = getUnifRand(1_IK, 5_IK)
662 400 : call setResized(mean_ref, ndim)
663 400 : call setResized(var_ref, ndim)
664 400 : call setResized(mean, ndim)
665 400 : call setResized(var, ndim)
666 400 : if (dim == 2) then
667 3919 : sample = getUnifRand(ONE, TWO, ndim, nsam)
668 993 : meang = sample(:,1)
669 : else
670 3755 : sample = getUnifRand(ONE, TWO, nsam, ndim)
671 1024 : meang = sample(1,:)
672 : end if
673 2570 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
674 2570 : iweight = getUnifRand(1_IK, 9_IK, nsam)
675 2170 : iweisum_ref = sum(iweight)
676 2170 : rweisum_ref = sum(rweight)
677 400 : iweisum = 0
678 400 : rweisum = 0
679 :
680 : ! integer weighted
681 :
682 2017 : var_ref = getVar(sample, dim, iweight)
683 2017 : mean_ref = getMean(sample, dim, iweight)
684 : call setVarMean(var, mean, sample, dim, iweight, iweisum, meang)
685 400 : call setAssertedSum(__LINE__, SK_"integer-weighted", real(iweisum, TKC), real(iweisum_ref, TKC))
686 400 : call setAssertedAvg(__LINE__, SK_"integer-weighted")
687 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
688 :
689 : ! real weighted
690 :
691 2017 : var_ref = getVar(sample, dim, rweight)
692 2017 : mean_ref = getMean(sample, dim, rweight)
693 200 : call setVarMean(var, mean, sample, dim, rweight, rweisum, meang)
694 400 : call setAssertedSum(__LINE__, SK_"real-weighted", rweisum, rweisum_ref)
695 400 : call setAssertedAvg(__LINE__, SK_"real-weighted")
696 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
697 :
698 : ! unweighted
699 :
700 2017 : var_ref = getVar(sample, dim)
701 2017 : mean_ref = getMean(sample, dim)
702 400 : call setVarMean(var, mean, sample, dim, meang)
703 400 : call setAssertedAvg(__LINE__, SK_"unweighted")
704 400 : call setAssertedVar(__LINE__, SK_"unweighted")
705 :
706 : end block
707 :
708 : ! test sample D2 ALL interface.
709 :
710 : block
711 :
712 400 : dim = 0
713 400 : ndim = getUnifRand(1_IK, 5_IK)
714 400 : call setResized(mean_ref, 1_IK)
715 400 : call setResized(var_ref, 1_IK)
716 400 : call setResized(mean, 1_IK)
717 400 : call setResized(var, 1_IK)
718 400 : if (getUnifRand()) then
719 3645 : sample = getUnifRand(ONE, TWO, ndim, nsam)
720 929 : meang = sample(:,1)
721 : else
722 3910 : sample = getUnifRand(ONE, TWO, nsam, ndim)
723 1067 : meang = sample(1,:)
724 : end if
725 6096 : rweight = getUnifRand(1._TKC, 9._TKC, nsam * ndim)
726 6096 : iweight = getUnifRand(1_IK, 9_IK, nsam * ndim)
727 5696 : iweisum_ref = sum(iweight)
728 5696 : rweisum_ref = sum(rweight)
729 400 : iweisum = 0
730 400 : rweisum = 0
731 :
732 : ! integer weighted
733 :
734 800 : var_ref = getVar(sample, iweight)
735 800 : mean_ref = getMean(sample, iweight)
736 400 : call setVarMean(var(1), mean(1), sample, iweight, iweisum, meang(1))
737 400 : call setAssertedSum(__LINE__, SK_"integer-weighted", real(iweisum, TKC), real(iweisum_ref, TKC))
738 400 : call setAssertedAvg(__LINE__, SK_"integer-weighted")
739 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
740 :
741 : ! real weighted
742 :
743 800 : var_ref = getVar(sample, rweight)
744 800 : mean_ref = getMean(sample, rweight)
745 400 : call setVarMean(var(1), mean(1), sample, rweight, rweisum, meang(1))
746 400 : call setAssertedSum(__LINE__, SK_"real-weighted", rweisum, rweisum_ref)
747 400 : call setAssertedAvg(__LINE__, SK_"real-weighted")
748 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
749 :
750 : ! unweighted
751 :
752 800 : var_ref = getVar(sample)
753 800 : mean_ref = getMean(sample)
754 400 : call setVarMean(var(1), mean(1), sample, meang(1))
755 400 : call setAssertedAvg(__LINE__, SK_"unweighted")
756 400 : call setAssertedVar(__LINE__, SK_"unweighted")
757 :
758 : end block
759 :
760 : ! test sample D1 DIM interface.
761 :
762 : block
763 :
764 400 : dim = 1_IK
765 400 : ndim = 1_IK
766 400 : call setResized(var, 1_IK)
767 400 : call setResized(mean, 1_IK)
768 400 : call setResized(var_ref, 1_IK)
769 400 : call setResized(mean_ref, 1_IK)
770 2970 : sample = getUnifRand(ONE, TWO, nsam, ndim)
771 1200 : meang = sample(1,:)
772 2570 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
773 2570 : iweight = getUnifRand(1_IK, 9_IK, nsam)
774 2170 : iweisum_ref = sum(iweight)
775 2170 : rweisum_ref = sum(rweight)
776 400 : iweisum = 0
777 400 : rweisum = 0
778 :
779 : ! integer weighted
780 :
781 800 : var_ref = getVar(sample(:,1), dim, iweight)
782 800 : mean_ref = getMean(sample(:,1), dim, iweight)
783 400 : call setVarMean(var(1), mean(1), sample(:,1), dim, iweight, iweisum, meang(1))
784 400 : call setAssertedSum(__LINE__, SK_"integer-weighted", real(iweisum, TKC), real(iweisum_ref, TKC))
785 400 : call setAssertedAvg(__LINE__, SK_"integer-weighted")
786 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
787 :
788 : ! real weighted
789 :
790 800 : var_ref = getVar(sample(:,1), dim, rweight)
791 800 : mean_ref = getMean(sample(:,1), dim, rweight)
792 400 : call setVarMean(var(1), mean(1), sample(:,1), dim, rweight, rweisum, meang(1))
793 400 : call setAssertedSum(__LINE__, SK_"real-weighted", rweisum, rweisum_ref)
794 400 : call setAssertedAvg(__LINE__, SK_"real-weighted")
795 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
796 :
797 : ! unweighted
798 :
799 800 : var_ref = getVar(sample(:,1), dim)
800 800 : mean_ref = getMean(sample(:,1), dim)
801 400 : call setVarMean(var(1), mean(1), sample(:,1), dim, meang(1))
802 400 : call setAssertedAvg(__LINE__, SK_"unweighted")
803 400 : call setAssertedVar(__LINE__, SK_"unweighted")
804 :
805 : end block
806 :
807 : ! test sample D1 ALL interface.
808 :
809 1608 : block
810 :
811 400 : dim = 0_IK
812 400 : ndim = 1_IK
813 400 : call setResized(var, 1_IK)
814 400 : call setResized(mean, 1_IK)
815 400 : call setResized(var_ref, 1_IK)
816 400 : call setResized(mean_ref, 1_IK)
817 2970 : sample = getUnifRand(ONE, TWO, nsam, ndim)
818 1200 : meang = sample(1,:)
819 2570 : rweight = getUnifRand(1._TKC, 9._TKC, nsam)
820 2570 : iweight = getUnifRand(1_IK, 9_IK, nsam)
821 2170 : iweisum_ref = sum(iweight)
822 2170 : rweisum_ref = sum(rweight)
823 400 : iweisum = 0
824 400 : rweisum = 0
825 :
826 : ! integer weighted
827 :
828 800 : var_ref = getVar(sample(:,1), iweight)
829 800 : mean_ref = getMean(sample(:,1), iweight)
830 400 : call setVarMean(var(1), mean(1), sample(:,1), iweight, iweisum, meang(1))
831 400 : call setAssertedSum(__LINE__, SK_"integer-weighted", real(iweisum, TKC), real(iweisum_ref, TKC))
832 400 : call setAssertedAvg(__LINE__, SK_"integer-weighted")
833 400 : call setAssertedVar(__LINE__, SK_"integer-weighted")
834 :
835 : ! real weighted
836 :
837 800 : var_ref = getVar(sample(:,1), rweight)
838 800 : mean_ref = getMean(sample(:,1), rweight)
839 400 : call setVarMean(var(1), mean(1), sample(:,1), rweight, rweisum, meang(1))
840 400 : call setAssertedSum(__LINE__, SK_"real-weighted", rweisum, rweisum_ref)
841 400 : call setAssertedAvg(__LINE__, SK_"real-weighted")
842 400 : call setAssertedVar(__LINE__, SK_"real-weighted")
843 :
844 : ! unweighted
845 :
846 800 : var_ref = getVar(sample(:,1))
847 800 : mean_ref = getMean(sample(:,1))
848 400 : call setVarMean(var(1), mean(1), sample(:,1), meang(1))
849 400 : call setAssertedAvg(__LINE__, SK_"unweighted")
850 400 : call setAssertedVar(__LINE__, SK_"unweighted")
851 :
852 : end block
853 :
854 : end do
855 :
856 : contains
857 :
858 3200 : subroutine setAssertedSum(line, this, weisum, weisum_ref)
859 : real(TKC), intent(in) :: weisum, weisum_ref
860 : character(*, SK), intent(in) :: this
861 : integer, intent(in) :: line
862 3200 : assertion = assertion .and. abs(weisum - weisum_ref) < rtol * weisum_ref
863 3200 : call report(line, SK_"The `weisum` must be computed correctly for "//this//SK_" sample.")
864 3200 : end subroutine
865 :
866 4800 : subroutine setAssertedAvg(line, this)
867 : character(*, SK), intent(in) :: this
868 : integer, intent(in) :: line
869 20517 : mdiff = abs(mean - mean_ref)
870 12051 : assertion = assertion .and. all(mdiff < ctol)
871 4800 : call report(line, SK_"The `mean` must be computed correctly for "//this//SK_" sample.")
872 4800 : end subroutine
873 :
874 4800 : subroutine setAssertedVar(line, this)
875 : character(*, SK), intent(in) :: this
876 : integer, intent(in) :: line
877 16851 : vdiff = abs(var - var_ref)
878 12051 : assertion = assertion .and. all(vdiff < rtol)
879 4800 : call report(line, SK_"The `var` must be computed correctly for "//this//SK_" sample.")
880 4800 : end subroutine
881 :
882 12800 : subroutine report(line, msg)
883 : integer, intent(in) :: line
884 : character(*, SK), intent(in) :: msg
885 12800 : if (test%traceable .and. .not. assertion) then
886 : ! LCOV_EXCL_START
887 : call test%disp%skip()
888 : call test%disp%show("[ndim, nsam, dim, shape(sample, IK)]")
889 : call test%disp%show( [ndim, nsam, dim, shape(sample, IK)] )
890 : call test%disp%show("sample")
891 : call test%disp%show( sample )
892 : call test%disp%show("meang")
893 : call test%disp%show( meang )
894 : call test%disp%show("iweight")
895 : call test%disp%show( iweight )
896 : call test%disp%show("rweight")
897 : call test%disp%show( rweight )
898 : call test%disp%show("iweisum")
899 : call test%disp%show( iweisum )
900 : call test%disp%show("iweisum_ref")
901 : call test%disp%show( iweisum_ref )
902 : call test%disp%show("rweisum")
903 : call test%disp%show( rweisum )
904 : call test%disp%show("rweisum_ref")
905 : call test%disp%show( rweisum_ref )
906 : call test%disp%show("var_ref")
907 : call test%disp%show( var_ref )
908 : call test%disp%show("var")
909 : call test%disp%show( var )
910 : call test%disp%show("vdiff")
911 : call test%disp%show( vdiff )
912 : call test%disp%show("mean_ref")
913 : call test%disp%show( mean_ref )
914 : call test%disp%show("mean")
915 : call test%disp%show( mean )
916 : call test%disp%show("mdiff")
917 : call test%disp%show( mdiff )
918 : call test%disp%skip()
919 : ! LCOV_EXCL_STOP
920 : end if
921 12800 : call test%assert(assertion, msg, int(line, IK))
922 12800 : end subroutine
923 :
924 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
925 : #elif getVarMerged_ENABLED || setVarMerged_ENABLED
926 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
927 :
928 : real(TKC) :: fracA
929 : integer(IK) :: itry, ndim, dim
930 : integer(IK) :: nsam, nsamA, nsamB
931 : TYPE_OF_SAMPLE, allocatable :: sample(:,:), meanDiff(:)
932 : real(TKC), allocatable :: varA(:), varB(:), var(:), var_ref(:), vdiff(:)
933 16 : assertion = .true._LK
934 16 : dim = 2_IK
935 :
936 816 : do itry = 1, 50
937 :
938 800 : nsamA = getUnifRand(2_IK, 5_IK)
939 800 : nsamB = getUnifRand(2_IK, 5_IK)
940 800 : nsam = nsamA + nsamB
941 800 : fracA = real(nsamA, TKC) / real(nsam, TKC)
942 :
943 : ! test D2 interface.
944 :
945 : block
946 :
947 800 : ndim = getUnifRand(1_IK, 5_IK)
948 800 : call setResized(var, ndim)
949 24708 : sample = getUnifRand(ONE, TWO, ndim, nsam)
950 :
951 4065 : var_ref = getVar(sample, dim)
952 4065 : varA = getVar(sample(:,1:nsamA), dim)
953 4065 : varB = getVar(sample(:,nsamA+1:), dim)
954 7330 : meanDiff = getMean(sample(:,1:nsamA), dim) - getMean(sample(:,nsamA+1:), dim)
955 :
956 : #if getVarMerged_ENABLED
957 2014 : var = getVarMerged(varB, varA, meanDiff, fracA)
958 2014 : vdiff = abs(var - var_ref)
959 1614 : assertion = assertion .and. all(vdiff < rtol)
960 400 : call report(__LINE__, SK_"The new `varMerged` must be computed correctly.")
961 : #elif setVarMerged_ENABLED
962 : ! new
963 400 : call setVarMerged(var, varB, varA, meanDiff, fracA)
964 2051 : vdiff = abs(var - var_ref)
965 1651 : assertion = assertion .and. all(vdiff < rtol)
966 400 : call report(__LINE__, SK_"The new `varMerged` must be computed correctly.")
967 : ! old
968 2051 : var = varB
969 400 : call setVarMerged(var, varA, meanDiff, fracA)
970 2051 : vdiff = abs(var - var_ref)
971 1651 : assertion = assertion .and. all(vdiff < rtol)
972 400 : call report(__LINE__, SK_"The in-place `varMerged` must be computed correctly.")
973 : #endif
974 :
975 : end block
976 :
977 : ! test D1 interface.
978 :
979 16 : block
980 :
981 800 : ndim = 1_IK
982 800 : call setResized(var, ndim)
983 12946 : sample = getUnifRand(ONE, TWO, ndim, nsam)
984 :
985 2400 : var_ref = getVar(sample, dim)
986 2400 : varA = getVar(sample(:,1:nsamA), dim)
987 2400 : varB = getVar(sample(:,nsamA+1:), dim)
988 4000 : meanDiff = getMean(sample(:,1:nsamA), dim) - getMean(sample(:,nsamA+1:), dim)
989 :
990 : #if getVarMerged_ENABLED
991 400 : var(1) = getVarMerged(varB(1), varA(1), meanDiff(1), fracA)
992 1200 : vdiff = abs(var - var_ref)
993 800 : assertion = assertion .and. all(vdiff < rtol)
994 400 : call report(__LINE__, SK_"The new `varMerged` must be computed correctly.")
995 : #elif setVarMerged_ENABLED
996 : ! new
997 400 : call setVarMerged(var(1), varB(1), varA(1), meanDiff(1), fracA)
998 1200 : vdiff = abs(var - var_ref)
999 800 : assertion = assertion .and. all(vdiff < rtol)
1000 400 : call report(__LINE__, SK_"The new `varMerged` must be computed correctly.")
1001 : ! old
1002 1200 : var = varB
1003 400 : call setVarMerged(var(1), varA(1), meanDiff(1), fracA)
1004 1200 : vdiff = abs(var - var_ref)
1005 800 : assertion = assertion .and. all(vdiff < rtol)
1006 400 : call report(__LINE__, SK_"The in-place `varMerged` must be computed correctly.")
1007 : #endif
1008 : end block
1009 :
1010 : end do
1011 :
1012 : contains
1013 :
1014 2400 : subroutine report(line, msg)
1015 : integer, intent(in) :: line
1016 : character(*, SK), intent(in) :: msg
1017 2400 : if (test%traceable .and. .not. assertion) then
1018 : ! LCOV_EXCL_START
1019 : call test%disp%skip()
1020 : call test%disp%show("[ndim, nsamA, nsamB, nsam, dim, shape(sample, IK)]")
1021 : call test%disp%show( [ndim, nsamA, nsamB, nsam, dim, shape(sample, IK)] )
1022 : call test%disp%show("sample")
1023 : call test%disp%show( sample )
1024 : call test%disp%show("meanDiff")
1025 : call test%disp%show( meanDiff )
1026 : call test%disp%show("varA")
1027 : call test%disp%show( varA )
1028 : call test%disp%show("varB")
1029 : call test%disp%show( varB )
1030 : call test%disp%show("var_ref")
1031 : call test%disp%show( var_ref )
1032 : call test%disp%show("var")
1033 : call test%disp%show( var )
1034 : call test%disp%show("vdiff")
1035 : call test%disp%show( vdiff )
1036 : call test%disp%skip()
1037 : ! LCOV_EXCL_STOP
1038 : end if
1039 2400 : call test%assert(assertion, msg, int(line, IK))
1040 2400 : end subroutine
1041 :
1042 : !%%%%%%%%%%%%%%%%%%%%%%%
1043 : #elif setVarMeanMerged_ENABLED
1044 : !%%%%%%%%%%%%%%%%%%%%%%%
1045 :
1046 : real(TKC) :: fracA
1047 : integer(IK) :: itry, ndim, dim
1048 : integer(IK) :: nsam, nsamA, nsamB
1049 : TYPE_OF_SAMPLE, allocatable :: sample(:,:), mean(:), mean_ref(:), meanA(:), meanB(:), mdiff(:)
1050 : real(TKC), allocatable :: varA(:), varB(:), var(:), var_ref(:), vdiff(:)
1051 8 : assertion = .true._LK
1052 8 : dim = 2_IK
1053 :
1054 408 : do itry = 1, 50
1055 :
1056 400 : nsamA = getUnifRand(2_IK, 5_IK)
1057 400 : nsamB = getUnifRand(2_IK, 5_IK)
1058 400 : nsam = nsamA + nsamB
1059 400 : fracA = real(nsamA, TKC) / real(nsam, TKC)
1060 :
1061 : ! test D2 interface.
1062 :
1063 : block
1064 :
1065 400 : ndim = getUnifRand(1_IK, 5_IK)
1066 400 : call setResized(var, ndim)
1067 400 : call setResized(mean, ndim)
1068 12191 : sample = getUnifRand(ONE, TWO, ndim, nsam)
1069 :
1070 2015 : var_ref = getVar(sample, dim)
1071 2015 : mean_ref = getMean(sample, dim)
1072 2015 : varA = getVar(sample(:,1:nsamA), dim)
1073 2015 : varB = getVar(sample(:,nsamA+1:), dim)
1074 2015 : meanA = getMean(sample(:,1:nsamA), dim)
1075 2015 : meanB = getMean(sample(:,nsamA+1:), dim)
1076 :
1077 : ! new
1078 400 : call setVarMeanMerged(var, mean, varB, meanB, varA, meanA, fracA)
1079 400 : call setAssertedAvg(__LINE__, SK_"new")
1080 400 : call setAssertedVar(__LINE__, SK_"new")
1081 :
1082 : ! old
1083 2015 : var = varB
1084 2015 : mean = meanB
1085 400 : call setVarMeanMerged(var, mean, varA, meanA, fracA)
1086 400 : call setAssertedAvg(__LINE__, SK_"in-place")
1087 400 : call setAssertedVar(__LINE__, SK_"in-place")
1088 :
1089 : end block
1090 :
1091 : ! test D1 interface.
1092 :
1093 8 : block
1094 :
1095 400 : ndim = 1_IK
1096 400 : call setResized(var, ndim)
1097 400 : call setResized(mean, ndim)
1098 6440 : sample = getUnifRand(ONE, TWO, ndim, nsam)
1099 :
1100 1200 : var_ref = getVar(sample, dim)
1101 1200 : mean_ref = getMean(sample, dim)
1102 1200 : varA = getVar(sample(:,1:nsamA), dim)
1103 1200 : varB = getVar(sample(:,nsamA+1:), dim)
1104 1200 : meanA = getMean(sample(:,1:nsamA), dim)
1105 1200 : meanB = getMean(sample(:,nsamA+1:), dim)
1106 :
1107 : ! new
1108 400 : call setVarMeanMerged(var(1), mean(1), varB(1), meanB(1), varA(1), meanA(1), fracA)
1109 400 : call setAssertedAvg(__LINE__, SK_"new")
1110 400 : call setAssertedVar(__LINE__, SK_"new")
1111 :
1112 : ! old
1113 1200 : var = varB
1114 1200 : mean = meanB
1115 400 : call setVarMeanMerged(var(1), mean(1), varA(1), meanA(1), fracA)
1116 400 : call setAssertedAvg(__LINE__, SK_"in-place")
1117 400 : call setAssertedVar(__LINE__, SK_"in-place")
1118 :
1119 : end block
1120 :
1121 : end do
1122 :
1123 : contains
1124 :
1125 1600 : subroutine setAssertedAvg(line, specific)
1126 : character(*, SK), intent(in) :: specific
1127 : integer, intent(in) :: line
1128 8018 : mdiff = abs(mean - mean_ref)
1129 4830 : assertion = assertion .and. all(mdiff < ctol)
1130 1600 : call report(line, SK_"The "//specific//SK_" `meanMerged` must be computed correctly.")
1131 1600 : end subroutine
1132 :
1133 1600 : subroutine setAssertedVar(line, specific)
1134 : character(*, SK), intent(in) :: specific
1135 : integer, intent(in) :: line
1136 6430 : vdiff = abs(var - var_ref)
1137 4830 : assertion = assertion .and. all(vdiff < rtol)
1138 1600 : call report(line, SK_"The "//specific//SK_" `varMerged` must be computed correctly.")
1139 1600 : end subroutine
1140 :
1141 3200 : subroutine report(line, msg)
1142 : integer, intent(in) :: line
1143 : character(*, SK), intent(in) :: msg
1144 3200 : if (test%traceable .and. .not. assertion) then
1145 : ! LCOV_EXCL_START
1146 : call test%disp%skip()
1147 : call test%disp%show("[ndim, nsamA, nsamB, nsam, dim, shape(sample, IK)]")
1148 : call test%disp%show( [ndim, nsamA, nsamB, nsam, dim, shape(sample, IK)] )
1149 : call test%disp%show("sample")
1150 : call test%disp%show( sample )
1151 : call test%disp%show("meanA")
1152 : call test%disp%show( meanA )
1153 : call test%disp%show("meanB")
1154 : call test%disp%show( meanB )
1155 : call test%disp%show("mean_ref")
1156 : call test%disp%show( mean_ref )
1157 : call test%disp%show("mean")
1158 : call test%disp%show( mean )
1159 : call test%disp%show("mdiff")
1160 : call test%disp%show( mdiff )
1161 : call test%disp%show("varA")
1162 : call test%disp%show( varA )
1163 : call test%disp%show("varB")
1164 : call test%disp%show( varB )
1165 : call test%disp%show("var_ref")
1166 : call test%disp%show( var_ref )
1167 : call test%disp%show("var")
1168 : call test%disp%show( var )
1169 : call test%disp%show("vdiff")
1170 : call test%disp%show( vdiff )
1171 : call test%disp%skip()
1172 : ! LCOV_EXCL_STOP
1173 : end if
1174 3200 : call test%assert(assertion, msg, int(line, IK))
1175 3200 : end subroutine
1176 :
1177 : #else
1178 : !%%%%%%%%%%%%%%%%%%%%%%%%
1179 : #error "Unrecognized interface."
1180 : !%%%%%%%%%%%%%%%%%%%%%%%%
1181 : #endif
1182 : #undef TYPE_OF_SAMPLE
|
