Line data Source code
1 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3 : !!!! !!!!
4 : !!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5 : !!!! !!!!
6 : !!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7 : !!!! !!!!
8 : !!!! This file is part of the ParaMonte library. !!!!
9 : !!!! !!!!
10 : !!!! LICENSE !!!!
11 : !!!! !!!!
12 : !!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13 : !!!! !!!!
14 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16 :
17 : !> \brief
18 : !> This include file contains the implementations of the tests of procedures of [pm_matrixChol](@ref pm_matrixChol).
19 : !>
20 : !> \fintest
21 : !>
22 : !> \author
23 : !> \AmirShahmoradi, Sunday 4:33 PM, September 19, 2021, Dallas, TX
24 :
25 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
26 :
27 : real(TKC), parameter :: RTOL = sqrt(epsilon(0._TKC))
28 : #if CK_ENABLED
29 : #define TYPE_OF cmplx
30 : #define GET_CONJG(X) conjg(X)
31 : #define TYPE_OF_MAT complex(TKC)
32 : complex(TKC), parameter :: ZERO = (0._TKC, 0._TKC), ONE = (1._TKC, 0._TKC), LB = (1._TKC, -1._TKC), UB = (2._TKC, 1._TKC), TOL = (RTOL, RTOL)
33 : #elif RK_ENABLED
34 : #define TYPE_OF real
35 : #define GET_CONJG(X) X
36 : #define TYPE_OF_MAT real(TKC)
37 : real(TKC), parameter :: ZERO = 0._TKC, ONE = 1._TKC, LB = 1._TKC, UB = 2._TKC, TOL = RTOL
38 : #else
39 : #error "Unrecognized interface."
40 : #endif
41 :
42 : !%%%%%%%%%%%%%%%%
43 : #if setChoLow_ENABLED
44 : !%%%%%%%%%%%%%%%%
45 :
46 : integer(IK), parameter :: ntry = 100_IK
47 : integer(IK) :: ndim, info, info_ref, itry
48 : TYPE_OF_MAT, allocatable :: chol_ref(:,:), chol(:,:), mat(:,:), diff(:,:), vdia(:), vdia_ref(:)
49 :
50 4 : assertion = .true._LK
51 :
52 404 : do itry = 1, ntry
53 :
54 : ! Set the matrix rank.
55 400 : ndim = getUnifRand(1_IK, 7_IK)
56 10842 : chol_ref = getFilled(ZERO, ndim, ndim)
57 : ! Generate random upper-triangular matrix.
58 10842 : mat = getCovRand(chol_ref(1,1), ndim)
59 400 : call setMatInit(mat, low, ZERO)
60 2449 : vdia = getFilled(ZERO, ndim)
61 400 : info = 0_IK
62 :
63 400 : call setMatChol(mat, uppDia, info_ref, chol_ref, transHerm)
64 : if (info_ref /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
65 4098 : vdia_ref = getMatCopy(lfpack, chol_ref, rdpack, dia)
66 400 : call setMatCopy(chol_ref, rdpack, mat, rdpack, uppDia)
67 11540 : chol = mat
68 400 : call setChoLow(chol, vdia, ndim)
69 10842 : diff = chol - chol_ref
70 10442 : assertion = assertion .and. all(diff < TOL)
71 2049 : assertion = assertion .and. all(abs(vdia - vdia_ref) < TOL)
72 400 : call report(__LINE__)
73 :
74 400 : info_ref = getUnifRand(1_IK, ndim)
75 400 : mat(info_ref, info_ref) = -mat(info_ref, info_ref)
76 400 : call setMatChol(mat, uppDia, info_ref, chol_ref, transHerm)
77 4098 : vdia_ref = getMatCopy(lfpack, chol_ref, rdpack, dia)
78 400 : call setMatCopy(chol_ref, rdpack, mat, rdpack, uppDia)
79 10842 : chol = mat
80 400 : call setChoLow(chol, vdia, ndim)
81 400 : info = int(-vdia(1), IK)
82 400 : assertion = assertion .and. info == info_ref
83 404 : call report(__LINE__)
84 :
85 : end do
86 :
87 : contains
88 :
89 800 : subroutine report(line)
90 : integer, intent(in) :: line
91 800 : if (test%traceable .and. .not. assertion) then
92 : ! LCOV_EXCL_START
93 : call test%disp%skip
94 : call test%disp%show("ndim")
95 : call test%disp%show( ndim )
96 : call test%disp%show("mat")
97 : call test%disp%show( mat )
98 : call test%disp%show("chol_ref")
99 : call test%disp%show( chol_ref )
100 : call test%disp%show("chol")
101 : call test%disp%show( chol )
102 : call test%disp%show("diff")
103 : call test%disp%show( diff )
104 : call test%disp%show("vdia_ref")
105 : call test%disp%show( vdia_ref )
106 : call test%disp%show("vdia")
107 : call test%disp%show( vdia )
108 : call test%disp%show("vdia - vdia_ref")
109 : call test%disp%show( vdia - vdia_ref )
110 : call test%disp%show("info_ref")
111 : call test%disp%show( info_ref )
112 : call test%disp%show("info")
113 : call test%disp%show( info )
114 : ! LCOV_EXCL_STOP
115 : end if
116 800 : call test%assert(assertion, SK_"The Cholesky factorization must not fail.", int(line, IK))
117 800 : end subroutine
118 :
119 : !%%%%%%%%%%%%%%%%%
120 : #elif getMatChol_ENABLED
121 : !%%%%%%%%%%%%%%%%%
122 :
123 : integer(IK) :: ndim, info, itry
124 : integer(IK), parameter :: ntry = 100_IK
125 : TYPE_OF_MAT, allocatable :: chol_ref(:,:), chol(:,:), mat(:,:), diff(:,:)
126 :
127 8 : assertion = .true._LK
128 :
129 808 : do itry = 1, ntry
130 :
131 : ! Set the matrix rank.
132 800 : ndim = getUnifRand(1_IK, 7_IK)
133 20480 : chol_ref = getFilled(ZERO, ndim, ndim)
134 : ! Generate random lower-triangular matrix.
135 20480 : mat = getCovRand(chol_ref(1,1), ndim)
136 800 : call setMatInit(mat, upp, ZERO)
137 :
138 : ! lowDia
139 :
140 : block
141 :
142 19680 : chol_ref = ZERO
143 800 : call setMatChol(mat, lowDia, info, chol_ref, nothing)
144 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
145 39360 : chol = getMatChol(mat, lowDia)
146 800 : call report(__LINE__, "lowDia")
147 :
148 19680 : chol_ref = ZERO
149 800 : call setMatChol(mat, lowDia, info, chol_ref, nothing)
150 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
151 39360 : chol = getMatChol(mat, lowDia, nothing)
152 800 : call report(__LINE__, "lowDia", "nothing")
153 :
154 19680 : chol_ref = ZERO
155 800 : call setMatChol(mat, lowDia, info, chol_ref, transHerm)
156 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
157 39360 : chol = getMatChol(mat, lowDia, transHerm)
158 800 : call report(__LINE__, "lowDia", "transHerm")
159 :
160 : end block
161 :
162 39360 : mat = transpose(GET_CONJG(mat))
163 :
164 : ! uppDia
165 :
166 8 : block
167 :
168 19680 : chol_ref = ZERO
169 800 : call setMatChol(mat, uppDia, info, chol_ref, nothing)
170 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
171 39360 : chol = getMatChol(mat, uppDia)
172 800 : call report(__LINE__, "uppDia")
173 :
174 19680 : chol_ref = ZERO
175 800 : call setMatChol(mat, uppDia, info, chol_ref, nothing)
176 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
177 39360 : chol = getMatChol(mat, uppDia, nothing)
178 800 : call report(__LINE__, "uppDia", "nothing")
179 :
180 19680 : chol_ref = ZERO
181 800 : call setMatChol(mat, uppDia, info, chol_ref, transHerm)
182 : if (info /= 0_IK) error stop getFine(__FILE__, __LINE__)//SK_": setMatChol() failed." ! LCOV_EXCL_LINE
183 39360 : chol = getMatChol(mat, uppDia, transHerm)
184 800 : call report(__LINE__, "uppDia", "transHerm")
185 :
186 : end block
187 :
188 : end do
189 :
190 : contains
191 :
192 4800 : subroutine report(line, subset, operation)
193 : integer, intent(in) :: line
194 : character(*, SK), intent(in) :: subset
195 : character(*, SK), intent(in), optional :: operation
196 122880 : diff = chol - chol_ref
197 118080 : assertion = assertion .and. all(diff < TOL)
198 4800 : if (test%traceable .and. .not. assertion) then
199 : ! LCOV_EXCL_START
200 : call test%disp%skip
201 : call test%disp%show("ndim")
202 : call test%disp%show( ndim )
203 : call test%disp%show("mat")
204 : call test%disp%show( mat )
205 : call test%disp%show("subset")
206 : call test%disp%show( subset )
207 : call test%disp%show("present(operation)")
208 : call test%disp%show( present(operation) )
209 : if (present(operation)) then
210 : call test%disp%show("operation")
211 : call test%disp%show( operation )
212 : end if
213 : call test%disp%show("chol_ref")
214 : call test%disp%show( chol_ref )
215 : call test%disp%show("chol")
216 : call test%disp%show( chol )
217 : call test%disp%show("diff")
218 : call test%disp%show( diff )
219 : ! LCOV_EXCL_STOP
220 : end if
221 4800 : call test%assert(assertion, SK_"The Cholesky factorization must not fail.", int(line, IK))
222 4800 : end subroutine
223 :
224 : !%%%%%%%%%%%%%%%%%
225 : #elif setMatChol_ENABLED
226 : !%%%%%%%%%%%%%%%%%
227 :
228 : integer(IK) :: ndim, info, info_def, itry
229 : integer(IK), parameter :: ntry = 100_IK
230 : TYPE_OF_MAT, allocatable :: choUpp_ref(:,:), choLow_ref(:,:), mat_ref(:,:)
231 : TYPE_OF_MAT, allocatable :: matUpp(:,:), matLow(:,:), choUpp(:,:), choLow(:,:), diff(:,:)
232 :
233 8 : info_def = 0_IK
234 8 : assertion = .true._LK
235 :
236 808 : do itry = 1, ntry
237 :
238 : ! Set the matrix rank.
239 800 : ndim = getUnifRand(1_IK, 7_IK)
240 : ! Generate random triangular matrix.
241 20966 : choLow_ref = getUnifRand(LB, UB, ndim, ndim)
242 40332 : choLow_ref = getMatCopy(rdpack, choLow_ref, rdpack, lowDia, init = ZERO)
243 : ! set the diagonals to positive real values.
244 4030 : call setMatInit(choLow_ref, dia, TYPE_OF(getUnifRand(1._TKC, 2._TKC, ndim), kind = TKC))
245 : ! Create the reference upper Cholesky matrix.
246 20966 : choUpp_ref = transpose(GET_CONJG(choLow_ref))
247 : ! Create the positive definite matrix.
248 582488 : mat_ref = matmul(choLow_ref, choUpp_ref)
249 40332 : matUpp = getMatCopy(rdpack, mat_ref, rdpack, uppDia, init = ZERO)
250 40332 : matLow = getMatCopy(rdpack, mat_ref, rdpack, lowDia, init = ZERO)
251 :
252 : ! Test the routines.
253 :
254 : internal_interface: block
255 :
256 : ! internal paramonte interface, no overwrite.
257 :
258 20966 : choUpp = getFilled(ZERO, ndim, ndim)
259 20966 : choLow = getFilled(ZERO, ndim, ndim)
260 :
261 20966 : choUpp = getFilled(ZERO, ndim, ndim)
262 800 : call setMatChol(matUpp, uppDia, info, choUpp, nothing)
263 800 : assertion = assertion .and. info == 0_IK
264 800 : call report(__LINE__, "uppDia", "nothing", .false._LK, info)
265 20966 : diff = abs(choUpp - choUpp_ref)
266 20166 : assertion = assertion .and. all(diff < TOL)
267 800 : call report(__LINE__, "uppDia", "nothing", .false._LK, diff = diff)
268 :
269 20966 : choLow = getFilled(ZERO, ndim, ndim)
270 800 : call setMatChol(matLow, lowDia, info, choLow, nothing)
271 800 : assertion = assertion .and. info == 0_IK
272 800 : call report(__LINE__, "lowDia", "nothing", .false._LK, info)
273 20966 : diff = abs(choLow - choLow_ref)
274 20166 : assertion = assertion .and. all(diff < TOL)
275 800 : call report(__LINE__, "lowDia", "nothing", .false._LK, diff = diff)
276 :
277 20966 : choLow = getFilled(ZERO, ndim, ndim)
278 800 : call setMatChol(matUpp, uppDia, info, choLow, transHerm)
279 800 : assertion = assertion .and. info == 0_IK
280 800 : call report(__LINE__, "uppDia", "tranHerm", .false._LK, info)
281 20966 : diff = abs(choLow - choLow_ref)
282 20166 : assertion = assertion .and. all(diff < TOL)
283 800 : call report(__LINE__, "uppDia", "tranHerm", .false._LK, diff = diff)
284 :
285 20966 : choUpp = getFilled(ZERO, ndim, ndim)
286 800 : call setMatChol(matLow, lowDia, info, choUpp, transHerm)
287 800 : assertion = assertion .and. info == 0_IK
288 800 : call report(__LINE__, "lowDia", "transHerm", .false._LK, info)
289 20966 : diff = abs(choUpp - choUpp_ref)
290 20166 : assertion = assertion .and. all(diff < TOL)
291 800 : call report(__LINE__, "lowDia", "transHerm", .false._LK, diff = diff)
292 :
293 : ! internal paramonte interface, with overwrite.
294 :
295 24996 : choUpp = getFilled(ZERO, ndim, ndim + 1)
296 24996 : choLow = getFilled(ZERO, ndim, ndim + 1)
297 :
298 20166 : choUpp(:,1:ndim) = matUpp
299 800 : call setMatChol(choUpp(:,1:ndim), uppDia, info, choUpp(:,1:ndim), nothing)
300 800 : assertion = assertion .and. info == 0_IK
301 800 : call report(__LINE__, "uppDia", "nothing", .true._LK, info)
302 800 : diff = getDiff(choUpp(:,1:ndim), choUpp_ref, uppDia)
303 20166 : assertion = assertion .and. all(diff < TOL)
304 800 : call report(__LINE__, "uppDia", "nothing", .true._LK, diff = diff)
305 :
306 20166 : choLow(:,1:ndim) = matLow
307 800 : call setMatChol(choLow(:,1:ndim), lowDia, info, choLow(:,1:ndim), nothing)
308 800 : assertion = assertion .and. info == 0_IK
309 800 : call report(__LINE__, "lowDia", "nothing", .true._LK, info)
310 800 : diff = getDiff(choLow(:,1:ndim), choLow_ref, lowDia)
311 20166 : assertion = assertion .and. all(diff < TOL)
312 800 : call report(__LINE__, "lowDia", "nothing", .true._LK, diff = diff)
313 :
314 20166 : choUpp(:,1:ndim) = matLow
315 800 : call setMatChol(choUpp(:,1:ndim), lowDia, info, choUpp(:,2:ndim+1), transHerm)
316 800 : assertion = assertion .and. info == 0_IK
317 800 : call report(__LINE__, "uppDia", "transHerm", .true._LK, info)
318 800 : diff = getDiff(choUpp(:,2:ndim+1), choUpp_ref, uppDia)
319 20166 : assertion = assertion .and. all(diff < TOL)
320 800 : call report(__LINE__, "uppDia", "transHerm", .true._LK, diff = diff)
321 :
322 20166 : choLow(:,2:ndim+1) = matUpp
323 800 : call setMatChol(choLow(:,2:ndim+1), uppDia, info, choLow(:,1:ndim), transHerm)
324 800 : assertion = assertion .and. info == 0_IK
325 800 : call report(__LINE__, "lowDia", "transHerm", .true._LK, info)
326 800 : diff = getDiff(choLow(:,1:ndim), choLow_ref, lowDia)
327 20166 : assertion = assertion .and. all(diff < TOL)
328 800 : call report(__LINE__, "lowDia", "transHerm", .true._LK, diff = diff)
329 :
330 : end block internal_interface
331 :
332 : recursion_interface: block
333 :
334 : ! implicit interface.
335 :
336 20966 : choUpp = getFilled(ZERO, ndim, ndim)
337 20966 : choLow = getFilled(ZERO, ndim, ndim)
338 :
339 40332 : choUpp = getMatCopy(rdpack, matUpp, rdpack, uppDia, init = ZERO)
340 800 : call setMatChol(choUpp, uppDia, info, recursion)
341 800 : assertion = assertion .and. info == 0_IK
342 800 : call report(__LINE__, "uppDia", "recursion", .true._LK, info)
343 20966 : diff = abs(choUpp - choUpp_ref)
344 20166 : assertion = assertion .and. all(diff < TOL)
345 800 : call report(__LINE__, "uppDia", "recursion", .true._LK, diff = diff)
346 :
347 40332 : choLow = getMatCopy(rdpack, matLow, rdpack, lowDia, init = ZERO)
348 800 : call setMatChol(choLow, lowDia, info, recursion)
349 800 : assertion = assertion .and. info == 0_IK
350 800 : call report(__LINE__, "lowDia", "recursion", .true._LK, info)
351 20966 : diff = abs(choLow - choLow_ref)
352 20166 : assertion = assertion .and. all(diff < TOL)
353 800 : call report(__LINE__, "lowDia", "recursion", .true._LK, diff = diff)
354 :
355 : end block recursion_interface
356 :
357 : iteration_interface: block
358 :
359 : ! implicit interface.
360 :
361 20966 : choUpp = getFilled(ZERO, ndim, ndim)
362 20966 : choLow = getFilled(ZERO, ndim, ndim)
363 :
364 40332 : choUpp = getMatCopy(rdpack, matUpp, rdpack, uppDia, init = ZERO)
365 800 : call setMatChol(choUpp, uppDia, info, iteration)
366 800 : assertion = assertion .and. info == 0_IK
367 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, info)
368 20966 : diff = abs(choUpp - choUpp_ref)
369 20166 : assertion = assertion .and. all(diff < TOL)
370 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, diff = diff)
371 :
372 40332 : choLow = getMatCopy(rdpack, matLow, rdpack, lowDia, init = ZERO)
373 800 : call setMatChol(choLow, lowDia, info, iteration)
374 800 : assertion = assertion .and. info == 0_IK
375 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, info)
376 20966 : diff = abs(choLow - choLow_ref)
377 20166 : assertion = assertion .and. all(diff < TOL)
378 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, diff = diff)
379 :
380 : end block iteration_interface
381 :
382 : iteration_bdim_interface: block
383 :
384 : ! implicit interface.
385 :
386 : integer(IK) :: bdim
387 :
388 800 : bdim = getUnifRand(2_IK, 2_IK * ndim + 1_IK)
389 20966 : choUpp = getFilled(ZERO, ndim, ndim)
390 20966 : choLow = getFilled(ZERO, ndim, ndim)
391 :
392 40332 : choUpp = getMatCopy(rdpack, matUpp, rdpack, uppDia, init = ZERO)
393 800 : call setMatChol(choUpp, uppDia, info, iteration, bdim = bdim)
394 800 : assertion = assertion .and. info == 0_IK
395 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, info, bdim = bdim)
396 20966 : diff = abs(choUpp - choUpp_ref)
397 20166 : assertion = assertion .and. all(diff < TOL)
398 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, diff = diff, bdim = bdim)
399 :
400 40332 : choLow = getMatCopy(rdpack, matLow, rdpack, lowDia, init = ZERO)
401 800 : call setMatChol(choLow, lowDia, info, iteration, bdim = bdim)
402 800 : assertion = assertion .and. info == 0_IK
403 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, info, bdim = bdim)
404 20966 : diff = abs(choLow - choLow_ref)
405 20166 : assertion = assertion .and. all(diff < TOL)
406 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, diff = diff, bdim = bdim)
407 :
408 : end block iteration_bdim_interface
409 :
410 8 : nonposdef_interface: block
411 :
412 : ! implicit interface.
413 :
414 : integer(IK) :: bdim
415 :
416 800 : info_def = getUnifRand(1_IK, ndim)
417 800 : bdim = getUnifRand(2_IK, 2_IK * ndim + 1_IK)
418 41932 : mat_ref = getMatInit([ndim, ndim], uppLowDia, ZERO, ZERO, getUnifRand(LB, UB, ndim))
419 800 : mat_ref(info_def, info_def) = -ONE
420 :
421 : ! default
422 :
423 20966 : choUpp = mat_ref
424 800 : call setMatChol(mat_ref, uppDia, info, choUpp, nothing)
425 800 : assertion = assertion .and. info == info_def
426 800 : call report(__LINE__, "uppDia", "nothing", .false._LK, info)
427 :
428 20966 : choLow = mat_ref
429 800 : call setMatChol(mat_ref, lowDia, info, choLow, nothing)
430 800 : assertion = assertion .and. info == info_def
431 800 : call report(__LINE__, "lowDia", "nothing", .false._LK, info)
432 :
433 20966 : choUpp = mat_ref
434 800 : call setMatChol(mat_ref, uppDia, info, choUpp, transHerm)
435 800 : assertion = assertion .and. info == info_def
436 800 : call report(__LINE__, "uppDia", "transHerm", .false._LK, info)
437 :
438 20966 : choLow = mat_ref
439 800 : call setMatChol(mat_ref, lowDia, info, choLow, transHerm)
440 800 : assertion = assertion .and. info == info_def
441 800 : call report(__LINE__, "lowDia", "transHerm", .false._LK, info)
442 :
443 20966 : choUpp = mat_ref
444 800 : call setMatChol(choUpp, uppDia, info)
445 800 : assertion = assertion .and. info == info_def
446 800 : call report(__LINE__, "uppDia", "nothing", .true._LK, info)
447 :
448 20966 : choLow = mat_ref
449 800 : call setMatChol(choLow, lowDia, info)
450 800 : assertion = assertion .and. info == info_def
451 800 : call report(__LINE__, "lowDia", "nothing", .true._LK, info)
452 :
453 : ! recursion
454 :
455 20966 : choUpp = mat_ref
456 800 : call setMatChol(choUpp, uppDia, info, recursion)
457 800 : assertion = assertion .and. info == info_def
458 800 : call report(__LINE__, "uppDia", "recursion", .true._LK, info)
459 :
460 20966 : choLow = mat_ref
461 800 : call setMatChol(choLow, lowDia, info, recursion)
462 800 : assertion = assertion .and. info == info_def
463 800 : call report(__LINE__, "lowDia", "recursion", .true._LK, info)
464 :
465 : ! iteration
466 :
467 20966 : choUpp = mat_ref
468 800 : call setMatChol(choUpp, uppDia, info, iteration)
469 800 : assertion = assertion .and. info == info_def
470 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, info)
471 :
472 20966 : choLow = mat_ref
473 800 : call setMatChol(choLow, lowDia, info, iteration)
474 800 : assertion = assertion .and. info == info_def
475 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, info)
476 :
477 20966 : choUpp = mat_ref
478 800 : call setMatChol(choUpp, uppDia, info, iteration, bdim = bdim)
479 800 : assertion = assertion .and. info == info_def
480 800 : call report(__LINE__, "uppDia", "iteration", .true._LK, info, bdim = bdim)
481 :
482 20966 : choLow = mat_ref
483 800 : call setMatChol(choLow, lowDia, info, iteration, bdim = bdim)
484 800 : assertion = assertion .and. info == info_def
485 800 : call report(__LINE__, "lowDia", "iteration", .true._LK, info, bdim = bdim)
486 :
487 : end block nonposdef_interface
488 :
489 : end do
490 :
491 : contains
492 :
493 3200 : pure function getDiff(mat, ref, subset) result(diff)
494 : TYPE_OF_MAT, intent(in) :: mat(:,:), ref(:,:)
495 : TYPE_OF_MAT :: diff(size(ref, 1, IK), size(ref, 2, IK))
496 : class(subset_type), intent(in) :: subset
497 : integer(IK) :: idim
498 80664 : diff = 0._TKC
499 16120 : do idim = 1, size(ref, 1, IK)
500 16120 : if (same_type_as(subset, uppDia)) then
501 25826 : diff(idim, idim :) = abs(mat(idim, idim :) - ref(idim, idim :))
502 6460 : elseif (same_type_as(subset, lowDia)) then
503 25826 : diff(idim :, idim) = abs(mat(idim :, idim) - ref(idim :, idim))
504 : else
505 0 : error stop "Internal library error: Unrecognized `subset` type. Please report this error to the developers."
506 : end if
507 : end do
508 3200 : end function
509 :
510 32000 : subroutine report(line, subset, operation, overwrite, info, diff, bdim)
511 : integer, intent(in) :: line
512 : integer, intent(in), optional :: info
513 : integer, intent(in), optional :: bdim
514 : character(*, SK), intent(in) :: subset, operation
515 : TYPE_OF_MAT, intent(in), optional :: diff(:,:)
516 : logical(LK), intent(in) :: overwrite
517 32000 : if (test%traceable .and. .not. assertion) then
518 : ! LCOV_EXCL_START
519 : call test%disp%skip
520 : call test%disp%show("ndim")
521 : call test%disp%show( ndim )
522 : call test%disp%show("choUpp_ref")
523 : call test%disp%show( choUpp_ref )
524 : call test%disp%show("choUpp")
525 : call test%disp%show( choUpp )
526 : call test%disp%show("choUpp - choUpp_ref")
527 : call test%disp%show( choUpp - choUpp_ref )
528 : call test%disp%show("choLow_ref")
529 : call test%disp%show( choLow_ref )
530 : call test%disp%show("choLow")
531 : call test%disp%show( choLow )
532 : call test%disp%show("choLow - choLow_ref")
533 : call test%disp%show( choLow - choLow_ref )
534 : call test%disp%show("mat_ref")
535 : call test%disp%show( mat_ref )
536 : call test%disp%show("subset")
537 : call test%disp%show( subset )
538 : call test%disp%show("operation")
539 : call test%disp%show( operation )
540 : call test%disp%show("overwrite")
541 : call test%disp%show( overwrite )
542 : call test%disp%show("present(info)")
543 : call test%disp%show( present(info) )
544 : if (present(info)) then
545 : call test%disp%show("info")
546 : call test%disp%show( info )
547 : end if
548 : call test%disp%show("info_def")
549 : call test%disp%show( info_def )
550 : if (present(diff)) then
551 : call test%disp%show("diff")
552 : call test%disp%show( diff )
553 : end if
554 : call test%disp%show("present(bdim)")
555 : call test%disp%show( present(bdim) )
556 : if (present(bdim)) then
557 : call test%disp%show("bdim")
558 : call test%disp%show( bdim )
559 : end if
560 : ! LCOV_EXCL_STOP
561 : end if
562 32000 : call test%assert(assertion, SK_"The Cholesky factorization must not fail.", int(line, IK))
563 32000 : end subroutine
564 :
565 : #else
566 : !%%%%%%%%%%%%%%%%%%%%%%%%%
567 : #error "Unrecognized interrface."
568 : !%%%%%%%%%%%%%%%%%%%%%%%%%
569 : #endif
570 : #undef TYPE_OF_MAT
571 : #undef GET_CONJG
572 : #undef TYPE_OF
|
