test.F90 17.9 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
!    This file is part of ELPA.
!
!    The ELPA library was originally created by the ELPA consortium,
!    consisting of the following organizations:
!
!    - Max Planck Computing and Data Facility (MPCDF), formerly known as
!      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!    - Bergische Universität Wuppertal, Lehrstuhl für angewandte
!      Informatik,
!    - Technische Universität München, Lehrstuhl für Informatik mit
!      Schwerpunkt Wissenschaftliches Rechnen ,
!    - Fritz-Haber-Institut, Berlin, Abt. Theorie,
!    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
!      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!      and
!    - IBM Deutschland GmbH
!
!
!    More information can be found here:
!    http://elpa.mpcdf.mpg.de/
!
!    ELPA is free software: you can redistribute it and/or modify
!    it under the terms of the version 3 of the license of the
!    GNU Lesser General Public License as published by the Free
!    Software Foundation.
!
!    ELPA is distributed in the hope that it will be useful,
!    but WITHOUT ANY WARRANTY; without even the implied warranty of
!    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!    GNU Lesser General Public License for more details.
!
!    You should have received a copy of the GNU Lesser General Public License
!    along with ELPA.  If not, see <http://www.gnu.org/licenses/>
!
!    ELPA reflects a substantial effort on the part of the original
!    ELPA consortium, and we ask you to respect the spirit of the
!    license that we chose: i.e., please contribute any changes you
!    may have back to the original ELPA library distribution, and keep
!    any derivatives of ELPA under the same license that we chose for
!    the original distribution, the GNU Lesser General Public License.
!
!
#include "config-f90.h"

! Define one of TEST_REAL or TEST_COMPLEX
! Define one of TEST_SINGLE or TEST_DOUBLE
! Define one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE
! Define TEST_GPU \in [0, 1]
49
! Define either TEST_ALL_KERNELS or a TEST_KERNEL \in [any valid kernel]
50

51
52
#if !(defined(TEST_REAL) ^ defined(TEST_COMPLEX))
error: define exactly one of TEST_REAL or TEST_COMPLEX
53
54
55
56
57
58
#endif

#if !(defined(TEST_SINGLE) ^ defined(TEST_DOUBLE))
error: define exactly one of TEST_SINGLE or TEST_DOUBLE
#endif

Pavel Kus's avatar
Pavel Kus committed
59
#if !(defined(TEST_SOLVER_1STAGE) ^ defined(TEST_SOLVER_2STAGE) ^ defined(TEST_SCALAPACK_ALL))
60
error: define exactly one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE or TEST_SCALAPACK_ALL
61
62
#endif

63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
#ifdef TEST_SOLVER_1STAGE
#ifdef TEST_ALL_KERNELS
error: TEST_ALL_KERNELS cannot be defined for TEST_SOLVER_1STAGE
#endif
#ifdef TEST_KERNEL
error: TEST_KERNEL cannot be defined for TEST_SOLVER_1STAGE
#endif
#endif

#ifdef TEST_SOLVER_2STAGE
#if !(defined(TEST_KERNEL) ^ defined(TEST_ALL_KERNELS))
error: define either TEST_ALL_KERNELS or a valid TEST_KERNEL
#endif
#endif


79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
#ifdef TEST_SINGLE
#  define EV_TYPE real(kind=C_FLOAT)
#  ifdef TEST_REAL
#    define MATRIX_TYPE real(kind=C_FLOAT)
#  else
#    define MATRIX_TYPE complex(kind=C_FLOAT_COMPLEX)
#  endif
#else
#  define EV_TYPE real(kind=C_DOUBLE)
#  ifdef TEST_REAL
#    define MATRIX_TYPE real(kind=C_DOUBLE)
#  else
#    define MATRIX_TYPE complex(kind=C_DOUBLE_COMPLEX)
#  endif
#endif

95
96
97
98
99
100
101
#ifdef TEST_REAL
#define KERNEL_KEY "real_kernel"
#endif
#ifdef TEST_COMPLEX
#define KERNEL_KEY "complex_kernel"
#endif

102
103
104
105
#include "assert.h"

program test
   use elpa
106
107
108
109
110
111
112

   use test_util
   use test_setup_mpi
   use test_prepare_matrix
   use test_read_input_parameters
   use test_blacs_infrastructure
   use test_check_correctness
Pavel Kus's avatar
Pavel Kus committed
113
   use test_analytic
114
#ifdef WITH_SCALAPACK_TESTS
Pavel Kus's avatar
Pavel Kus committed
115
   use test_scalapack
116
#endif
117

118
119
120
#ifdef HAVE_REDIRECT
   use test_redirect
#endif
121
122
123
   implicit none

   ! matrix dimensions
124
   integer                     :: na, nev, nblk
125
126

   ! mpi
127
128
129
130
131
   integer                     :: myid, nprocs
   integer                     :: na_cols, na_rows  ! local matrix size
   integer                     :: np_cols, np_rows  ! number of MPI processes per column/row
   integer                     :: my_prow, my_pcol  ! local MPI task position (my_prow, my_pcol) in the grid (0..np_cols -1, 0..np_rows -1)
   integer                     :: mpierr
132
133

   ! blacs
134
   integer                     :: my_blacs_ctxt, sc_desc(9), info, nprow, npcol
135
136

   ! The Matrix
137
   MATRIX_TYPE, allocatable    :: a(:,:), as(:,:)
138
#if defined(TEST_HERMITIAN_MULTIPLY)
139
   MATRIX_TYPE, allocatable    :: b(:,:), c(:,:)
140
#endif
141
   ! eigenvectors
142
   MATRIX_TYPE, allocatable    :: z(:,:)
143
   ! eigenvalues
144
   EV_TYPE, allocatable        :: ev(:), ev_analytic(:)
145
146

#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || defined(TEST_HERMITIAN_MULTIPLY)
147
148
   EV_TYPE, allocatable        :: d(:), sd(:), ds(:), sds(:)
   EV_TYPE                     :: diagonalELement, subdiagonalElement
149
#endif
150
#if defined(TEST_CHOLESKY)
151
152
   MATRIX_TYPE, allocatable    :: d(:), sd(:), ds(:), sds(:)
   MATRIX_TYPE                 :: diagonalELement, subdiagonalElement
153
#endif
154

155
   integer                     :: error, status
156

157
158
   type(output_t)              :: write_to_file
   class(elpa_t), pointer      :: e
159
#ifdef TEST_ALL_KERNELS
160
   integer                     :: i
161
162
163
#endif
#ifdef TEST_ALL_LAYOUTS
   character(len=1), parameter :: layouts(2) = [ 'C', 'R' ]
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
   integer                     :: i_layout
#endif
   integer                     :: kernel
   character(len=1)            :: layout
#ifdef TEST_COMPLEX
   EV_TYPE                     :: norm, normmax
   MATRIX_TYPE, allocatable    :: tmp1(:,:), tmp2(:,:)
#ifdef TEST_DOUBLE
   MATRIX_TYPE, parameter      :: CONE = (1.0_c_double, 0.0_c_double), &
                                  CZERO = (0.0_c_double, 0.0_c_double)
   EV_TYPE :: pzlange, zlange
#else
   MATRIX_TYPE, parameter      :: CONE = (1.0_c_float, 0.0_c_float), &
                                  CZERO = (0.0_c_float, 0.0_c_float)
   EV_TYPE :: pclange, clange
#endif
#endif
181
   call read_input_parameters_traditional(na, nev, nblk, write_to_file)
182
   call setup_mpi(myid, nprocs)
183
184
185
186
#ifdef HAVE_REDIRECT
#ifdef WITH_MPI
     call MPI_BARRIER(MPI_COMM_WORLD, mpierr)
     call redirect_stdout(myid)
187
#endif
188
#endif
189

190
191
192
193
194
   if (elpa_init(CURRENT_API_VERSION) /= ELPA_OK) then
     print *, "ELPA API version not supported"
     stop 1
   endif

195
196
197
198
199
   if (myid == 0) then
     print '((a,i0))', 'Program ' // TEST_CASE
     print *, ""
   endif

200
#ifdef TEST_ALL_LAYOUTS
201
   do i_layout = 1, size(layouts)               ! layouts
202
     layout = layouts(i_layout)
203
     do np_cols = 1, nprocs                     ! factors
204
205
206
207
208
       if (mod(nprocs,np_cols) /= 0 ) then
         cycle
       endif
#else
   layout = 'C'
209
210
211
   do np_cols = NINT(SQRT(REAL(nprocs))),2,-1
      if(mod(nprocs,np_cols) == 0 ) exit
   enddo
212
#endif
213
214

   np_rows = nprocs/np_cols
215
   assert(nprocs == np_rows * np_cols)
216

