elpa_impl.F90 35.7 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
!> \brief Fortran module which provides the implementation of the API
51
module elpa_impl
52
  use elpa_abstract_impl
53
  use, intrinsic :: iso_c_binding
54
  implicit none
55

56
57
  private
  public :: elpa_impl_allocate
58

59
!> \brief Definition of the extended elpa_impl_t type
60
  type, extends(elpa_abstract_impl_t) :: elpa_impl_t
Andreas Marek's avatar
Andreas Marek committed
61
   private
62
   type(c_ptr)         :: index = C_NULL_PTR
63

64
   !> \brief methods available with the elpa_impl_t type
65
   contains
66
     !> \brief the puplic methods
67
     ! con-/destructor
68
69
     procedure, public :: setup => elpa_setup                   !< a setup method: implemented in elpa_setup
     procedure, public :: destroy => elpa_destroy               !< a destroy method: implemented in elpa_destroy
70

71
     ! KV store
72
73
74
75
76
77
     procedure, public :: get => elpa_get_integer               !< a get method for integer key/values: implemented in elpa_get_integer
     procedure, public :: get_double => elpa_get_double         !< a get method for double key/values: implemented in elpa_get_double
     procedure, public :: is_set => elpa_is_set                 !< a method to check whether a key/value pair has been set : implemented
                                                                !< in elpa_is_set
     procedure, public :: can_set => elpa_can_set               !< a method to check whether a key/value pair can be set : implemented
                                                                !< in elpa_can_set
78

79
80
81
82
83
84

     ! timer
     procedure, public :: get_time => elpa_get_time
     procedure, public :: print_times => elpa_print_times


85
     !> \brief the private methods
86

87
     procedure, private :: elpa_set_integer                     !< private methods to implement the setting of an integer/double key/value pair
88
     procedure, private :: elpa_set_double
89

90
91
     procedure, private :: elpa_solve_real_double               !< private methods to implement the solve step for real/complex
                                                                !< double/single matrices
92
     procedure, private :: elpa_solve_real_single
93
     procedure, private :: elpa_solve_complex_double
94
     procedure, private :: elpa_solve_complex_single
95

96
97
98
     procedure, private :: elpa_multiply_at_b_double            !< private methods to implement a "hermitian" multiplication of matrices a and b
     procedure, private :: elpa_multiply_at_b_single            !< for real valued matrices:   a**T * b
     procedure, private :: elpa_multiply_ah_b_double            !< for complex valued matrices:   a**H * b
99
     procedure, private :: elpa_multiply_ah_b_single
100

101
102
     procedure, private :: elpa_cholesky_double_real            !< private methods to implement the cholesky factorisation of
                                                                !< real/complex double/single matrices
103
104
105
     procedure, private :: elpa_cholesky_single_real
     procedure, private :: elpa_cholesky_double_complex
     procedure, private :: elpa_cholesky_single_complex
106

107
108
     procedure, private :: elpa_invert_trm_double_real          !< private methods to implement the inversion of a triangular
                                                                !< real/complex double/single matrix
109
110
111
     procedure, private :: elpa_invert_trm_single_real
     procedure, private :: elpa_invert_trm_double_complex
     procedure, private :: elpa_invert_trm_single_complex
112

113
114
     procedure, private :: elpa_solve_tridi_double_real         !< private methods to implement the solve step for a real valued
     procedure, private :: elpa_solve_tridi_single_real         !< double/single tridiagonal matrix
115

116
     procedure, private :: associate_int => elpa_associate_int  !< private method to set some pointers
117

118
  end type elpa_impl_t
119

120
  !> \brief the implementation of the private methods
121
  contains
122
123
124
125
    !> \brief function to allocate an ELPA object
    !> Parameters
    !> \param   error      integer, optional to get an error code
    !> \result  obj        class(elpa_impl_t) allocated ELPA object
126
    function elpa_impl_allocate(error) result(obj)
Andreas Marek's avatar
Andreas Marek committed
127
128
      use precision
      use elpa_utilities, only : error_unit
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
129
      use elpa_generated_fortran_interfaces
Andreas Marek's avatar
Andreas Marek committed
130

131
132
133
134
      type(elpa_impl_t), pointer   :: obj
      integer, optional            :: error

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

