elpa_impl.F90 47.3 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 63 64 65
   contains
     procedure, public :: setup => elpa_setup

     procedure, public :: get => elpa_get_integer
66 67 68 69 70 71
     procedure, public :: get_double => elpa_get_double
     procedure, public :: is_set => elpa_is_set
     procedure, public :: print_options => elpa_print_options

     procedure, public :: get_real_kernel => elpa_get_real_kernel
     procedure, public :: get_complex_kernel => elpa_get_complex_kernel
72 73 74 75 76

     procedure, public :: destroy => elpa_destroy

     ! privates:
     procedure, private :: elpa_set_integer
77
     procedure, private :: elpa_set_double
78

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

84 85 86 87
     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
88

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

     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
98 99 100

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

     procedure, private :: associate_int => elpa_associate_int

104
  end type elpa_impl_t
105 106 107

  contains

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

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

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

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

127 128 129
      obj%index = elpa_index_instance()

      ! Associate some important integer pointers for convenience
130 131 132 133 134 135 136 137
      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
138 139
      endif
    end function
Andreas Marek's avatar
Andreas Marek committed
140

141

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

148
      success = ELPA_ERROR
149

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

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

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

        success = ELPA_OK
      endif
166

167 168 169
      if (self%is_set("mpi_comm_rows") == ELPA_OK .and. self%is_set("mpi_comm_cols") == ELPA_OK) then
        success = ELPA_OK
      endif
170 171

    end function
172 173

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