217
218
219
220
   if (myid == 0) then
     print '((a,i0))', 'Matrix size: ', na
     print '((a,i0))', 'Num eigenvectors: ', nev
     print '((a,i0))', 'Blocksize: ', nblk
221
#ifdef WITH_MPI
222
223
     print '((a,i0))', 'Num MPI proc: ', nprocs
     print '(3(a,i0))','Number of processor rows=',np_rows,', cols=',np_cols,', total=',nprocs
224
     print '(a)',      'Process layout: ' // layout
225
#endif
226
227
228
     print *,''
   endif

229
230
   call set_up_blacsgrid(mpi_comm_world, np_rows, np_cols, layout, &
                         my_blacs_ctxt, my_prow, my_pcol)
231
232
233
234

   call set_up_blacs_descriptor(na, nblk, my_prow, my_pcol, np_rows, np_cols, &
                                na_rows, na_cols, sc_desc, my_blacs_ctxt, info)

Pavel Kus's avatar
Pavel Kus committed
235
236
   allocate(a (na_rows,na_cols))
   allocate(as(na_rows,na_cols))
237
238
239
   allocate(z (na_rows,na_cols))
   allocate(ev(na))

240
241
242
243
244
#ifdef TEST_HERMITIAN_MULTIPLY
   allocate(b (na_rows,na_cols))
   allocate(c (na_rows,na_cols))
#endif

245
#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || defined(TEST_CHOLESKY)
246
247
248
249
250
   allocate(d (na), ds(na))
   allocate(sd (na), sds(na))
   allocate(ev_analytic(na))
#endif

251
252
253
254
   a(:,:) = 0.0
   z(:,:) = 0.0
   ev(:) = 0.0

255
#if defined(TEST_EIGENVECTORS) || defined(TEST_HERMITIAN_MULTIPLY)
Pavel Kus's avatar
Pavel Kus committed
256
257
#ifdef TEST_MATRIX_ANALYTIC
   call prepare_matrix_analytic(na, a, nblk, myid, np_rows, np_cols, my_prow, my_pcol)
258
   as(:,:) = a
259
260
261
262
263
264
265
266
#else
   if (nev .ge. 1) then
     call prepare_matrix(na, myid, sc_desc, a, z, as)
   else
     ! zero eigenvectors and not analytic test => toeplitz matrix
#ifdef TEST_SINGLE
     diagonalElement = 0.45_c_float
     subdiagonalElement =  0.78_c_float
Pavel Kus's avatar
Pavel Kus committed
267
#else
268
269
270
271
272
273
274
     diagonalElement = 0.45_c_double
     subdiagonalElement =  0.78_c_double
#endif
     call prepare_toeplitz_matrix(na, diagonalElement, subdiagonalElement, &
                                  d, sd, ds, sds, a, as, nblk, np_rows, &
                                  np_cols, my_prow, my_pcol)
   endif
275
276
277
278
279
280
281
282
283
284

#ifdef TEST_HERMITIAN_MULTIPLY
#if REALCASE == 1

#ifdef DOUBLE_PRECISION_REAL
   b(:,:) = 2.0_rk8 * a(:,:)
   c(:,:) = 0.0_rk8
#else
   b(:,:) = 2.0_rk4 * a(:,:)
   c(:,:) = 0.0_rk4
285
#endif
286

Pavel Kus's avatar
Pavel Kus committed
287
#endif
288

289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
#if COMPLEXCASE == 1

#ifdef DOUBLE_PRECISION_COMPLEX
   b(:,:) = 2.0_ck8 * a(:,:)
   c(:,:) = 0.0_ck8
#else
   b(:,:) = 2.0_ck4 * a(:,:)
   c(:,:) = 0.0_ck4
#endif

#endif

#endif /* TEST_HERMITIAN_MULTIPLY */

#endif /* TEST_MATRIX_ANALYTIC */
#endif /* defined(TEST_EIGENVECTORS) || defined(TEST_HERMITIAN_MULTIPLY) */

306
#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL)
Andreas Marek's avatar
Andreas Marek committed
307
308
309
310
311
312
313
314
315
316
317

#ifdef TEST_SINGLE
   diagonalElement = 0.45_c_float
   subdiagonalElement =  0.78_c_float
#else
   diagonalElement = 0.45_c_double
   subdiagonalElement =  0.78_c_double
#endif
   call prepare_toeplitz_matrix(na, diagonalElement, subdiagonalElement, &
                                d, sd, ds, sds, a, as, nblk, np_rows, &
                                np_cols, my_prow, my_pcol)