Andreas Marek's avatar
Andreas Marek committed
136
      ! check whether init has ever been called
137
      if ( elpa_initialized() .ne. ELPA_OK) then
138
        write(error_unit, *) "elpa_allocate(): you must call elpa_init() once before creating instances of ELPA"
139
140
        if(present(error)) then
          error = ELPA_ERROR
141
        endif
Andreas Marek's avatar
Andreas Marek committed
142
143
        return
      endif
Andreas Marek's avatar
Andreas Marek committed
144

145
      obj%index = elpa_index_instance_c()
146
147

      ! Associate some important integer pointers for convenience
148
149
150
151
152
153
154
155
      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
156
157
      endif
    end function
Andreas Marek's avatar
Andreas Marek committed
158

159
160
161
162
    !> \brief function to setup an ELPA object and to store the MPI communicators internally
    !> Parameters
    !> \param   self       class(elpa_impl_t), the allocated ELPA object
    !> \result  error      integer, the error code
163
    function elpa_setup(self) result(error)
164
      use elpa1_impl, only : elpa_get_communicators_impl
165
      class(elpa_impl_t), intent(inout) :: self
166
      integer :: error, error2
167
      integer :: mpi_comm_rows, mpi_comm_cols, mpierr
168

169
      error = ELPA_ERROR
170

171
172
173
      if (self%is_set("mpi_comm_parent") == 1 .and. &
          self%is_set("process_row") == 1 .and. &
          self%is_set("process_col") == 1) then
174

175
176
177
178
179
180
        mpierr = elpa_get_communicators_impl(&
                        self%get("mpi_comm_parent"), &
                        self%get("process_row"), &
                        self%get("process_col"), &
                        mpi_comm_rows, &
                        mpi_comm_cols)
181

182
183
184
        call self%set("mpi_comm_rows", mpi_comm_rows)
        call self%set("mpi_comm_cols", mpi_comm_cols)

185
        error = ELPA_OK
186
      endif
187

188
189
      if (self%is_set("mpi_comm_rows") == 1 .and. self%is_set("mpi_comm_cols") == 1) then
        error = ELPA_OK
190
      endif
191

192
193
194
195
      if (self%get("timings") == 1) then
        call self%timer%enable()
      endif

196
    end function
197

198
199
200
201
202
203
    !> \brief subroutine to set an integer key/value pair
    !> Parameters
    !> \param   self       class(elpa_impl_t) the allocated ELPA object
    !> \param   name       string, the key
    !> \param   value      integer, the value to be set
    !> \result  error      integer, the error code
204
    subroutine elpa_set_integer(self, name, value, error)
205
      use iso_c_binding
206
207
      use elpa_generated_fortran_interfaces
      use elpa_utilities, only : error_unit
208
      class(elpa_impl_t)              :: self
209
210
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
211
212
      integer, optional               :: error
      integer                         :: actual_error
213

214
      actual_error = elpa_index_set_int_value_c(self%index, name // c_null_char, value, 0)
215

216
217
      if (present(error)) then
        error = actual_error
218

219
      else if (actual_error /= ELPA_OK) then
220
221
        write(error_unit,'(a,i0,a)') "ELPA: Error setting option '" // name // "' to value ", value, &
                " (got: " // elpa_strerr(actual_error) // ") and you did not check for errors!"
222
      end if
223
224
    end subroutine

225
226
227
228
229
230
    !> \brief function to get an integer key/value pair
    !> Parameters
    !> \param   self       class(elpa_impl_t) the allocated ELPA object
    !> \param   name       string, the key
    !> \param   error      integer, optional, to store an error code
    !> \result  value      integer, the value of the key/vaue pair
231
    function elpa_get_integer(self, name, error) result(value)
232
      use iso_c_binding
233
      use elpa_generated_fortran_interfaces
234
      use elpa_utilities, only : error_unit
235
      class(elpa_impl_t)             :: self
236
237
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
238
      integer, intent(out), optional :: error
239
      integer                        :: actual_error
240

241
242
243
244
245
246
247
      value = elpa_index_get_int_value_c(self%index, name // c_null_char, actual_error)
      if (present(error)) then
        error = actual_error
      else if (actual_error /= ELPA_OK) then
        write(error_unit,'(a)') "ELPA: Error getting option '" // name // "'" // &
                " (got: " // elpa_strerr(actual_error) // ") and you did not check for errors!"
      end if
248
    end function
Andreas Marek's avatar
Andreas Marek committed
249

250
251
252
253
254
    !> \brief function to check whether a key/value pair is set
    !> Parameters
    !> \param   self       class(elpa_impl_t) the allocated ELPA object
    !> \param   name       string, the key
    !> \result  state      integer, the state of the key/value pair
255
    function elpa_is_set(self, name) result(state)
256
257
      use iso_c_binding
      use elpa_generated_fortran_interfaces
258
      class(elpa_impl_t)       :: self
259
      character(*), intent(in) :: name
260
      integer                  :: state
261

262
      state = elpa_index_value_is_set_c(self%index, name // c_null_char)
263
264
    end function

265
266
267
268
269
270
    !> \brief function to check whether a key/value pair can be set
    !> Parameters
    !> \param   self       class(elpa_impl_t) the allocated ELPA object
    !> \param   name       string, the key
    !> \param   value      integer, value
    !> \result  error      integer, error code
271
272
273
274
275
276
277
278
279
280
281
282
283
    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)
284
285
286
      use elpa_generated_fortran_interfaces
      class(elpa_impl_t), intent(in) :: self
      character(kind=c_char, len=*), intent(in) :: option_name
287
288
289
290
      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
291

292
293
294
295
296
297
298
299
      nullify(string)

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

302
303
304
305
      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
306

307
308
309
310
      if (present(error)) then
        error = actual_error
      endif
    end function
311

312
313

    subroutine elpa_set_double(self, name, value, error)
Andreas Marek's avatar
Andreas Marek committed
314
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
315
      use elpa_generated_fortran_interfaces
316
      use elpa_utilities, only : error_unit
317
      class(elpa_impl_t)              :: self
318
      character(*), intent(in)        :: name
319
      real(kind=c_double), intent(in) :: value
320
321
      integer, optional               :: error
      integer                         :: actual_error
Andreas Marek's avatar
Andreas Marek committed
322

323
      actual_error = elpa_index_set_double_value_c(self%index, name // c_null_char, value, 0)
Andreas Marek's avatar
Andreas Marek committed
324

325
326
327
      if (present(error)) then
        error = actual_error
      else if (actual_error /= ELPA_OK) then
328
329
        write(error_unit,'(a,es12.5,a)') "ELPA: Error setting option '" // name // "' to value ", value, &
                " (got: " // elpa_strerr(actual_error) // ") and you did not check for errors!"
330
331
      end if
    end subroutine
Andreas Marek's avatar
Andreas Marek committed
332
333


334
    function elpa_get_double(self, name, error) result(value)
Andreas Marek's avatar
Andreas Marek committed
335
      use iso_c_binding
Lorenz Huedepohl's avatar
Lorenz Huedepohl committed
336
      use elpa_generated_fortran_interfaces
337
      use elpa_utilities, only : error_unit
338
      class(elpa_impl_t)             :: self
339
      character(*), intent(in)       :: name
340
      real(kind=c_double)            :: value
341
      integer, intent(out), optional :: error
342
      integer                        :: actual_error
343

344
345
346
347
348
349
350
      value = elpa_index_get_double_value_c(self%index, name // c_null_char, actual_error)
      if (present(error)) then
        error = actual_error
      else if (actual_error /= ELPA_OK) then
        write(error_unit,'(a)') "ELPA: Error getting option '" // name // "'" // &
                " (got: " // elpa_strerr(actual_error) // ") and you did not check for errors!"
      end if
351
    end function
Andreas Marek's avatar
Andreas Marek committed
352
353


354
    function elpa_associate_int(self, name) result(value)
Andreas Marek's avatar
Andreas Marek committed
355
      use iso_c_binding
356
      use elpa_generated_fortran_interfaces
357
358
      use elpa_utilities, only : error_unit
      class(elpa_impl_t)             :: self
359
360
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
Andreas Marek's avatar
Andreas Marek committed
361

362
363
      type(c_ptr)                    :: value_p

364
      value_p = elpa_index_get_int_loc_c(self%index, name // c_null_char)
365
366
367
      if (.not. c_associated(value_p)) then
        write(error_unit, '(a,a,a)') "ELPA: Warning, received NULL pointer for entry '", name, "'"
      endif
368
369
      call c_f_pointer(value_p, value)
    end function
Andreas Marek's avatar
Andreas Marek committed
370

371

372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
    function elpa_get_time(self, name1, name2, name3, name4, name5, name6) result(s)
      class(elpa_impl_t), intent(in) :: self
      ! this is clunky, but what can you do..
      character(len=*), intent(in), optional :: name1, name2, name3, name4, name5, name6
      real(kind=c_double) :: s

      s = self%timer%get(name1, name2, name3, name4, name5, name6)
    end function


    subroutine elpa_print_times(self)
      class(elpa_impl_t), intent(in) :: self
      call self%timer%print()
    end subroutine


388
    subroutine elpa_solve_real_double(self, a, ev, q, error)
389
390
      use elpa2_impl
      use elpa1_impl
391
      use elpa_utilities, only : error_unit
392
      use precision
Andreas Marek's avatar
Andreas Marek committed
393
      use iso_c_binding
394
      class(elpa_impl_t)  :: self
Andreas Marek's avatar
Andreas Marek committed
395

396
397
398
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
399
      real(kind=c_double) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
400
#endif
401
      real(kind=c_double) :: ev(self%na)
402

403
404
      integer, optional   :: error
      integer(kind=c_int) :: error_actual
405
      logical             :: success_l
406

407

408
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
409
        success_l = elpa_solve_evp_real_1stage_double_impl(self, a, ev, q)
410

411
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
412
        success_l = elpa_solve_evp_real_2stage_double_impl(self, a, ev, q)
413
414
415
416
      else
        print *,"unknown solver"
        stop
      endif
417

418
      if (present(error)) then
419
        if (success_l) then
420
          error = ELPA_OK
421
        else
422
          error = ELPA_ERROR
423
424
425
426
427
428
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
    end subroutine

429

430
    subroutine elpa_solve_real_single(self, a, ev, q, error)
431
432
      use elpa2_impl
      use elpa1_impl
433
      use elpa_utilities, only : error_unit
434
      use precision
435
      use iso_c_binding
436
      class(elpa_impl_t)  :: self
437
438
439
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)  :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
440
      real(kind=c_float)  :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
441
#endif
442
      real(kind=c_float)  :: ev(self%na)
443

444
445
      integer, optional   :: error
      integer(kind=c_int) :: error_actual
446
      logical             :: success_l
447

448
#ifdef WANT_SINGLE_PRECISION_REAL
449

450
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
451
        success_l = elpa_solve_evp_real_1stage_single_impl(self, a, ev, q)
452

453
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
454
        success_l = elpa_solve_evp_real_2stage_single_impl(self, a, ev, q)
455
456
457
458
      else
        print *,"unknown solver"
        stop
      endif
459

460
      if (present(error)) then
461
        if (success_l) then
462
          error = ELPA_OK
463
        else
464
          error = ELPA_ERROR
465
466
467
468
469
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
#else
470
      print *,"This installation of the ELPA library has not been build with single-precision support"
471
      error = ELPA_ERROR
472
473
474
#endif
    end subroutine

475

476
    subroutine elpa_solve_complex_double(self, a, ev, q, error)
477
478
      use elpa2_impl
      use elpa1_impl
479
      use elpa_utilities, only : error_unit
480
      use precision
481
      use iso_c_binding
482
      class(elpa_impl_t)             :: self
483

484
485
486
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
487
      complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
488
#endif
489
      real(kind=c_double)            :: ev(self%na)
490

491
492
      integer, optional              :: error
      integer(kind=c_int)            :: error_actual
493
      logical                        :: success_l
494

495
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
496
        success_l = elpa_solve_evp_complex_1stage_double_impl(self, a, ev, q)
497

498
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
499
        success_l = elpa_solve_evp_complex_2stage_double_impl(self,  a, ev, q)
500
501
502
503
      else
        print *,"unknown solver"
        stop
      endif
504

505
      if (present(error)) then
506
        if (success_l) then
507
          error = ELPA_OK
508
        else
509
          error = ELPA_ERROR
510
511
512
513
514
515
516
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
    end subroutine


517
    subroutine elpa_solve_complex_single(self, a, ev, q, error)
518
519
      use elpa2_impl
      use elpa1_impl
520
521
522
      use elpa_utilities, only : error_unit

      use iso_c_binding
523
      use precision
524
      class(elpa_impl_t)            :: self
525
526
527
528
529
530
#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)
531

532
533
      integer, optional             :: error
      integer(kind=c_int)           :: error_actual
534
      logical                       :: success_l
535
536

#ifdef WANT_SINGLE_PRECISION_COMPLEX
537

538
      if (self%get("solver") .eq. ELPA_SOLVER_1STAGE) then
539
        success_l = elpa_solve_evp_complex_1stage_single_impl(self, a, ev, q)
540

541
      else if (self%get("solver") .eq. ELPA_SOLVER_2STAGE) then
542
        success_l = elpa_solve_evp_complex_2stage_single_impl(self,  a, ev, q)
543
544
545
546
      else
        print *,"unknown solver"
        stop
      endif
547

548
      if (present(error)) then
549
        if (success_l) then
550
          error = ELPA_OK
551
        else
552
          error = ELPA_ERROR
553
554
555
556
557
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
      endif
#else
558
      print *,"This installation of the ELPA library has not been build with single-precision support"
559
      error = ELPA_ERROR
560
561
562
#endif
    end subroutine

563

564
    subroutine elpa_multiply_at_b_double (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
565
                                          c, ldc, ldcCols, error)
566
      use iso_c_binding
567
      use elpa1_auxiliary_impl
568
      use precision
569
      class(elpa_impl_t)              :: self
570
      character*1                     :: uplo_a, uplo_c
571
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
572
573
574
#ifdef USE_ASSUMED_SIZE
      real(kind=rk8)                  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
575
      real(kind=rk8)                  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
576
#endif
577
      integer, optional               :: error
578
579
      logical                         :: success_l

580
581
      success_l = elpa_mult_at_b_real_double_impl(self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                                  c, ldc, ldcCols)
582
      if (present(error)) then
583
        if (success_l) then
584
          error = ELPA_OK
585
        else
586
          error = ELPA_ERROR
587
588
589
590
591
592
        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

593

594
    subroutine elpa_multiply_at_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
595
                                          c, ldc, ldcCols, error)
596
      use iso_c_binding
597
      use elpa1_auxiliary_impl
598
      use precision
599
      class(elpa_impl_t)              :: self
600
      character*1                     :: uplo_a, uplo_c
601
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
602
603
604
#ifdef USE_ASSUMED_SIZE
      real(kind=rk4)                  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
605
      real(kind=rk4)                  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
606
#endif
607
      integer, optional               :: error
608
609
      logical                         :: success_l
#ifdef WANT_SINGLE_PRECISION_REAL
610
611
      success_l = elpa_mult_at_b_real_single_impl(self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                                  c, ldc, ldcCols)
612
      if (present(error)) then
613
        if (success_l) then
614
          error = ELPA_OK
615
        else
616
          error = ELPA_ERROR
617
618
619
620
        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
621
622
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
623
      error = ELPA_ERROR
624
625
626
#endif
    end subroutine

627

628
    subroutine elpa_multiply_ah_b_double (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
629
                                          c, ldc, ldcCols, error)
630
      use iso_c_binding
631
      use elpa1_auxiliary_impl
632
      use precision
633
      class(elpa_impl_t)              :: self
634
      character*1                     :: uplo_a, uplo_c
635
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
636
637
638
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck8)               :: a(lda,*), b(ldb,*), c(ldc,*)
#else
639
      complex(kind=ck8)               :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
640
#endif
641
      integer, optional               :: error
642
643
      logical                         :: success_l

644
645
      success_l = elpa_mult_ah_b_complex_double_impl(self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                                     c, ldc, ldcCols)
646
      if (present(error)) then
647
        if (success_l) then
648
          error = ELPA_OK
649
        else
650
          error = ELPA_ERROR
651
652
653
654
655
656
        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

657

658
    subroutine elpa_multiply_ah_b_single (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
659
                                          c, ldc, ldcCols, error)
660
      use iso_c_binding
661
      use elpa1_auxiliary_impl
662
      use precision
663
      class(elpa_impl_t)              :: self
664
      character*1                     :: uplo_a, uplo_c
665
      integer(kind=ik), intent(in)    :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
666
667
668
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck4)               :: a(lda,*), b(ldb,*), c(ldc,*)
#else
669
      complex(kind=ck4)               :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
670
#endif
671
      integer, optional               :: error
672
673
674
      logical                         :: success_l

#ifdef WANT_SINGLE_PRECISION_COMPLEX
675
676
      success_l = elpa_mult_ah_b_complex_single_impl(self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
                                                     c, ldc, ldcCols)
677
      if (present(error)) then
678
        if (success_l) then
679
          error = ELPA_OK
680
        else
681
          error = ELPA_ERROR
682
683
684
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in multiply_a_b() and you did not check for errors!"
685
      endif 
686
687
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
688
      error = ELPA_ERROR
689
690
691
#endif
    end subroutine

692

693
    subroutine elpa_cholesky_double_real (self, a, error)
694
      use iso_c_binding
695
      use elpa1_auxiliary_impl
696
      use precision
697
      class(elpa_impl_t)              :: self
698
699
700
#ifdef USE_ASSUMED_SIZE
      real(kind=rk8)                  :: a(self%local_nrows,*)
#else
701
      real(kind=rk8)                  :: a(self%local_nrows,self%local_ncols)
702
#endif
703
      integer, optional               :: error
704
      logical                         :: success_l
705
      integer(kind=c_int)             :: error_actual
706

707
      success_l = elpa_cholesky_real_double_impl (self, a)
708
      if (present(error)) then
709
        if (success_l) then
710
          error = ELPA_OK
711
        else
712
          error = ELPA_ERROR
713
714
715
716
717
718
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
      endif
    end subroutine

719

720
    subroutine elpa_cholesky_single_real(self, a, error)
721
      use iso_c_binding
722
      use elpa1_auxiliary_impl
723
      use precision
724
      class(elpa_impl_t)              :: self
725
726
727
#ifdef USE_ASSUMED_SIZE
      real(kind=rk4)                  :: a(self%local_nrows,*)
#else
728
      real(kind=rk4)                  :: a(self%local_nrows,self%local_ncols)
729
#endif
730
      integer, optional               :: error
731
      logical                         :: success_l
732
      integer(kind=c_int)             :: error_actual
733
734

#if WANT_SINGLE_PRECISION_REAL
735
      success_l = elpa_cholesky_real_single_impl (self, a)
736
737
#else
      print *,"This installation of the ELPA library has not been build with single-precision support"
738
      error = ELPA_ERROR
739
#endif
740
      if (present(error)) then
741
        if (success_l) then
742
          error = ELPA_OK
743
        else
744
          error = ELPA_ERROR
745
746
747
748
749
750
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
      endif
    end subroutine

751

752
    subroutine elpa_cholesky_double_complex (self, a, error)
753
      use iso_c_binding
754
      use elpa1_auxiliary_impl
755
      use precision
756
      class(elpa_impl_t)              :: self
757
758
759
760
761
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck8)               :: a(self%local_nrows,*)
#else
      complex(kind=ck8)               :: a(self%local_nrows,self%local_ncols)
#endif
762
      integer, optional               :: error
763
      logical                         :: success_l
764
      integer(kind=c_int)             :: error_actual
765

766
      success_l = elpa_cholesky_complex_double_impl (self, a)
767
      if (present(error)) then
768
        if (success_l) then
769
          error = ELPA_OK
770
        else
771
          error = ELPA_ERROR
772
773
774
775
776
777
        endif
      else if (.not. success_l) then
        write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
      endif
    end subroutine

778

779
    subroutine elpa_cholesky_single_complex (self, a, error)
780
      use iso_c_binding
781
      use elpa1_auxiliary_impl
782
      use precision
783
      class(elpa_impl_t)              :: self
784
785
786
787
788
#ifdef USE_ASSUMED_SIZE
      complex(kind=ck4)               :: a(self%local_nrows,*)
#else
      complex(kind=ck4)               :: a(self%local_nrows,self%local_ncols)
#endif
789
      integer, optional               :: error
790
      logical                         :: success_l
791
      integer(kind=c_int)             :: error_actual
792

793
#if WANT_SINGLE_PRECISION_COMPLEX
794
      success_l = elpa_cholesky_complex_single_impl (self, a)
795
796
#else