Line data Source code
1 : ! This file is part of ELPA.
2 : !
3 : ! The ELPA library was originally created by the ELPA consortium,
4 : ! consisting of the following organizations:
5 : !
6 : ! - Max Planck Computing and Data Facility (MPCDF), formerly known as
7 : ! Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
8 : ! - Bergische Universität Wuppertal, Lehrstuhl für angewandte
9 : ! Informatik,
10 : ! - Technische Universität München, Lehrstuhl für Informatik mit
11 : ! Schwerpunkt Wissenschaftliches Rechnen ,
12 : ! - Fritz-Haber-Institut, Berlin, Abt. Theorie,
13 : ! - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
14 : ! Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
15 : ! and
16 : ! - IBM Deutschland GmbH
17 : !
18 : !
19 : ! More information can be found here:
20 : ! http://elpa.mpcdf.mpg.de/
21 : !
22 : ! ELPA is free software: you can redistribute it and/or modify
23 : ! it under the terms of the version 3 of the license of the
24 : ! GNU Lesser General Public License as published by the Free
25 : ! Software Foundation.
26 : !
27 : ! ELPA is distributed in the hope that it will be useful,
28 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : ! GNU Lesser General Public License for more details.
31 : !
32 : ! You should have received a copy of the GNU Lesser General Public License
33 : ! along with ELPA. If not, see <http://www.gnu.org/licenses/>
34 : !
35 : ! ELPA reflects a substantial effort on the part of the original
36 : ! ELPA consortium, and we ask you to respect the spirit of the
37 : ! license that we chose: i.e., please contribute any changes you
38 : ! may have back to the original ELPA library distribution, and keep
39 : ! any derivatives of ELPA under the same license that we chose for
40 : ! the original distribution, the GNU Lesser General Public License.
41 : !
42 : !
43 : #include "config-f90.h"
44 :
45 : ! Define one of TEST_REAL or TEST_COMPLEX
46 : ! Define one of TEST_SINGLE or TEST_DOUBLE
47 : ! Define one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE
48 : ! Define TEST_GPU \in [0, 1]
49 : ! Define either TEST_ALL_KERNELS or a TEST_KERNEL \in [any valid kernel]
50 :
51 : #if !(defined(TEST_REAL) ^ defined(TEST_COMPLEX))
52 : error: define exactly one of TEST_REAL or TEST_COMPLEX
53 : #endif
54 :
55 : #if !(defined(TEST_SINGLE) ^ defined(TEST_DOUBLE))
56 : error: define exactly one of TEST_SINGLE or TEST_DOUBLE
57 : #endif
58 :
59 : #if !(defined(TEST_SOLVER_1STAGE) ^ defined(TEST_SOLVER_2STAGE) ^ defined(TEST_SCALAPACK_ALL) ^ defined(TEST_SCALAPACK_PART))
60 : error: define exactly one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE or TEST_SCALAPACK_ALL or TEST_SCALAPACK_PART
61 : #endif
62 :
63 : #ifdef TEST_SOLVER_1STAGE
64 : #ifdef TEST_ALL_KERNELS
65 : error: TEST_ALL_KERNELS cannot be defined for TEST_SOLVER_1STAGE
66 : #endif
67 : #ifdef TEST_KERNEL
68 : error: TEST_KERNEL cannot be defined for TEST_SOLVER_1STAGE
69 : #endif
70 : #endif
71 :
72 : #ifdef TEST_SOLVER_2STAGE
73 : #if !(defined(TEST_KERNEL) ^ defined(TEST_ALL_KERNELS))
74 : error: define either TEST_ALL_KERNELS or a valid TEST_KERNEL
75 : #endif
76 : #endif
77 :
78 :
79 : #ifdef TEST_SINGLE
80 : # define EV_TYPE real(kind=C_FLOAT)
81 : # ifdef TEST_REAL
82 : # define MATRIX_TYPE real(kind=C_FLOAT)
83 : # else
84 : # define MATRIX_TYPE complex(kind=C_FLOAT_COMPLEX)
85 : # endif
86 : #else
87 : # define EV_TYPE real(kind=C_DOUBLE)
88 : # ifdef TEST_REAL
89 : # define MATRIX_TYPE real(kind=C_DOUBLE)
90 : # else
91 : # define MATRIX_TYPE complex(kind=C_DOUBLE_COMPLEX)
92 : # endif
93 : #endif
94 :
95 : #ifdef TEST_REAL
96 : #define KERNEL_KEY "real_kernel"
97 : #endif
98 : #ifdef TEST_COMPLEX
99 : #define KERNEL_KEY "complex_kernel"
100 : #endif
101 :
102 : #include "assert.h"
103 :
104 4944 : program test
105 4944 : use elpa
106 :
107 : use test_util
108 : use test_setup_mpi
109 : use test_prepare_matrix
110 : use test_read_input_parameters
111 : use test_blacs_infrastructure
112 : use test_check_correctness
113 : use test_analytic
114 : #ifdef WITH_SCALAPACK_TESTS
115 : use test_scalapack
116 : #endif
117 :
118 : #ifdef HAVE_REDIRECT
119 : use test_redirect
120 : #endif
121 : implicit none
122 :
123 : ! matrix dimensions
124 : integer :: na, nev, nblk
125 :
126 : ! mpi
127 : integer :: myid, nprocs
128 : integer :: na_cols, na_rows ! local matrix size
129 : integer :: np_cols, np_rows ! number of MPI processes per column/row
130 : integer :: my_prow, my_pcol ! local MPI task position (my_prow, my_pcol) in the grid (0..np_cols -1, 0..np_rows -1)
131 : integer :: mpierr
132 :
133 : ! blacs
134 : integer :: my_blacs_ctxt, sc_desc(9), info, nprow, npcol
135 :
136 : ! The Matrix
137 4944 : MATRIX_TYPE, allocatable :: a(:,:), as(:,:)
138 : #if defined(TEST_HERMITIAN_MULTIPLY)
139 384 : MATRIX_TYPE, allocatable :: b(:,:), c(:,:)
140 : #endif
141 : ! eigenvectors
142 4944 : MATRIX_TYPE, allocatable :: z(:,:)
143 : ! eigenvalues
144 4944 : EV_TYPE, allocatable :: ev(:), ev_analytic(:)
145 :
146 : logical :: check_all_evals
147 :
148 : #if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || TEST_QR_DECOMPOSITION == 1 || defined(TEST_HERMITIAN_MULTIPLY)
149 8544 : EV_TYPE, allocatable :: d(:), sd(:), ds(:), sds(:)
150 : EV_TYPE :: diagonalELement, subdiagonalElement
151 : #endif
152 : #if defined(TEST_CHOLESKY)
153 576 : MATRIX_TYPE, allocatable :: d(:), sd(:), ds(:), sds(:)
154 : MATRIX_TYPE :: diagonalELement, subdiagonalElement
155 : #endif
156 :
157 : integer :: error, status
158 :
159 : type(output_t) :: write_to_file
160 : class(elpa_t), pointer :: e
161 : #ifdef TEST_ALL_KERNELS
162 : integer :: i
163 : #endif
164 : #ifdef TEST_ALL_LAYOUTS
165 : character(len=1), parameter :: layouts(2) = [ 'C', 'R' ]
166 : integer :: i_layout
167 : #endif
168 : integer :: kernel
169 : character(len=1) :: layout
170 : logical :: do_test_numeric_residual, do_test_analytic_eigenvalues, &
171 : do_test_analytic_eigenvalues_eigenvectors, &
172 : do_test_frank_eigenvalues, &
173 : do_test_toeplitz_eigenvalues, do_test_cholesky, &
174 : do_test_hermitian_multiply
175 :
176 4944 : call read_input_parameters_traditional(na, nev, nblk, write_to_file)
177 4944 : call setup_mpi(myid, nprocs)
178 : #ifdef HAVE_REDIRECT
179 : #ifdef WITH_MPI
180 : call MPI_BARRIER(MPI_COMM_WORLD, mpierr)
181 : call redirect_stdout(myid)
182 : #endif
183 : #endif
184 :
185 4944 : check_all_evals = .true.
186 :
187 :
188 4944 : do_test_numeric_residual = .false.
189 4944 : do_test_analytic_eigenvalues = .false.
190 4944 : do_test_analytic_eigenvalues_eigenvectors = .false.
191 4944 : do_test_frank_eigenvalues = .false.
192 4944 : do_test_toeplitz_eigenvalues = .false.
193 :
194 4944 : do_test_cholesky = .false.
195 : #if defined(TEST_CHOLESKY)
196 288 : do_test_cholesky = .true.
197 : #endif
198 4944 : do_test_hermitian_multiply = .false.
199 : #if defined(TEST_HERMITIAN_MULTIPLY)
200 384 : do_test_hermitian_multiply = .true.
201 : #endif
202 :
203 4944 : status = 0
204 4944 : if (elpa_init(CURRENT_API_VERSION) /= ELPA_OK) then
205 0 : print *, "ELPA API version not supported"
206 0 : stop 1
207 : endif
208 :
209 4944 : if (myid == 0) then
210 3264 : print '((a,i0))', 'Program ' // TEST_CASE
211 3264 : print *, ""
212 : endif
213 :
214 : #ifdef TEST_ALL_LAYOUTS
215 : do i_layout = 1, size(layouts) ! layouts
216 : layout = layouts(i_layout)
217 : do np_cols = 1, nprocs ! factors
218 : if (mod(nprocs,np_cols) /= 0 ) then
219 : cycle
220 : endif
221 : #else
222 4944 : layout = 'C'
223 4944 : do np_cols = NINT(SQRT(REAL(nprocs))),2,-1
224 0 : if(mod(nprocs,np_cols) == 0 ) exit
225 : enddo
226 : #endif
227 :
228 4944 : np_rows = nprocs/np_cols
229 4944 : assert(nprocs == np_rows * np_cols)
230 :
231 4944 : if (myid == 0) then
232 3264 : print '((a,i0))', 'Matrix size: ', na
233 3264 : print '((a,i0))', 'Num eigenvectors: ', nev
234 3264 : print '((a,i0))', 'Blocksize: ', nblk
235 : #ifdef WITH_MPI
236 1680 : print '((a,i0))', 'Num MPI proc: ', nprocs
237 1680 : print '(3(a,i0))','Number of processor rows=',np_rows,', cols=',np_cols,', total=',nprocs
238 1680 : print '(a)', 'Process layout: ' // layout
239 : #endif
240 3264 : print *,''
241 : endif
242 :
243 : #if TEST_QR_DECOMPOSITION == 1
244 :
245 : #if TEST_GPU == 1
246 : #ifdef WITH_MPI
247 : call mpi_finalize(mpierr)
248 : #endif
249 : stop 77
250 : #endif /* TEST_GPU */
251 288 : if (nblk .lt. 64) then
252 288 : if (myid .eq. 0) then
253 192 : print *,"At the moment QR decomposition need blocksize of at least 64"
254 : endif
255 288 : if ((na .lt. 64) .and. (myid .eq. 0)) then
256 0 : print *,"This is why the matrix size must also be at least 64 or only 1 MPI task can be used"
257 : endif
258 :
259 : #ifdef WITH_MPI
260 192 : call mpi_finalize(mpierr)
261 : #endif
262 288 : stop 77
263 : endif
264 : #endif /* TEST_QR_DECOMPOSITION */
265 :
266 :
267 : call set_up_blacsgrid(mpi_comm_world, np_rows, np_cols, layout, &
268 4656 : my_blacs_ctxt, my_prow, my_pcol)
269 :
270 : call set_up_blacs_descriptor(na, nblk, my_prow, my_pcol, np_rows, np_cols, &
271 4656 : na_rows, na_cols, sc_desc, my_blacs_ctxt, info)
272 :
273 4656 : allocate(a (na_rows,na_cols))
274 4656 : allocate(as(na_rows,na_cols))
275 4656 : allocate(z (na_rows,na_cols))
276 4656 : allocate(ev(na))
277 :
278 : #ifdef TEST_HERMITIAN_MULTIPLY
279 384 : allocate(b (na_rows,na_cols))
280 384 : allocate(c (na_rows,na_cols))
281 : #endif
282 :
283 : #if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || TEST_QR_DECOMPOSITION == 1|| defined(TEST_CHOLESKY)
284 4272 : allocate(d (na), ds(na))
285 4272 : allocate(sd (na), sds(na))
286 4272 : allocate(ev_analytic(na))
287 : #endif
288 :
289 4656 : a(:,:) = 0.0
290 4656 : z(:,:) = 0.0
291 4656 : ev(:) = 0.0
292 :
293 : #if defined(TEST_MATRIX_RANDOM) && !defined(TEST_SOLVE_TRIDIAGONAL) && !defined(TEST_CHOLESKY) && !defined(TEST_EIGENVALUES)
294 : ! the random matrix can be used in allmost all tests; but for some no
295 : ! correctness checks have been implemented; do not allow these
296 : ! combinations
297 : ! RANDOM + TEST_SOLVE_TRIDIAGONAL: we need a TOEPLITZ MATRIX
298 : ! RANDOM + TEST_CHOLESKY: no correctness check yet implemented
299 : ! RANDOM + TEST_EIGENVALUES: no correctness check known
300 :
301 : ! We also have to take care of special case in TEST_EIGENVECTORS
302 : #if !defined(TEST_EIGENVECTORS)
303 288 : call prepare_matrix_random(na, myid, sc_desc, a, z, as)
304 :
305 288 : do_test_analytic_eigenvalues = .false.
306 288 : do_test_analytic_eigenvalues_eigenvectors = .false.
307 288 : do_test_frank_eigenvalues = .false.
308 288 : do_test_toeplitz_eigenvalues = .false.
309 :
310 : #else /* TEST_EIGENVECTORS */
311 :
312 864 : if (nev .ge. 1) then
313 864 : call prepare_matrix_random(na, myid, sc_desc, a, z, as)
314 :
315 864 : do_test_analytic_eigenvalues = .false.
316 864 : do_test_analytic_eigenvalues_eigenvectors = .false.
317 864 : do_test_frank_eigenvalues = .false.
318 864 : do_test_toeplitz_eigenvalues = .false.
319 : #ifndef TEST_HERMITIAN_MULTIPLY
320 864 : do_test_numeric_residual = .true.
321 : #endif
322 : else
323 0 : if (myid .eq. 0) then
324 0 : print *,"At the moment with the random matrix you need nev >=1"
325 : endif
326 : #ifdef WITH_MPI
327 0 : call mpi_finalize(mpierr)
328 : #endif
329 0 : stop 77
330 :
331 : endif
332 :
333 : #endif /* TEST_EIGENVECTORS */
334 : #endif /* (TEST_MATRIX_RANDOM) */
335 : #if defined(TEST_MATRIX_RANDOM) && (defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_CHOLESKY) || defined(TEST_EIGENVALUES))
336 : #error "Random matrix is not allowed in this configuration"
337 : #endif
338 :
339 : #if defined(TEST_MATRIX_ANALYTIC) && !defined(TEST_SOLVE_TRIDIAGONAL) && !defined(TEST_CHOLESKY)
340 : ! the analytic matrix can be used in allmost all tests; but for some no
341 : ! correctness checks have been implemented; do not allow these
342 : ! combinations
343 : ! ANALYTIC + TEST_SOLVE_TRIDIAGONAL: we need a TOEPLITZ MATRIX
344 : ! ANALTIC + TEST_CHOLESKY: no correctness check yet implemented
345 :
346 1056 : call prepare_matrix_analytic(na, a, nblk, myid, np_rows, np_cols, my_prow, my_pcol)
347 1056 : as(:,:) = a
348 :
349 1056 : do_test_numeric_residual = .false.
350 1056 : do_test_analytic_eigenvalues_eigenvectors = .false.
351 : #ifndef TEST_HERMITIAN_MULTIPLY
352 1056 : do_test_analytic_eigenvalues = .true.
353 : #endif
354 : #if defined(TEST_EIGENVECTORS)
355 1056 : if (nev .ge. 1) then
356 1056 : do_test_analytic_eigenvalues_eigenvectors = .true.
357 1056 : do_test_numeric_residual = .true.
358 : else
359 0 : do_test_analytic_eigenvalues_eigenvectors = .false.
360 0 : do_test_numeric_residual = .false.
361 : endif
362 : #endif
363 1056 : do_test_frank_eigenvalues = .false.
364 1056 : do_test_toeplitz_eigenvalues = .false.
365 : #endif /* TEST_MATRIX_ANALYTIC */
366 : #if defined(TEST_MATRIX_ANALYTIC) && (defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_CHOLESKY))
367 : #error "Analytic matrix is not allowd in this configuration"
368 : #endif
369 :
370 : #if defined(TEST_MATRIX_TOEPLITZ)
371 : ! The Toeplitz matrix works in each test
372 : #ifdef TEST_SINGLE
373 624 : diagonalElement = 0.45_c_float
374 624 : subdiagonalElement = 0.78_c_float
375 : #else
376 1248 : diagonalElement = 0.45_c_double
377 1248 : subdiagonalElement = 0.78_c_double
378 : #endif
379 :
380 : #if defined(TEST_CHOLESKY)
381 : #ifdef TEST_SINGLE
382 96 : diagonalElement = (2.546_c_float, 0.0_c_float)
383 96 : subdiagonalElement = (0.0_c_float, 0.0_c_float)
384 : #else
385 192 : diagonalElement = (2.546_c_double, 0.0_c_double)
386 192 : subdiagonalElement = (0.0_c_double, 0.0_c_double)
387 : #endif
388 : #endif /* TEST_CHOLESKY */
389 :
390 : call prepare_matrix_toeplitz(na, diagonalElement, subdiagonalElement, &
391 : d, sd, ds, sds, a, as, nblk, np_rows, &
392 1872 : np_cols, my_prow, my_pcol)
393 :
394 :
395 1872 : do_test_numeric_residual = .false.
396 : #if defined(TEST_EIGENVECTORS)
397 864 : if (nev .ge. 1) then
398 864 : do_test_numeric_residual = .true.
399 : else
400 0 : do_test_numeric_residual = .false.
401 : endif
402 : #endif
403 :
404 1872 : do_test_analytic_eigenvalues = .false.
405 1872 : do_test_analytic_eigenvalues_eigenvectors = .false.
406 1872 : do_test_frank_eigenvalues = .false.
407 : #if defined(TEST_CHOLESKY)
408 288 : do_test_toeplitz_eigenvalues = .false.
409 : #else
410 1584 : do_test_toeplitz_eigenvalues = .true.
411 : #endif
412 : #endif /* TEST_MATRIX_TOEPLITZ */
413 :
414 :
415 : #if defined(TEST_MATRIX_FRANK) && !defined(TEST_SOLVE_TRIDIAGONAL) && !defined(TEST_CHOLESKY)
416 : ! the random matrix can be used in allmost all tests; but for some no
417 : ! correctness checks have been implemented; do not allow these
418 : ! combinations
419 : ! FRANK + TEST_SOLVE_TRIDIAGONAL: we need a TOEPLITZ MATRIX
420 : ! FRANK + TEST_CHOLESKY: no correctness check yet implemented
421 :
422 : ! We also have to take care of special case in TEST_EIGENVECTORS
423 : #if !defined(TEST_EIGENVECTORS)
424 288 : call prepare_matrix_frank(na, a, z, as, nblk, np_rows, np_cols, my_prow, my_pcol)
425 :
426 288 : do_test_analytic_eigenvalues = .false.
427 288 : do_test_analytic_eigenvalues_eigenvectors = .false.
428 : #ifndef TEST_HERMITIAN_MULTIPLY
429 192 : do_test_frank_eigenvalues = .true.
430 : #endif
431 288 : do_test_toeplitz_eigenvalues = .false.
432 :
433 : #else /* TEST_EIGENVECTORS */
434 :
435 288 : if (nev .ge. 1) then
436 288 : call prepare_matrix_frank(na, a, z, as, nblk, np_rows, np_cols, my_prow, my_pcol)
437 :
438 288 : do_test_analytic_eigenvalues = .false.
439 288 : do_test_analytic_eigenvalues_eigenvectors = .false.
440 : #ifndef TEST_HERMITIAN_MULTIPLY
441 288 : do_test_frank_eigenvalues = .true.
442 : #endif
443 288 : do_test_toeplitz_eigenvalues = .false.
444 288 : do_test_numeric_residual = .false.
445 : else
446 0 : do_test_analytic_eigenvalues = .false.
447 0 : do_test_analytic_eigenvalues_eigenvectors = .false.
448 : #ifndef TEST_HERMITIAN_MULTIPLY
449 0 : do_test_frank_eigenvalues = .true.
450 : #endif
451 0 : do_test_toeplitz_eigenvalues = .false.
452 0 : do_test_numeric_residual = .false.
453 :
454 : endif
455 :
456 : #endif /* TEST_EIGENVECTORS */
457 : #endif /* (TEST_MATRIX_FRANK) */
458 : #if defined(TEST_MATRIX_FRANK) && (defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_CHOLESKY))
459 : #error "FRANK matrix is not allowed in this configuration"
460 : #endif
461 :
462 :
463 : #ifdef TEST_HERMITIAN_MULTIPLY
464 : #ifdef TEST_REAL
465 :
466 : #ifdef TEST_DOUBLE
467 192 : b(:,:) = 2.0_c_double * a(:,:)
468 192 : c(:,:) = 0.0_c_double
469 : #else
470 48 : b(:,:) = 2.0_c_float * a(:,:)
471 48 : c(:,:) = 0.0_c_float
472 : #endif
473 :
474 : #endif /* TEST_REAL */
475 :
476 : #ifdef TEST_COMPLEX
477 :
478 : #ifdef TEST_DOUBLE
479 96 : b(:,:) = 2.0_c_double * a(:,:)
480 96 : c(:,:) = (0.0_c_double, 0.0_c_double)
481 : #else
482 48 : b(:,:) = 2.0_c_float * a(:,:)
483 48 : c(:,:) = (0.0_c_float, 0.0_c_float)
484 : #endif
485 :
486 : #endif /* TEST_COMPLEX */
487 :
488 : #endif /* TEST_HERMITIAN_MULTIPLY */
489 :
490 4656 : e => elpa_allocate()
491 :
492 4656 : call e%set("na", na, error)
493 4656 : assert_elpa_ok(error)
494 4656 : call e%set("nev", nev, error)
495 4656 : assert_elpa_ok(error)
496 4656 : call e%set("local_nrows", na_rows, error)
497 4656 : assert_elpa_ok(error)
498 4656 : call e%set("local_ncols", na_cols, error)
499 4656 : assert_elpa_ok(error)
500 4656 : call e%set("nblk", nblk, error)
501 4656 : assert_elpa_ok(error)
502 :
503 : #ifdef WITH_MPI
504 3168 : call e%set("mpi_comm_parent", MPI_COMM_WORLD, error)
505 3168 : assert_elpa_ok(error)
506 3168 : call e%set("process_row", my_prow, error)
507 3168 : assert_elpa_ok(error)
508 3168 : call e%set("process_col", my_pcol, error)
509 3168 : assert_elpa_ok(error)
510 : #endif
511 4656 : call e%set("timings",1,error)
512 4656 : assert_elpa_ok(e%setup())
513 :
514 : #ifdef TEST_SOLVER_1STAGE
515 2160 : call e%set("solver", ELPA_SOLVER_1STAGE,error)
516 : #else
517 2496 : call e%set("solver", ELPA_SOLVER_2STAGE,error)
518 : #endif
519 4656 : assert_elpa_ok(error)
520 :
521 4656 : call e%set("gpu", TEST_GPU, error)
522 4656 : assert_elpa_ok(error)
523 :
524 : #if TEST_QR_DECOMPOSITION == 1
525 0 : call e%set("qr", 1, error)
526 0 : assert_elpa_ok(error)
527 : #endif
528 :
529 4656 : if (myid == 0) print *, ""
530 :
531 : #ifdef TEST_ALL_KERNELS
532 20640 : do i = 0, elpa_option_cardinality(KERNEL_KEY) ! kernels
533 19680 : kernel = elpa_option_enumerate(KERNEL_KEY, i)
534 : #endif
535 : #ifdef TEST_KERNEL
536 1344 : kernel = TEST_KERNEL
537 : #endif
538 :
539 : #ifdef TEST_SOLVER_2STAGE
540 21024 : call e%set(KERNEL_KEY, kernel, error)
541 : #ifdef TEST_KERNEL
542 1344 : assert_elpa_ok(error)
543 : #else
544 19680 : if (error /= ELPA_OK) then
545 8232 : cycle
546 : endif
547 : ! actually used kernel might be different if forced via environment variables
548 11448 : call e%get(KERNEL_KEY, kernel, error)
549 : #endif
550 12792 : if (myid == 0) then
551 8528 : print *, elpa_int_value_to_string(KERNEL_KEY, kernel) // " kernel"
552 : endif
553 : #endif
554 :
555 : #ifdef TEST_ALL_KERNELS
556 11448 : call e%timer_start(elpa_int_value_to_string(KERNEL_KEY, kernel))
557 : #endif
558 :
559 : ! The actual solve step
560 : #if defined(TEST_EIGENVECTORS)
561 : #if TEST_QR_DECOMPOSITION == 1
562 0 : call e%timer_start("e%eigenvectors_qr()")
563 : #else
564 13560 : call e%timer_start("e%eigenvectors()")
565 : #endif
566 : #ifdef TEST_SCALAPACK_ALL
567 96 : call solve_scalapack_all(na, a, sc_desc, ev, z)
568 : #elif TEST_SCALAPACK_PART
569 96 : call solve_scalapack_part(na, a, sc_desc, nev, ev, z)
570 96 : check_all_evals = .false. ! scalapack does not compute all eigenvectors
571 : #else
572 13368 : call e%eigenvectors(a, ev, z, error)
573 : #endif
574 : #if TEST_QR_DECOMPOSITION == 1
575 0 : call e%timer_stop("e%eigenvectors_qr()")
576 : #else
577 13560 : call e%timer_stop("e%eigenvectors()")
578 : #endif
579 : #endif /* TEST_EIGENVECTORS */
580 :
581 : #ifdef TEST_EIGENVALUES
582 768 : call e%timer_start("e%eigenvalues()")
583 768 : call e%eigenvalues(a, ev, error)
584 768 : call e%timer_stop("e%eigenvalues()")
585 : #endif
586 :
587 : #if defined(TEST_SOLVE_TRIDIAGONAL)
588 144 : call e%timer_start("e%solve_tridiagonal()")
589 144 : call e%solve_tridiagonal(d, sd, z, error)
590 144 : call e%timer_stop("e%solve_tridiagonal()")
591 144 : ev(:) = d(:)
592 : #endif
593 :
594 : #if defined(TEST_CHOLESKY)
595 288 : call e%timer_start("e%cholesky()")
596 288 : call e%cholesky(a, error)
597 288 : call e%timer_stop("e%cholesky()")
598 : #endif
599 :
600 : #if defined(TEST_HERMITIAN_MULTIPLY)
601 384 : call e%timer_start("e%hermitian_multiply()")
602 384 : call e%hermitian_multiply('F','F', na, a, b, na_rows, na_cols, c, na_rows, na_cols, error)
603 384 : call e%timer_stop("e%hermitian_multiply()")
604 : #endif
605 :
606 15144 : assert_elpa_ok(error)
607 :
608 : #ifdef TEST_ALL_KERNELS
609 11448 : call e%timer_stop(elpa_int_value_to_string(KERNEL_KEY, kernel))
610 : #endif
611 :
612 15144 : if (myid .eq. 0) then
613 : #ifdef TEST_ALL_KERNELS
614 7632 : call e%print_times(elpa_int_value_to_string(KERNEL_KEY, kernel))
615 : #else /* TEST_ALL_KERNELS */
616 :
617 : #if defined(TEST_EIGENVECTORS)
618 : #if TEST_QR_DECOMPOSITION == 1
619 0 : call e%print_times("e%eigenvectors_qr()")
620 : #else
621 1376 : call e%print_times("e%eigenvectors()")
622 : #endif
623 : #endif
624 : #ifdef TEST_EIGENVALUES
625 512 : call e%print_times("e%eigenvalues()")
626 : #endif
627 : #ifdef TEST_SOLVE_TRIDIAGONAL
628 96 : call e%print_times("e%solve_tridiagonal()")
629 : #endif
630 : #ifdef TEST_CHOLESKY
631 192 : call e%print_times("e%cholesky()")
632 : #endif
633 : #ifdef TEST_HERMITIAN_MULTIPLY
634 256 : call e%print_times("e%hermitian_multiply()")
635 : #endif
636 : #endif /* TEST_ALL_KERNELS */
637 : endif
638 :
639 15144 : if (do_test_analytic_eigenvalues) then
640 4152 : status = check_correctness_analytic(na, nev, ev, z, nblk, myid, np_rows, np_cols, my_prow, my_pcol, check_all_evals, .false.)
641 4152 : call check_status(status, myid)
642 : endif
643 :
644 15144 : if (do_test_analytic_eigenvalues_eigenvectors) then
645 4152 : status = check_correctness_analytic(na, nev, ev, z, nblk, myid, np_rows, np_cols, my_prow, my_pcol, check_all_evals, .true.)
646 4152 : call check_status(status, myid)
647 : endif
648 15144 : if(do_test_numeric_residual) then
649 12072 : status = check_correctness_evp_numeric_residuals(na, nev, as, z, ev, sc_desc, nblk, myid, np_rows,np_cols, my_prow, my_pcol)
650 12072 : call check_status(status, myid)
651 : endif
652 :
653 15144 : if (do_test_frank_eigenvalues) then
654 1680 : status = check_correctness_eigenvalues_frank(na, ev, z, myid)
655 1680 : call check_status(status, myid)
656 : endif
657 :
658 912 : if (do_test_toeplitz_eigenvalues) then
659 : #if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL)
660 : status = check_correctness_eigenvalues_toeplitz(na, diagonalElement, &
661 720 : subdiagonalElement, ev, z, myid)
662 720 : call check_status(status, myid)
663 : #endif
664 : endif
665 :
666 15144 : if (do_test_cholesky) then
667 288 : status = check_correctness_cholesky(na, a, as, na_rows, sc_desc, myid )
668 288 : call check_status(status, myid)
669 : endif
670 :
671 : #ifdef TEST_HERMITIAN_MULTIPLY
672 384 : if (do_test_hermitian_multiply) then
673 384 : status = check_correctness_hermitian_multiply(na, a, b, c, na_rows, sc_desc, myid )
674 384 : call check_status(status, myid)
675 : endif
676 : #endif
677 :
678 15144 : if (myid == 0) then
679 10064 : print *, ""
680 : endif
681 :
682 : #ifdef TEST_ALL_KERNELS
683 11448 : a(:,:) = as(:,:)
684 : #if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || TEST_QR_DECOMPOSITION == 1 || defined(TEST_CHOLESKY)
685 11448 : d = ds
686 11448 : sd = sds
687 : #endif
688 : end do ! kernels
689 : #endif
690 :
691 4656 : call elpa_deallocate(e)
692 :
693 4656 : deallocate(a)
694 4656 : deallocate(as)
695 4656 : deallocate(z)
696 4656 : deallocate(ev)
697 : #ifdef TEST_HERMITIAN_MULTIPLY
698 384 : deallocate(b)
699 384 : deallocate(c)
700 : #endif
701 : #if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || TEST_QR_DECOMPOSITION == 1 || defined(TEST_CHOLESKY)
702 4272 : deallocate(d, ds)
703 4272 : deallocate(sd, sds)
704 4272 : deallocate(ev_analytic)
705 : #endif
706 :
707 : #ifdef TEST_ALL_LAYOUTS
708 : end do ! factors
709 : end do ! layouts
710 : #endif
711 4656 : call elpa_uninit()
712 :
713 : #ifdef WITH_MPI
714 3168 : call blacs_gridexit(my_blacs_ctxt)
715 3168 : call mpi_finalize(mpierr)
716 : #endif
717 4656 : call exit(status)
718 :
719 : contains
720 :
721 23448 : subroutine check_status(status, myid)
722 : implicit none
723 : integer, intent(in) :: status, myid
724 : integer :: mpierr
725 23448 : if (status /= 0) then
726 0 : if (myid == 0) print *, "Result incorrect!"
727 : #ifdef WITH_MPI
728 0 : call mpi_finalize(mpierr)
729 : #endif
730 0 : call exit(status)
731 : endif
732 23448 : end subroutine
733 :
734 : end program
|