183
      actual_success = elpa_index_set_int_value(self%index, name // c_null_char, value)
184

185 186 187 188
      if (present(success)) then
        success = actual_success

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

      end if
193 194
    end subroutine

195

196 197 198 199 200 201 202 203 204 205 206 207 208
    subroutine elpa_set_integer_by_string(self, name, value, success)
      use iso_c_binding
      use elpa_generated_fortran_interfaces
      use elpa_utilities, only : error_unit
      class(elpa_impl_t)              :: self
      character(*), intent(in)        :: name
      character(*), intent(in)        :: value
      integer, optional               :: success
      integer                         :: actual_success

    end subroutine


209
    function elpa_get_integer(self, name, success) result(value)
210
      use iso_c_binding
211
      use elpa_generated_fortran_interfaces
212
      class(elpa_impl_t)             :: self
213 214 215
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
      integer, intent(out), optional :: success
216

217
      value = elpa_index_get_int_value(self%index, name // c_null_char, success)
218

219
    end function
Andreas Marek's avatar
Andreas Marek committed
220

221 222

    function elpa_is_set(self, name) result(state)
223 224
      use iso_c_binding
      use elpa_generated_fortran_interfaces
225
      class(elpa_impl_t)       :: self
226
      character(*), intent(in) :: name
227
      integer                  :: state
228

229
      state = elpa_index_value_is_set(self%index, name // c_null_char)
230 231
    end function

232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255

    subroutine elpa_print_options(self, option_name, unit)
      use elpa_utilities, only : output_unit
      use elpa_generated_fortran_interfaces
      class(elpa_impl_t), intent(in) :: self
      character(kind=c_char, len=*), intent(in) :: option_name
      integer, intent(in), optional :: unit
      integer :: unit_actual, i, val

      if (present(unit)) then
        unit_actual = unit
      else
        unit_actual = output_unit
      endif

      do i = 0, elpa_index_cardinality_c(self%index, option_name // C_NULL_CHAR)
        val = elpa_index_enumerate_c(self%index, option_name // C_NULL_CHAR, i)
        if (elpa_index_int_is_valid_c(self%index, option_name // C_NULL_CHAR, val) == ELPA_OK) then
          write(unit_actual, '(a)') " " // elpa_c_string(elpa_index_int_value_to_string_helper_c(option_name // C_NULL_CHAR, val))
        end if
      end do
    end subroutine


256
    subroutine elpa_set_double(self, name, value, success)
Andreas Marek's avatar
Andreas Marek committed
257
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
258
      use elpa_generated_fortran_interfaces
259
      use elpa_utilities, only : error_unit
260
      class(elpa_impl_t)              :: self
261
      character(*), intent(in)        :: name
262
      real(kind=c_double), intent(in) :: value
263 264
      integer, optional               :: success
      integer                         :: actual_success
Andreas Marek's avatar
Andreas Marek committed
265

266
      actual_success = elpa_index_set_double_value(self%index, name // c_null_char, value)
Andreas Marek's avatar
Andreas Marek committed
267

268 269 270 271
      if (present(success)) then
        success = actual_success

      else if (actual_success /= ELPA_OK) then
272
        write(error_unit,'(a,a,es12.5,a)') "ELPA: Error setting option '", name, "' to value ", value, &
273
                " and you did not check for errors!"
274 275 276

      end if
    end subroutine
Andreas Marek's avatar
Andreas Marek committed
277 278


279
    function elpa_get_double(self, name, success) result(value)
Andreas Marek's avatar
Andreas Marek committed
280
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
281
      use elpa_generated_fortran_interfaces
282
      class(elpa_impl_t)             :: self
283
      character(*), intent(in)       :: name
284
      real(kind=c_double)            :: value
285 286 287
      integer, intent(out), optional :: success
      type(c_ptr) :: c_success_ptr

288
      value = elpa_index_get_double_value(self%index, name // c_null_char, success)
289 290

    end function
Andreas Marek's avatar
Andreas Marek committed
291 292


293
    function elpa_associate_int(self, name) result(value)
Andreas Marek's avatar
Andreas Marek committed
294
      use iso_c_binding
295
      use elpa_generated_fortran_interfaces
296 297
      use elpa_utilities, only : error_unit
      class(elpa_impl_t)             :: self
298 299
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
Andreas Marek's avatar
Andreas Marek committed
300

301 302 303
      type(c_ptr)                    :: value_p

      value_p = elpa_index_get_int_loc(self%index, name // c_null_char)
304 305 306
      if (.not. c_associated(value_p)) then
        write(error_unit, '(a,a,a)') "ELPA: Warning, received NULL pointer for entry '", name, "'"
      endif
307 308
      call c_f_pointer(value_p, value)
    end function
Andreas Marek's avatar
Andreas Marek committed
309

310 311

    subroutine elpa_solve_real_double(self, a, ev, q, success)
312 313
      use elpa2_impl
      use elpa1_impl
314
      use elpa_utilities, only : error_unit
315
      use precision
Andreas Marek's avatar
Andreas Marek committed
316
      use iso_c_binding
317
      class(elpa_impl_t)  :: self
Andreas Marek's avatar
Andreas Marek committed
318

319 320 321
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
322
      real(kind=c_double) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
323
#endif
324
      real(kind=c_double) :: ev(self%na)
325 326

      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
327 328
      integer, optional   :: success
      integer(kind=c_int) :: success_internal
329
      logical             :: success_l, summary_timings
330

331
      logical             :: useGPU, useQR
332
      integer(kind=c_int) :: kernel
333

334 335 336 337 338 339 340 341 342 343 344 345
      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

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


346 347 348 349 350 351 352 353 354 355 356
      if (self%get("gpu",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry gpu"
          stop
        endif

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

357
      useQR = self%get("qr") == 1
358

359 360
      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,  &
361 362 363 364
                                                           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)
365

366 367 368
      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,  &
369 370 371 372 373
                                                           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)
374 375 376 377
      else
        print *,"unknown solver"
        stop
      endif
378 379 380 381 382 383 384 385 386 387

      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        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
388

389 390 391 392 393 394
      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

395 396 397
        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)
398 399
      else

400 401 402
        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)
403
      endif
404 405
    end subroutine

406

407
    subroutine elpa_solve_real_single(self, a, ev, q, success)
408 409
      use elpa2_impl
      use elpa1_impl
410
      use elpa_utilities, only : error_unit
411
      use precision
412
      use iso_c_binding
413
      class(elpa_impl_t)  :: self
414 415 416
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)  :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
417
      real(kind=c_float)  :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
418
#endif
419
      real(kind=c_float)  :: ev(self%na)
420 421

      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
422 423
      integer, optional   :: success
      integer(kind=c_int) :: success_internal
424
      logical             :: success_l, summary_timings
425

426
      logical             :: useGPU, useQR
427
      integer(kind=c_int) :: kernel
428

429
#ifdef WANT_SINGLE_PRECISION_REAL
430 431 432 433 434 435 436 437 438 439
      if (self%get("timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif
440 441 442 443 444 445 446 447 448 449 450 451

      if (self%get("gpu",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry gpu"
          stop
        endif

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

452
      useQR = self%get("qr") == 1
453

454 455
      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,  &
456
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
457 458
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
459
                                                          time_evp_solve, time_evp_back, summary_timings)
460

461 462 463
      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,  &
464
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
465 466
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), time_evp_fwd,     &
467
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
468
                                                          kernel, useQR)
469 470 471 472
      else
        print *,"unknown solver"
        stop
      endif
473 474 475 476 477 478 479 480 481 482

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


      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

491 492 493
        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)
494 495
      else

496 497 498
        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)
499
      endif
500
#else
501
      print *,"This installation of the ELPA library has not been build with single-precision support"
502 503 504 505
      success = ELPA_ERROR
#endif
    end subroutine

506 507

    subroutine elpa_solve_complex_double(self, a, ev, q, success)
508 509
      use elpa2_impl
      use elpa1_impl
510
      use elpa_utilities, only : error_unit
511
      use precision
512
      use iso_c_binding
513
      class(elpa_impl_t)             :: self
514

515 516 517
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
518
      complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
519
#endif
520
      real(kind=c_double)            :: ev(self%na)
521

522 523
      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back

524 525
      integer, optional              :: success
      integer(kind=c_int)            :: success_internal
526
      logical                        :: success_l, summary_timings
527

528
      logical                        :: useGPU
529
      integer(kind=c_int) :: kernel
530 531 532 533 534 535 536 537 538 539
      if (self%get("timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

        summary_timings = .true.
      else
        summary_timings = .false.
      endif
540 541 542 543 544 545 546 547 548 549 550 551

      if (self%get("gpu",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry gpu"
          stop
        endif

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

552 553
      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,  &
554
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
555 556
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
557
                                                          time_evp_solve, time_evp_back, summary_timings)
558

559 560 561
      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,  &
562
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
563 564
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), time_evp_fwd,     &
565
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
566
                                                          kernel)
567 568 569 570
      else
        print *,"unknown solver"
        stop
      endif
571 572 573 574 575 576 577 578 579 580 581

      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif

582 583 584 585 586 587
      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

588 589 590
        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)
591 592
      else

593 594 595
        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)
596
      endif
597 598 599
    end subroutine


600
    subroutine elpa_solve_complex_single(self, a, ev, q, success)
601 602
      use elpa2_impl
      use elpa1_impl
603 604 605
      use elpa_utilities, only : error_unit

      use iso_c_binding
606
      use precision
607
      class(elpa_impl_t)            :: self
608 609 610 611 612 613
#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)
614

615
      real(kind=c_double) :: time_evp_fwd, time_evp_solve, time_evp_back
616 617
      integer, optional             :: success
      integer(kind=c_int)           :: success_internal
618
      logical                       :: success_l, summary_timings
619

620
      logical                       :: useGPU
621
      integer(kind=c_int) :: kernel
622
#ifdef WANT_SINGLE_PRECISION_COMPLEX
623 624 625 626 627 628 629 630 631 632 633 634

      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

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

635 636 637 638 639 640 641 642 643 644 645
      if (self%get("gpu",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry gpu"
          stop
        endif

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

646 647
      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,  &
648
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
649 650
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"), useGPU, time_evp_fwd,     &
651
                                                          time_evp_solve, time_evp_back, summary_timings)
652

653 654 655
      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,  &
656
                                                          self%local_nrows,  self%nblk, self%local_ncols, &
657 658
                                                          self%get("mpi_comm_rows"), self%get("mpi_comm_cols"),         &
                                                          self%get("mpi_comm_parent"),  time_evp_fwd,     &
659
                                                          time_evp_solve, time_evp_back, summary_timings, useGPU, &
660
                                                          kernel)
661 662 663 664
      else
        print *,"unknown solver"
        stop
      endif
665 666 667 668 669 670 671 672 673 674

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

      if (self%get("summary_timings",success_internal) .eq. 1) then
        if (success_internal .ne. ELPA_OK) then
          print *,"Could not querry summary timings"
          stop
        endif

682 683 684
        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)
685 686
      else

687 688 689
        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)
690 691
      endif

692
#else
693
      print *,"This installation of the ELPA library has not been build with single-precision support"
694 695 696 697
      success = ELPA_ERROR
#endif
    end subroutine

698

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

715
      success_l = elpa_mult_at_b_real_double_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
716
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
717 718 719 720 721 722 723 724 725 726 727
      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        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

728

729 730 731
    subroutine elpa_multiply_at_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                          c, ldc, ldcCols, success)
      use iso_c_binding
732
      use elpa1_auxiliary_impl
733
      use precision
734
      class(elpa_impl_t)              :: self
735
      character*1                     :: uplo_a, uplo_c
736
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
737 738 739
#ifdef USE_ASSUMED_SIZE
      real(kind=rk4)                  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
740
      real(kind=rk4)                  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
741
#endif
742 743 744
      integer, optional               :: success
      logical                         :: success_l
#ifdef WANT_SINGLE_PRECISION_REAL
745
      success_l = elpa_mult_at_b_real_single_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
746
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
747 748 749 750 751 752 753 754 755
      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        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
756 757 758
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
      success = ELPA_ERROR
759 760 761
#endif
    end subroutine

762

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

779
      success_l = elpa_mult_ah_b_complex_double_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
780
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
781 782 783 784 785 786 787 788 789 790 791
      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        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

792

793 794 795
    subroutine elpa_multiply_ah_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                          c, ldc, ldcCols, success)
      use iso_c_binding
796
      use elpa1_auxiliary_impl
797
      use precision
798
      class(elpa_impl_t)              :: self
799
      character*1                     :: uplo_a, uplo_c
800
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
801 802 803
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck4)               :: a(lda,*), b(ldb,*), c(ldc,*)
#else
804
      complex(kind=ck4)               :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
805
#endif
806 807 808 809
      integer, optional               :: success
      logical                         :: success_l

#ifdef WANT_SINGLE_PRECISION_COMPLEX
810
      success_l = elpa_mult_ah_b_complex_single_impl(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, self%nblk, &
811
                              self%get("mpi_comm_rows"), self%get("mpi_comm_cols"), c, ldc, ldcCols)
812 813 814 815 816 817 818 819 820
      if (present(success)) then
        if (success_l) then
          success = ELPA_OK
        else
          success = ELPA_ERROR
        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
821 822 823
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
      success = ELPA_ERROR
824 825 826
#endif
    end subroutine

827

828
    subroutine elpa_cholesky_double_real (self, a, success)
829
      use iso_c_binding
830
      use elpa1_auxiliary_impl
831
      use precision
832
      class(elpa_impl_t)              :: self
833 834 835
#ifdef USE_ASSUMED_SIZE
      real(kind=rk8)                  :: a(self%local_nrows,*)
#else
836
      real(kind=rk8)                  :: a(self%local_nrows,self%local_ncols)
837
#endif
838 839 840 841 842 843 844 845 846 847 848 849 850 851