elpa_impl.F90 46.8 KB
Newer Older
1
2
3
!
!    Copyright 2017, L. Hüdepohl and A. Marek, MPCDF
!
Andreas Marek's avatar
Andreas Marek committed
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
!    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
!
!    This particular source code file contains additions, changes and
!    enhancements authored by Intel Corporation which is not part of
!    the ELPA consortium.
!
!    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.
!
48
#include "config-f90.h"
49

50
51
module elpa_impl
  use elpa_api
52
  use, intrinsic :: iso_c_binding
53
  implicit none
54

55
56
  private
  public :: elpa_impl_allocate
57

58
  type, extends(elpa_t) :: elpa_impl_t
Andreas Marek's avatar
Andreas Marek committed
59
   private
60
   type(c_ptr)         :: index = C_NULL_PTR
61

62
   contains
63
     ! con-/destructor
64
     procedure, public :: setup => elpa_setup
65
     procedure, public :: destroy => elpa_destroy
66

67
     ! KV store
68
     procedure, public :: get => elpa_get_integer
69
70
     procedure, public :: get_double => elpa_get_double
     procedure, public :: is_set => elpa_is_set
71
     procedure, public :: can_set => elpa_can_set
72
73

     ! privates:
74
75
76
     procedure, private :: get_real_kernel => elpa_get_real_kernel
     procedure, private :: get_complex_kernel => elpa_get_complex_kernel

77
     procedure, private :: elpa_set_integer
78
     procedure, private :: elpa_set_double
79

80
     procedure, private :: elpa_solve_real_double
81
     procedure, private :: elpa_solve_real_single
82
     procedure, private :: elpa_solve_complex_double
83
     procedure, private :: elpa_solve_complex_single
84

85
86
87
88
     procedure, private :: elpa_multiply_at_b_double
     procedure, private :: elpa_multiply_at_b_single
     procedure, private :: elpa_multiply_ah_b_double
     procedure, private :: elpa_multiply_ah_b_single
89

90
91
92
93
     procedure, private :: elpa_cholesky_double_real
     procedure, private :: elpa_cholesky_single_real
     procedure, private :: elpa_cholesky_double_complex
     procedure, private :: elpa_cholesky_single_complex
94
95
96
97
98

     procedure, private :: elpa_invert_trm_double_real
     procedure, private :: elpa_invert_trm_single_real
     procedure, private :: elpa_invert_trm_double_complex
     procedure, private :: elpa_invert_trm_single_complex
99
100
101

     procedure, private :: elpa_solve_tridi_double_real
     procedure, private :: elpa_solve_tridi_single_real
102
103
104

     procedure, private :: associate_int => elpa_associate_int

105
  end type elpa_impl_t
106
107
108

  contains

109
    function elpa_impl_allocate(error) result(obj)
Andreas Marek's avatar
Andreas Marek committed
110
111
      use precision
      use elpa_utilities, only : error_unit
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
112
      use elpa_generated_fortran_interfaces
Andreas Marek's avatar
Andreas Marek committed
113

114
115
116
117
      type(elpa_impl_t), pointer   :: obj
      integer, optional            :: error

      allocate(obj)
Andreas Marek's avatar
Andreas Marek committed
118

Andreas Marek's avatar
Andreas Marek committed
119
120
      ! check whether init has ever been called
      if (.not.(elpa_initialized())) then
121
        write(error_unit, *) "elpa_allocate(): you must call elpa_init() once before creating instances of ELPA"
122
123
        if(present(error)) then
          error = ELPA_ERROR
124
        endif
Andreas Marek's avatar
Andreas Marek committed
125
126
        return
      endif
Andreas Marek's avatar
Andreas Marek committed
127

128
      obj%index = elpa_index_instance_c()
129
130

      ! Associate some important integer pointers for convenience