318
#endif /* EIGENVALUES OR TRIDIAGONAL */
319

320
321
322
#if defined(TEST_CHOLESKY)

#ifdef TEST_SINGLE
323
324
   diagonalElement = (2.546_c_float, 0.0_c_float)
   subdiagonalElement =  (0.0_c_float, 0.0_c_float)
325
#else
326
327
   diagonalElement = (2.546_c_double, 0.0_c_double)
   subdiagonalElement =  (0.0_c_double, 0.0_c_double)
328
329
330
331
#endif
   call prepare_toeplitz_matrix(na, diagonalElement, subdiagonalElement, &
                                d, sd, ds, sds, a, as, nblk, np_rows, &
                                np_cols, my_prow, my_pcol)
332
333


334
335
#endif /* TEST_CHOLESKY */

336
337
   e => elpa_allocate()

338
339
340
341
342
343
344
345
346
347
   call e%set("na", na, error)
   assert_elpa_ok(error)
   call e%set("nev", nev, error)
   assert_elpa_ok(error)
   call e%set("local_nrows", na_rows, error)
   assert_elpa_ok(error)
   call e%set("local_ncols", na_cols, error)
   assert_elpa_ok(error)
   call e%set("nblk", nblk, error)
   assert_elpa_ok(error)
348
349

#ifdef WITH_MPI
350
351
352
353
354
355
   call e%set("mpi_comm_parent", MPI_COMM_WORLD, error)
   assert_elpa_ok(error)
   call e%set("process_row", my_prow, error)
   assert_elpa_ok(error)
   call e%set("process_col", my_pcol, error)
   assert_elpa_ok(error)
356
#endif
357

Andreas Marek's avatar
Andreas Marek committed
358
359
   call e%set("timings",1)

360
361
362
363
364
365
366
   assert_elpa_ok(e%setup())

#ifdef TEST_SOLVER_1STAGE
   call e%set("solver", ELPA_SOLVER_1STAGE)
#else
   call e%set("solver", ELPA_SOLVER_2STAGE)
#endif
367
   assert_elpa_ok(error)
368

369
370
   call e%set("gpu", TEST_GPU, error)
   assert_elpa_ok(error)
371

372
373
   if (myid == 0) print *, ""

374
#ifdef TEST_ALL_KERNELS
375
   do i = 0, elpa_option_cardinality(KERNEL_KEY)  ! kernels
376
     kernel = elpa_option_enumerate(KERNEL_KEY, i)
377
#endif
378
#ifdef TEST_KERNEL
379
     kernel = TEST_KERNEL
380
#endif
381
382

#ifdef TEST_SOLVER_2STAGE
383
     call e%set(KERNEL_KEY, kernel, error)
384
385
386
#ifdef TEST_KERNEL
     assert_elpa_ok(error)
#else
387
388
389
     if (error /= ELPA_OK) then
       cycle
     endif
390
391
392
393
394
395
     ! actually used kernel might be different if forced via environment variables
     call e%get(KERNEL_KEY, kernel)
#endif
     if (myid == 0) then
       print *, elpa_int_value_to_string(KERNEL_KEY, kernel) // " kernel"
     endif
396
397
#endif

398
#ifdef TEST_ALL_KERNELS
399
     call e%timer_start(elpa_int_value_to_string(KERNEL_KEY, kernel))
400
#endif
401

402
     ! The actual solve step
403
404
#ifdef TEST_EIGENVECTORS
     call e%timer_start("e%eigenvectors()")
Pavel Kus's avatar
Pavel Kus committed
405
406
407
#ifdef TEST_SCALAPACK_ALL
     call solve_scalapack_all(na, a, sc_desc, ev, z)
#else
408
     call e%eigenvectors(a, ev, z, error)
Pavel Kus's avatar
Pavel Kus committed
409
#endif
410
     call e%timer_stop("e%eigenvectors()")
411
#endif /* TEST_EIGENVECTORS */
412
413
414

#ifdef TEST_EIGENVALUES
     call e%timer_start("e%eigenvalues()")
415
     call e%eigenvalues(a, ev, error)
416
     call e%timer_stop("e%eigenvalues()")
417
#endif
418
419
420

#if defined(TEST_SOLVE_TRIDIAGONAL)
     call e%timer_start("e%solve_tridiagonal()")
421
     call e%solve_tridiagonal(d, sd, z, error)
422
     call e%timer_stop("e%solve_tridiagonal()")
423
424
425
     ev(:) = d(:)
#endif

426
427
428
429
430
431
#if defined(TEST_CHOLESKY)
     call e%timer_start("e%cholesky()")
     call e%cholesky(a, error)
     call e%timer_stop("e%cholesky()")
#endif

432
433
434
435
436
#if defined(TEST_HERMITIAN_MULTIPLY)
     call e%timer_start("e%hermitian_multiply()")
     call e%hermitian_multiply('F','F', na, a, b, na_rows, na_cols, c, na_rows, na_cols, error)
     call e%timer_stop("e%hermitian_multiply()")
#endif
Pavel Kus's avatar
Pavel Kus committed
437

438
439
     assert_elpa_ok(error)

440
#ifdef TEST_ALL_KERNELS
441
442
443
     call e%timer_stop(elpa_int_value_to_string(KERNEL_KEY, kernel))
#endif

444
     if (myid .eq. 0) then
445
#ifdef TEST_ALL_KERNELS
446
       call e%print_times(elpa_int_value_to_string(KERNEL_KEY, kernel))
447
448
#else /* TEST_ALL_KERNELS */

449
#ifdef TEST_EIGENVECTORS
450
       call e%print_times("e%eigenvectors()")
451
#endif
452
#ifdef TEST_EIGENVALUES
453
454
       call e%print_times("e%eigenvalues()")
#endif
455
456
#ifdef TEST_SOLVE_TRIDIAGONAL
       call e%print_times("e%solve_tridiagonal()")
457
#endif
458
459
#ifdef TEST_CHOLESKY
       call e%print_times("e%cholesky()")
460
#endif
461
462
463
#ifdef TEST_HERMITIAN_MULTIPLY
       call e%print_times("e%hermitian_multiply()")
#endif
464
#endif /* TEST_ALL_KERNELS */
465
     endif
466

467
#ifdef TEST_EIGENVECTORS
Pavel Kus's avatar
Pavel Kus committed
468
469
470
#ifdef TEST_MATRIX_ANALYTIC
     status = check_correctness_analytic(na, nev, ev, z, nblk, myid, np_rows, np_cols, my_prow, my_pcol)
#else
471
472
473
474
475
476
477
     if (nev .ge. 1) then
       status = check_correctness(na, nev, as, z, ev, sc_desc, nblk, myid, np_rows,np_cols, my_prow, my_pcol)
     else
       ! zero eigenvectors and no analytic test => toeplitz
       status = check_correctness_eigenvalues_toeplitz(na, diagonalElement, &
         subdiagonalElement, ev, z, myid)
     endif
Andreas Marek's avatar
Andreas Marek committed
478
#endif
479
     call check_status(status, myid)
480
481
#endif

482
#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL)
Andreas Marek's avatar
Andreas Marek committed
483
484
     status = check_correctness_eigenvalues_toeplitz(na, diagonalElement, &
         subdiagonalElement, ev, z, myid)
485
     call check_status(status, myid)
486

487
#ifdef TEST_SOLVE_TRIDIAGONAL
488
     ! check eigenvectors
Pavel Kus's avatar
Pavel Kus committed
489
     status = check_correctness(na, nev, as, z, ev, sc_desc, nblk, myid, np_rows, np_cols, my_prow, my_pcol)
490
     call check_status(status, myid)
491
492
#endif
#endif
493

494
495
496
497
498
#if defined(TEST_CHOLESKY)
     status = check_correctness_cholesky(na, a, as, na_rows, sc_desc, myid )
     call check_status(status, myid)
#endif

499
#if defined(TEST_HERMITIAN_MULTIPLY)
500
501

#ifdef TEST_REAL
502
503
504
     status = check_correctness_hermitian_multiply(na, a, b, c, na_rows, sc_desc, myid )
     call check_status(status, myid)
#endif
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
#ifdef TEST_COMPLEX
   status = 0

   !-------------------------------------------------------------------------------
   ! Test correctness of result (using plain scalapack routines)
   allocate(tmp1(na_rows,na_cols))
   allocate(tmp2(na_rows,na_cols))
#ifdef TEST_DOUBLE
   tmp1(:,:) = 0.0_ck8
#else
   tmp1(:,:) = 0.0_ck4
#endif
   ! tmp1 = a**T
#ifdef WITH_MPI
#ifdef TEST_DOUBLE
   call pztranc(na, na, CONE, a, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#else
   call pctranc(na, na, CONE, a, 1, 1, sc_desc, CZERO, tmp1, 1, 1, sc_desc)
#endif
#else
   tmp1 = transpose(conjg(a))
#endif
   ! tmp2 = tmp1 * b
#ifdef TEST_DOUBLE
#ifdef WITH_MPI
   call pzgemm("N","N", na, na, na, CONE, tmp1, 1, 1, sc_desc, b, 1, 1, &
               sc_desc, CZERO, tmp2, 1, 1, sc_desc)
#else
   call zgemm("N","N", na, na, na, CONE, tmp1, na, b, na, CZERO, tmp2, na)
#endif
#else
#ifdef WITH_MPI
   call pcgemm("N","N", na, na, na, CONE, tmp1, 1, 1, sc_desc, b, 1, 1, &
               sc_desc, CZERO, tmp2, 1, 1, sc_desc)
#else
   call cgemm("N","N", na, na, na, CONE, tmp1, na, b, na, CZERO, tmp2, na)
#endif
#endif

   ! compare tmp2 with c
   tmp2(:,:) = tmp2(:,:) - c(:,:)
#ifdef TEST_DOUBLE
#ifdef WITH_MPI
   norm = pzlange("M",na, na, tmp2, 1, 1, sc_desc, tmp1)
#else
   norm = zlange("M",na, na, tmp2, na_rows, tmp1)
#endif
#else
#ifdef WITH_MPI
   norm = pclange("M",na, na, tmp2, 1, 1, sc_desc, tmp1)
#else
   norm = clange("M",na, na, tmp2, na_rows, tmp1)
#endif
#endif
#ifdef WITH_MPI
#ifdef TEST_DOUBLE
   call mpi_allreduce(norm,normmax,1,MPI_REAL8,MPI_MAX,MPI_COMM_WORLD,mpierr)
#else
   call mpi_allreduce(norm,normmax,1,MPI_REAL4,MPI_MAX,MPI_COMM_WORLD,mpierr)
#endif
#else
   normmax = norm
#endif
   if (myid .eq. 0) then
     print *," Maximum error of result: ", normmax
   endif

#ifdef TEST_DOUBLE
Andreas Marek's avatar
Andreas Marek committed
573
   if (normmax .gt. 5e-10_rk8) then
574
#else
Andreas Marek's avatar
Andreas Marek committed
575
   if (normmax .gt. 5e-3_rk4) then
576
577
578
579
580
581
582
583
584
585
#endif
        print *,"norm= ",normmax
        status = 1
   endif

   deallocate(tmp1)
   deallocate(tmp2)

#endif
#endif /* TEST_HERMITIAN_MULTIPLY */
586
587


588
589
590
591
     if (myid == 0) then
       print *, ""
     endif

592
593
#ifdef TEST_ALL_KERNELS
     a(:,:) = as(:,:)
594
#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || defined(TEST_CHOLESKY)
595
596
597
     d = ds
     sd = sds
#endif
598
   end do ! kernels
599
#endif
Andreas Marek's avatar
Andreas Marek committed
600

601
602
603
604
605
606
607
   call elpa_deallocate(e)

   deallocate(a)
   deallocate(as)
   deallocate(z)
   deallocate(ev)

608
609
610
611
612
#ifdef TEST_HERMITIAN_MULTIPLY
   deallocate(b)
   deallocate(c)
#endif

613
#if defined(TEST_EIGENVALUES) || defined(TEST_SOLVE_TRIDIAGONAL) || defined(TEST_EIGENVECTORS) || defined(TEST_CHOLESKY)
614
615
616
617
618
   deallocate(d, ds)
   deallocate(sd, sds)
   deallocate(ev_analytic)
#endif

619
620
621
622
623
624
625
#ifdef TEST_ALL_LAYOUTS
   end do ! factors
   end do ! layouts
#endif

   call elpa_uninit()

626
627
628
629
630
631
632
#ifdef WITH_MPI
   call blacs_gridexit(my_blacs_ctxt)
   call mpi_finalize(mpierr)
#endif

   call exit(status)

633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
   contains

     subroutine check_status(status, myid)
       implicit none
       integer, intent(in) :: status, myid
       integer :: mpierr
       if (status /= 0) then
         if (myid == 0) print *, "Result incorrect!"
#ifdef WITH_MPI
         call mpi_finalize(mpierr)
#endif
         call exit(status)
       endif
     end subroutine

648
end program