131
132
133
134
135
136
137
138
      obj%na => obj%associate_int("na")
      obj%nev => obj%associate_int("nev")
      obj%local_nrows => obj%associate_int("local_nrows")
      obj%local_ncols => obj%associate_int("local_ncols")
      obj%nblk => obj%associate_int("nblk")

      if(present(error)) then
        error = ELPA_OK
139
140
      endif
    end function
Andreas Marek's avatar
Andreas Marek committed
141

142

143
    function elpa_setup(self) result(error)
144
      use elpa1_impl, only : elpa_get_communicators_impl
145
      class(elpa_impl_t), intent(inout) :: self
146
      integer :: error, error2
147
      integer :: mpi_comm_rows, mpi_comm_cols, mpierr
148

149
      error = ELPA_ERROR
150

151
152
153
      if (self%is_set("mpi_comm_parent") == 1 .and. &
          self%is_set("process_row") == 1 .and. &
          self%is_set("process_col") == 1) then
154

155
156
157
158
159
160
        mpierr = elpa_get_communicators_impl(&
                        self%get("mpi_comm_parent"), &
                        self%get("process_row"), &
                        self%get("process_col"), &
                        mpi_comm_rows, &
                        mpi_comm_cols)
161

162
163
164
        call self%set("mpi_comm_rows", mpi_comm_rows)
        call self%set("mpi_comm_cols", mpi_comm_cols)

165
        error = ELPA_OK
166
      endif
167

168
169
      if (self%is_set("mpi_comm_rows") == 1 .and. self%is_set("mpi_comm_cols") == 1) then
        error = ELPA_OK
170
      endif
171
172

    end function
173

174
    subroutine elpa_set_integer(self, name, value, error)
175
      use iso_c_binding
176
177
      use elpa_generated_fortran_interfaces
      use elpa_utilities, only : error_unit
178
      class(elpa_impl_t)              :: self
179
180
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
181
182
      integer, optional               :: error
      integer                         :: actual_error
183

184
      actual_error = elpa_index_set_int_value_c(self%index, name // c_null_char, value)
185

186
187
      if (present(error)) then
        error = actual_error
188

189
      else if (actual_error /= ELPA_OK) then
190
191
        write(error_unit,'(a,a,a,i0,a)') "ELPA: Error setting option '", name, "' to value ", value, &
                " and you did not check for errors!"
192
193

      end if
194
195
    end subroutine

196

197
    function elpa_get_integer(self, name, error) result(value)
198
      use iso_c_binding
199
      use elpa_generated_fortran_interfaces
200
      class(elpa_impl_t)             :: self
201
202
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
203
      integer, intent(out), optional :: error
204

205
      value = elpa_index_get_int_value_c(self%index, name // c_null_char, error)
206

207
    end function
Andreas Marek's avatar
Andreas Marek committed
208

209
210

    function elpa_is_set(self, name) result(state)
211
212
      use iso_c_binding
      use elpa_generated_fortran_interfaces
213
      class(elpa_impl_t)       :: self
214
      character(*), intent(in) :: name
215
      integer                  :: state
216

217
      state = elpa_index_value_is_set_c(self%index, name // c_null_char)
218
219
    end function

220

221
222
223
224
225
226
227
228
229
230
231
232
233
    function elpa_can_set(self, name, value) result(error)
      use iso_c_binding
      use elpa_generated_fortran_interfaces
      class(elpa_impl_t)       :: self
      character(*), intent(in) :: name
      integer(kind=c_int), intent(in) :: value
      integer                  :: error

      error = elpa_index_int_is_valid_c(self%index, name // c_null_char, value)
    end function


    function elpa_value_to_string(self, option_name, error) result(string)
234
235
236
      use elpa_generated_fortran_interfaces
      class(elpa_impl_t), intent(in) :: self
      character(kind=c_char, len=*), intent(in) :: option_name
237
238
239
240
      type(c_ptr) :: ptr
      integer, intent(out), optional :: error
      integer :: val, actual_error
      character(kind=c_char, len=elpa_index_int_value_to_strlen_c(self%index, option_name // C_NULL_CHAR)), pointer :: string
241

242
243
244
245
246
247
248
249
      nullify(string)

      val = self%get(option_name, actual_error)
      if (actual_error /= ELPA_OK) then
        if (present(error)) then
          error = actual_error
        endif
        return
250
251
      endif

252
253
254
255
      actual_error = elpa_int_value_to_string_c(option_name // C_NULL_CHAR, val, ptr)
      if (c_associated(ptr)) then
        call c_f_pointer(ptr, string)
      endif
256

257
258
259
260
      if (present(error)) then
        error = actual_error
      endif
    end function
261

262
263

    subroutine elpa_set_double(self, name, value, error)
Andreas Marek's avatar
Andreas Marek committed
264
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
265
      use elpa_generated_fortran_interfaces
266
      use elpa_utilities, only : error_unit
267
      class(elpa_impl_t)              :: self
268
      character(*), intent(in)        :: name
269
      real(kind=c_double), intent(in) :: value
270
271
      integer, optional               :: error
      integer                         :: actual_error
Andreas Marek's avatar
Andreas Marek committed
272

273
      actual_error = elpa_index_set_double_value_c(self%index, name // c_null_char, value)
Andreas Marek's avatar
Andreas Marek committed
274

275
276
      if (present(error)) then
        error = actual_error
277

278
      else if (actual_error /= ELPA_OK) then
279
        write(error_unit,'(a,a,es12.5,a)') "ELPA: Error setting option '", name, "' to value ", value, &
280
                " and you did not check for errors!"
281
282
283

      end if
    end subroutine
Andreas Marek's avatar
Andreas Marek committed
284
285


286
    function elpa_get_double(self, name, error) result(value)
Andreas Marek's avatar
Andreas Marek committed
287
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
288
      use elpa_generated_fortran_interfaces
289
      class(elpa_impl_t)             :: self
290
      character(*), intent(in)       :: name
291
      real(kind=c_double)            :: value
292
      integer, intent(out), optional :: error
293

294
      value = elpa_index_get_double_value_c(self%index, name // c_null_char, error)
295
296

    end function
Andreas Marek's avatar
Andreas Marek committed
297
298


299
    function elpa_associate_int(self, name) result(value)
Andreas Marek's avatar
Andreas Marek committed
300
      use iso_c_binding
301
      use elpa_generated_fortran_interfaces
302
303
      use elpa_utilities, only : error_unit
      class(elpa_impl_t)             :: self
304
305
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
Andreas Marek's avatar
Andreas Marek committed
306

307
308
      type(c_ptr)                    :: value_p

309
      value_p = elpa_index_get_int_loc_c(self%index, name // c_null_char)
310
311
312
      if (.not. c_associated(value_p)) then
        write(error_unit, '(a,a,a)') "ELPA: Warning, received NULL pointer for entry '", name, "'"
      endif
313
314
      call c_f_pointer(value_p, value)
    end function
Andreas Marek's avatar
Andreas Marek committed
315

316

317
    subroutine elpa_solve_real_double(self, a, ev, q, error)
318
319
      use elpa2_impl
      use elpa1_impl
320
      use elpa_utilities, only : error_unit
321
      use precision
Andreas Marek's avatar
Andreas Marek committed
322
      use iso_c_binding
323
      class(elpa_impl_t)  :: self
Andreas Marek's avatar
Andreas Marek committed
324

325
326
327
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
328
      real(kind=c_double) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
329
#endif
330
      real(kind=c_double) :: ev(self%na)
331
332

      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
333
334
      integer, optional   :: error
      integer(kind=c_int) :: error_actual
335
      logical             :: success_l, summary_timings
336

337
      logical             :: useGPU, useQR
338
      integer(kind=c_int) :: kernel
339

340
341
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
342
343
344
345
346
347
348
349
350
351
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif


352
353
      if (self%get("gpu",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
354
355
356
357
358
359
360
361
362
          print *,"Could not querry gpu"
          stop
        endif

        useGPU = .true.
      else
        useGPU = .false.
      endif

363
      useQR = self%get("qr") == 1
364

365
366
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
        success_l = elpa_solve_evp_real_1stage_double_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
367
368
369
370
                                                           self%local_nrows,  self%nblk, self%local_ncols, &
                                                           self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), &
                                                           self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
                                                           time_evp_solve, time_evp_back, summary_timings)
371

372
373
374
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
        kernel = self%get_real_kernel()
        success_l = elpa_solve_evp_real_2stage_double_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
375
376
377
378
379
                                                           self%local_nrows,  self%nblk, self%local_ncols, &
                                                           self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), &
                                                           self%get("mpi_comm_parent"), time_evp_fwd,     &
                                                           time_evp_solve, time_evp_back, summary_timings, useGPU, &
                                                           kernel, useQR)
380
381
382
383
      else
        print *,"unknown solver"
        stop
      endif
384

385
      if (present(error)) then
386
        if (success_l) then
387
          error = ELPA_OK
388
        else
389
          error = ELPA_ERROR
390
391
392
393
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
Andreas Marek's avatar
Andreas Marek committed
394

395
396
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
397
398
399
400
          print *,"Could not querry summary timings"
          stop
        endif

401
402
403
        call self%set("time_evp_fwd", time_evp_fwd)
        call self%set("time_evp_solve", time_evp_solve)
        call self%set("time_evp_back", time_evp_back)
404
405
      else

406
407
408
        call self%set("time_evp_fwd", -1.0_rk8)
        call self%set("time_evp_solve", -1.0_rk8)
        call self%set("time_evp_back", -1.0_rk8)
409
      endif
410
411
    end subroutine

412

413
    subroutine elpa_solve_real_single(self, a, ev, q, error)
414
415
      use elpa2_impl
      use elpa1_impl
416
      use elpa_utilities, only : error_unit
417
      use precision
418
      use iso_c_binding
419
      class(elpa_impl_t)  :: self
420
421
422
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)  :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
423
      real(kind=c_float)  :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
424
#endif
425
      real(kind=c_float)  :: ev(self%na)
426
427

      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
428
429
      integer, optional   :: error
      integer(kind=c_int) :: error_actual
430
      logical             :: success_l, summary_timings
431

432
      logical             :: useGPU, useQR
433
      integer(kind=c_int) :: kernel
434

435
#ifdef WANT_SINGLE_PRECISION_REAL
436
437
      if (self%get("timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
438
439
440
441
442
443
444
445
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif
446

447
448
      if (self%get("gpu",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
449
450
451
452
453
454
455
456
457
          print *,"Could not querry gpu"
          stop
        endif

        useGPU = .true.
      else
        useGPU = .false.
      endif

458
      useQR = self%get("qr") == 1
459

460
461
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
        success_l = elpa_solve_evp_real_1stage_single_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
462
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
463
464
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
465
                                                          time_evp_solve, time_evp_back, summary_timings)
466

467
468
469
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
        kernel = self%get_real_kernel()
        success_l = elpa_solve_evp_real_2stage_single_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
470
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
471
472
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), time_evp_fwd,     &
473
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
474
                                                          kernel, useQR)
475
476
477
478
      else
        print *,"unknown solver"
        stop
      endif
479

480
      if (present(error)) then
481
        if (success_l) then
482
          error = ELPA_OK
483
        else
484
          error = ELPA_ERROR
485
486
487
488
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
489
490


491
492
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
493
494
495
496
          print *,"Could not querry summary timings"
          stop
        endif

497
498
499
        call self%set("time_evp_fwd", time_evp_fwd)
        call self%set("time_evp_solve", time_evp_solve)
        call self%set("time_evp_back", time_evp_back)
500
501
      else

502
503
504
        call self%set("time_evp_fwd", -1.0_rk8)
        call self%set("time_evp_solve", -1.0_rk8)
        call self%set("time_evp_back", -1.0_rk8)
505
      endif
506
#else
507
      print *,"This installation of the ELPA library has not been build with single-precision support"
508
      error = ELPA_ERROR
509
510
511
#endif
    end subroutine

512

513
    subroutine elpa_solve_complex_double(self, a, ev, q, error)
514
515
      use elpa2_impl
      use elpa1_impl
516
      use elpa_utilities, only : error_unit
517
      use precision
518
      use iso_c_binding
519
      class(elpa_impl_t)             :: self
520

521
522
523
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
524
      complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
525
#endif
526
      real(kind=c_double)            :: ev(self%na)
527

528
529
      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back

530
531
      integer, optional              :: error
      integer(kind=c_int)            :: error_actual
532
      logical                        :: success_l, summary_timings
533

534
      logical                        :: useGPU
535
      integer(kind=c_int) :: kernel
536
537
      if (self%get("timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
538
539
540
541
542
543
544
545
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif
546

547
548
      if (self%get("gpu",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
549
550
551
552
553
554
555
556
557
          print *,"Could not querry gpu"
          stop
        endif

        useGPU = .true.
      else
        useGPU = .false.
      endif

558
559
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
        success_l = elpa_solve_evp_complex_1stage_double_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
560
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
561
562
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
563
                                                          time_evp_solve, time_evp_back, summary_timings)
564

565
566
567
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
        kernel = self%get_complex_kernel()
        success_l = elpa_solve_evp_complex_2stage_double_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
568
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
569
570
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), time_evp_fwd,     &
571
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
572
                                                          kernel)
573
574
575
576
      else
        print *,"unknown solver"
        stop
      endif
577

578
      if (present(error)) then
579
        if (success_l) then
580
          error = ELPA_OK
581
        else
582
          error = ELPA_ERROR
583
584
585
586
587
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif

588
589
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
590
591
592
593
          print *,"Could not querry summary timings"
          stop
        endif

594
595
596
        call self%set("time_evp_fwd", time_evp_fwd)
        call self%set("time_evp_solve", time_evp_solve)
        call self%set("time_evp_back", time_evp_back)
597
598
      else

599
600
601
        call self%set("time_evp_fwd", -1.0_rk8)
        call self%set("time_evp_solve", -1.0_rk8)
        call self%set("time_evp_back", -1.0_rk8)
602
      endif
603
604
605
    end subroutine


606
    subroutine elpa_solve_complex_single(self, a, ev, q, error)
607
608
      use elpa2_impl
      use elpa1_impl
609
610
611
      use elpa_utilities, only : error_unit

      use iso_c_binding
612
      use precision
613
      class(elpa_impl_t)            :: self
614
615
616
617
618
619
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck4)             :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
      complex(kind=ck4)             :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
#endif
      real(kind=rk4)                :: ev(self%na)
620

621
      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
622
623
      integer, optional             :: error
      integer(kind=c_int)           :: error_actual
624
      logical                       :: success_l, summary_timings
625

626
      logical                       :: useGPU
627
      integer(kind=c_int) :: kernel
628
#ifdef WANT_SINGLE_PRECISION_COMPLEX
629

630
631
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
632
633
634
635
636
637
638
639
640
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif

641
642
      if (self%get("gpu",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
643
644
645
646
647
648
649
650
651
          print *,"Could not querry gpu"
          stop
        endif

        useGPU = .true.
      else
        useGPU = .false.
      endif

652
653
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
        success_l = elpa_solve_evp_complex_1stage_single_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
654
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
655
656
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
657
                                                          time_evp_solve, time_evp_back, summary_timings)
658

659
660
661
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
        kernel = self%get_complex_kernel()
        success_l = elpa_solve_evp_complex_2stage_single_impl(self, self%na, self%nev, a, self%local_nrows, ev, q,  &
662
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
663
664
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"),  time_evp_fwd,     &
665
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
666
                                                          kernel)
667
668
669
670
      else
        print *,"unknown solver"
        stop
      endif
671

672
      if (present(error)) then
673
        if (success_l) then
674
          error = ELPA_OK
675
        else
676
          error = ELPA_ERROR
677
678
679
680
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
681

682
683
      if (self%get("summary_timings",error_actual) .eq. 1) then
        if (error_actual .ne. ELPA_OK) then
684
685
686
687
          print *,"Could not querry summary timings"
          stop
        endif

688
689
690
        call self%set("time_evp_fwd", time_evp_fwd)
        call self%set("time_evp_solve", time_evp_solve)
        call self%set("time_evp_back", time_evp_back)
691
692
      else

693
694
695
        call self%set("time_evp_fwd", -1.0_rk8)
        call self%set("time_evp_solve", -1.0_rk8)
        call self%set("time_evp_back", -1.0_rk8)
696
697
      endif

698
#else
699
      print *,"This installation of the ELPA library has not been build with single-precision support"
700
      error = ELPA_ERROR
701
702
703
#endif
    end subroutine

704

705
    subroutine elpa_multiply_at_b_double (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
706
                                          c, ldc, ldcCols, error)
707
      use iso_c_binding
708
      use elpa1_auxiliary_impl
709
      use precision
710
      class(elpa_impl_t)              :: self
711
      character*1                     :: uplo_a, uplo_c
712
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
713
714
715
#ifdef USE_ASSUMED_SIZE
      real(kind=rk8)                  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
716
      real(kind=rk8)                  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
717
#endif
718
      integer, optional               :: error
719
720
      logical                         :: success_l

721
      success_l = elpa_mult_at_b_real_double_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
722
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
723
      if (present(error)) then
724
        if (success_l) then
725
          error = ELPA_OK
726
        else
727
          error = ELPA_ERROR
728
729
730
731
732
733
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in multiply_a_b() and you did not check for errors!"
      endif
    end subroutine

734

735
    subroutine elpa_multiply_at_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
736
                                          c, ldc, ldcCols, error)
737
      use iso_c_binding
738
      use elpa1_auxiliary_impl
739
      use precision
740
      class(elpa_impl_t)              :: self
741
      character*1                     :: uplo_a, uplo_c
742
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
743
744
745
#ifdef USE_ASSUMED_SIZE
      real(kind=rk4)                  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
746
      real(kind=rk4)                  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
747
#endif
748
      integer, optional               :: error
749
750
      logical                         :: success_l
#ifdef WANT_SINGLE_PRECISION_REAL
751
      success_l = elpa_mult_at_b_real_single_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
752
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
753
      if (present(error)) then
754
        if (success_l) then
755
          error = ELPA_OK
756
        else
757
          error = ELPA_ERROR
758
759
760
761
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in multiply_a_b() and you did not check for errors!"
      endif
762
763
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
764
      error = ELPA_ERROR
765
766
767
#endif
    end subroutine

768

769
    subroutine elpa_multiply_ah_b_double (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
770
                                          c, ldc, ldcCols, error)
771
      use iso_c_binding
772
      use elpa1_auxiliary_impl
773
      use precision
774
      class(elpa_impl_t)              :: self
775
      character*1                     :: uplo_a, uplo_c
776
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
777
778
779
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck8)               :: a(lda,*), b(ldb,*), c(ldc,*)
#else
780
      complex(kind=ck8)               :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
781
#endif
782
      integer, optional               :: error
783
784
      logical                         :: success_l

785
      success_l = elpa_mult_ah_b_complex_double_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
786
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
787
      if (present(error)) then
788
        if (success_l) then
789
          error = ELPA_OK
790
        else
791
          error = ELPA_ERROR
792
793
794
795
796
797
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in multiply_a_b() and you did not check for errors!"
      endif
    end subroutine

798

799
    subroutine elpa_multiply_ah_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
800
                                          c, ldc, ldcCols, error)
801
      use iso_c_binding
802
      use elpa1_auxiliary_impl
803
      use precision
804
      class(elpa_impl_t)              :: self
805
      character*1                     :: uplo_a, uplo_c
806
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols,