elpa_api.F90 44.1 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
!
!    Copyright 2017, L. Hüdepohl and A. Marek, MPCDF
!
!    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.
!
#include "config-f90.h"
49
50
51
!> \brief Fortran module which provides the definition of the ELPA API. Do not use directly! Use the module "elpa"


52
53
54
55
56
module elpa_api
  use elpa_constants
  use, intrinsic :: iso_c_binding
  implicit none

57
58
59
60
  integer, private, parameter :: earliest_api_version = EARLIEST_API_VERSION !< Definition of the earliest API version supported
                                                                             !< with the current release
  integer, private, parameter :: current_api_version  = CURRENT_API_VERSION  !< Definition of the current API version

61
62
63
64
65
66
67
  logical, private :: initDone = .false.

  public :: elpa_t, &
      c_int, &
      c_double, c_double_complex, &
      c_float, c_float_complex

68
  !> \brief Abstract definition of the elpa_t type
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115

  !>
  !> A typical usage of ELPA might look like this:
  !>
  !> Fortran synopsis
  !>
  !> \code{.f90}
  !>  use elpa
  !>  class(elpa_t), pointer :: elpa
  !>
  !>  if (elpa_init(20170403) /= ELPA_OK) then
  !>     print *, "ELPA API version not supported"
  !>     stop
  !>   endif
  !>   elpa => elpa_allocate()
  !>
  !>   call elpa%set("na", na, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("nev", nev, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("local_nrows", na_rows, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("local_ncols", na_cols, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("nblk", nblk, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("mpi_comm_parent", MPI_COMM_WORLD, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("process_row", my_prow, success)
  !>   assert_elpa_ok(success)
  !>   call elpa%set("process_col", my_pcol, success)
  !>   assert_elpa_ok(success)
  !>
  !>   assert(elpa%setup() .eq. ELPA_OK)
  !>
  !>   call e%set("solver", ELPA_SOLVER_2STAGE, success)
  !>   assert_elpa_ok(success)
  !> \endcode
  !>   ... set and get all other options that are desired
  !> \code{.f90}
  !>   call e%solve(a, ev, z, success)
  !>   assert_elpa_ok(success)
  !>
  !>   call elpa_deallocate(e)
  !>
  !>   call elpa_uninit()
  !> \endcode
116
117
118
  type, abstract :: elpa_t
    private

119
    !< these have to be public for proper bounds checking, sadly
120
121
122
123
124
125
126
    integer(kind=c_int), public, pointer :: na => NULL()
    integer(kind=c_int), public, pointer :: nev => NULL()
    integer(kind=c_int), public, pointer :: local_nrows => NULL()
    integer(kind=c_int), public, pointer :: local_ncols => NULL()
    integer(kind=c_int), public, pointer :: nblk => NULL()

    contains
127
      !> \brief methods available with the elpa_t type
128
      ! general
129
130
      procedure(elpa_setup_i),   deferred, public :: setup          !< export a setup method
      procedure(elpa_destroy_i), deferred, public :: destroy        !< export a destroy method
131

132
      ! key/value store
133
      generic, public :: set => &                                   !< export a method to set integer/double key/values
134
135
          elpa_set_integer, &
          elpa_set_double
136
137
      procedure(elpa_get_integer_i), deferred, public :: get        !< get method for integer key/values
      procedure(elpa_get_double_i),  deferred, public :: get_double !< get method for double key/values
138

139
140
      procedure(elpa_is_set_i),  deferred, public :: is_set         !< method to check whether key/value is set
      procedure(elpa_can_set_i), deferred, public :: can_set        !< method to check whether key/value can be set
141

142
      ! Timer
143
144
145
      procedure(elpa_get_time_i), deferred, public :: get_time
      procedure(elpa_print_times_i), deferred, public :: print_times

146
      ! Actual math routines
147
      generic, public :: solve => &                                 !< method solve for solving the eigenvalue problem
148
149
          elpa_solve_d, &                                           !< for symmetric real valued / hermitian complex valued
          elpa_solve_f, &                                           !< matrices
150
151
          elpa_solve_dc, &
          elpa_solve_fc
152

153
      generic, public :: hermitian_multiply => &                    !< method for a "hermitian" multiplication of matrices a and b
154
          elpa_hermitian_multiply_d, &                              !< for real valued matrices:   a**T * b
155
          elpa_hermitian_multiply_dc, &                             !< for complex valued matrices a**H * b
156
157
          elpa_hermitian_multiply_f, &
          elpa_hermitian_multiply_fc
158

159
      generic, public :: cholesky => &                              !< method for the cholesky factorisation of matrix a
160
161
162
163
          elpa_cholesky_d, &
          elpa_cholesky_f, &
          elpa_cholesky_dc, &
          elpa_cholesky_fc
164

Andreas Marek's avatar
Andreas Marek committed
165
      generic, public :: invert_triangular => &                     !< method to invert a upper triangular matrix a
166
167
168
169
          elpa_invert_trm_d, &
          elpa_invert_trm_f, &
          elpa_invert_trm_dc, &
          elpa_invert_trm_fc
170

171
      generic, public :: solve_tridi => &                           !< method to solve the eigenvalue problem for a tridiagonal
172
          elpa_solve_tridi_d, &                                     !< matrix
173
          elpa_solve_tridi_f
174
175


176
      !> \brief private methods of elpa_t type. NOT accessible for the user
177
178
179
180
      ! privates
      procedure(elpa_set_integer_i), deferred, private :: elpa_set_integer
      procedure(elpa_set_double_i),  deferred, private :: elpa_set_double

181
182
183
184
      procedure(elpa_solve_d_i),    deferred, private :: elpa_solve_d
      procedure(elpa_solve_f_i),    deferred, private :: elpa_solve_f
      procedure(elpa_solve_dc_i), deferred, private :: elpa_solve_dc
      procedure(elpa_solve_fc_i), deferred, private :: elpa_solve_fc
185

186
187
188
189
      procedure(elpa_hermitian_multiply_d_i),  deferred, private :: elpa_hermitian_multiply_d
      procedure(elpa_hermitian_multiply_f_i),  deferred, private :: elpa_hermitian_multiply_f
      procedure(elpa_hermitian_multiply_dc_i), deferred, private :: elpa_hermitian_multiply_dc
      procedure(elpa_hermitian_multiply_fc_i), deferred, private :: elpa_hermitian_multiply_fc
190

191
192
193
194
      procedure(elpa_cholesky_d_i),    deferred, private :: elpa_cholesky_d
      procedure(elpa_cholesky_f_i),    deferred, private :: elpa_cholesky_f
      procedure(elpa_cholesky_dc_i), deferred, private :: elpa_cholesky_dc
      procedure(elpa_cholesky_fc_i), deferred, private :: elpa_cholesky_fc
195

196
197
198
199
      procedure(elpa_invert_trm_d_i),    deferred, private :: elpa_invert_trm_d
      procedure(elpa_invert_trm_f_i),    deferred, private :: elpa_invert_trm_f
      procedure(elpa_invert_trm_dc_i), deferred, private :: elpa_invert_trm_dc
      procedure(elpa_invert_trm_fc_i), deferred, private :: elpa_invert_trm_fc
200

201
202
      procedure(elpa_solve_tridi_d_i), deferred, private :: elpa_solve_tridi_d
      procedure(elpa_solve_tridi_f_i), deferred, private :: elpa_solve_tridi_f
203
204
205
  end type elpa_t


206
207
208
209
210
  !> \brief definition of helper function to get C strlen
  !> Parameters
  !> \details
  !> \param   ptr         type(c_ptr) : pointer to string
  !> \result  size        integer(kind=c_size_t) : length of string
211
212
213
  interface
    pure function elpa_strlen_c(ptr) result(size) bind(c, name="strlen")
      use, intrinsic :: iso_c_binding
214
      implicit none
215
216
217
218
219
      type(c_ptr), intent(in), value :: ptr
      integer(kind=c_size_t) :: size
    end function
  end interface

220
221
222
223
224
  !> \brief abstract definition of setup method
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \result  error       integer : error code, which can be queried with elpa_strerr()
225
  abstract interface
226
    function elpa_setup_i(self) result(error)
227
      import elpa_t
228
      implicit none
229
      class(elpa_t), intent(inout) :: self
230
      integer :: error
231
232
233
    end function
  end interface

234
235
236
237
238
239
240
  !> \brief abstract definition of set method for integer values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \param   value       integer : the value to set for the key
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr()
241
  abstract interface
242
    subroutine elpa_set_integer_i(self, name, value, error)
243
244
      use iso_c_binding
      import elpa_t
245
      implicit none
246
247
248
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
249
      integer, optional               :: error
250
251
252
    end subroutine
  end interface

253
254
255
256
257
258
259
  !> \brief abstract definition of get method for integer values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr()
  !> \result  value       integer : the value corresponding to the key
260
  abstract interface
261
    function elpa_get_integer_i(self, name, error) result(value)
262
263
      use iso_c_binding
      import elpa_t
264
      implicit none
265
266
267
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
268
      integer, intent(out), optional :: error
269
270
271
    end function
  end interface

272
273
274
275
276
277
278
  !> \brief abstract definition of is_set method for integer values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \result  state       integer : 1 is set, 0 if not, else a negativ error code
  !>                                                    which can be queried with elpa_strerr
279
  abstract interface
280
    function elpa_is_set_i(self, name) result(state)
281
      import elpa_t
282
      implicit none
283
284
      class(elpa_t)            :: self
      character(*), intent(in) :: name
285
      integer                  :: state
286
287
288
    end function
  end interface

289
290
291
292
293
294
295
296
  !> \brief abstract definition of can_set method for integer values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \param   value       integer: the valye to associate with the key
  !> \result  state       integer : 1 is set, 0 if not, else a negativ error code
  !>                                                    which can be queried with elpa_strerr
297
  abstract interface
298
    function elpa_can_set_i(self, name, value) result(state)
299
      import elpa_t, c_int
300
      implicit none
301
302
303
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
304
      integer                         :: state
305
    end function
306
307
  end interface

308
309
310
311
312
313
314
  !> \brief abstract definition of set method for double values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !? \param   value       double: the value to associate with the key
  !> \param   error       integer. optional : error code, which can be queried with elpa_strerr
315
  abstract interface
316
    subroutine elpa_set_double_i(self, name, value, error)
317
318
      use iso_c_binding
      import elpa_t
319
      implicit none
320
321
322
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      real(kind=c_double), intent(in) :: value
323
      integer, optional               :: error
324
325
326
    end subroutine
  end interface

327
328
329
330
331
332
333
  !> \brief abstract definition of get method for double values
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \param   error       integer. optional : error code, which can be queried with elpa_strerr
  !> \result  value       double: the value associated with the key
334
  abstract interface
335
    function elpa_get_double_i(self, name, error) result(value)
336
337
      use iso_c_binding
      import elpa_t
338
      implicit none
339
340
341
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      real(kind=c_double)            :: value
342
      integer, intent(out), optional :: error
343
344
345
    end function
  end interface

346
347
348
349
350
351
  !> \brief abstract definition of associate method for integer pointers
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name        string: the name of the key
  !> \result  value       integer pointer: the value associated with the key
352
353
354
355
  abstract interface
    function elpa_associate_int_i(self, name) result(value)
      use iso_c_binding
      import elpa_t
356
      implicit none
357
358
359
360
361
362
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
    end function
  end interface

363
364
365

  ! Timer routines

366
367
368
369
370
371
  !> \brief abstract definition of get_time method to querry the timer
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
  !> \param   name1..6    string: the name of the timer entry, supports up to 6 levels
  !> \result  s           double: the time for the entry name1..6
372
373
374
  abstract interface
    function elpa_get_time_i(self, name1, name2, name3, name4, name5, name6) result(s)
      import elpa_t, c_double
375
      implicit none
376
377
378
379
380
381
382
      class(elpa_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
    end function
  end interface

383
384
385
386
  !> \brief abstract definition of print method for timer
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
387
388
389
  abstract interface
    subroutine elpa_print_times_i(self)
      import elpa_t
390
      implicit none
391
392
393
394
395
      class(elpa_t), intent(in) :: self
    end subroutine
  end interface


396
  ! Actual math routines
397

398
  !> \brief abstract definition of interface to solve double real eigenvalue problem
399
400
401
402
403
404
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix a: defines the problem to solve
  !> \param   ev          double real: on output stores the eigenvalues
  !> \param   q           double real matrix q: on output stores the eigenvalues
405
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
406
  abstract interface
407
    subroutine elpa_solve_d_i(self, a, ev, q, error)
408
409
      use iso_c_binding
      import elpa_t
410
      implicit none
411
412
413
414
415
416
417
418
      class(elpa_t)       :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
      real(kind=c_double) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
#endif
      real(kind=c_double) :: ev(self%na)

419
      integer, optional   :: error
420
421
422
    end subroutine
  end interface

423
  !> \brief abstract definition of interface to solve single real eigenvalue problem
424
425
426
427
428
429
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix a: defines the problem to solve
  !> \param   ev          single real: on output stores the eigenvalues
  !> \param   q           single real matrix q: on output stores the eigenvalues
430
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
431
  abstract interface
432
    subroutine elpa_solve_f_i(self, a, ev, q, error)
433
434
      use iso_c_binding
      import elpa_t
435
      implicit none
436
437
438
439
440
441
442
443
      class(elpa_t)       :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)  :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
      real(kind=c_float)  :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
#endif
      real(kind=c_float)  :: ev(self%na)

444
      integer, optional   :: error
445
446
447
    end subroutine
  end interface

448
  !> \brief abstract definition of interface to solve double complex eigenvalue problem
449
450
451
452
453
454
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix a: defines the problem to solve
  !> \param   ev          double real: on output stores the eigenvalues
  !> \param   q           double complex matrix q: on output stores the eigenvalues
455
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
456
  abstract interface
457
    subroutine elpa_solve_dc_i(self, a, ev, q, error)
458
459
      use iso_c_binding
      import elpa_t
460
      implicit none
461
462
463
464
465
466
467
468
469
      class(elpa_t)                  :: self

#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
      complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
#endif
      real(kind=c_double)            :: ev(self%na)

470
      integer, optional              :: error
471
472
473
    end subroutine
  end interface

474
  !> \brief abstract definition of interface to solve single complex eigenvalue problem
475
476
477
478
479
480
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix a: defines the problem to solve
  !> \param   ev          single real: on output stores the eigenvalues
  !> \param   q           single complex matrix q: on output stores the eigenvalues
481
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
482
  abstract interface
483
    subroutine elpa_solve_fc_i(self, a, ev, q, error)
484
485
      use iso_c_binding
      import elpa_t
486
      implicit none
487
488
489
490
491
492
493
494
      class(elpa_t)                 :: self
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_float_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
#else
      complex(kind=c_float_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
#endif
      real(kind=c_float)            :: ev(self%na)

495
      integer, optional             :: error
496
497
498
    end subroutine
  end interface

499
  !> \brief abstract definition of interface to compute C : = A**T * B
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   uplo_a      'U' if A is upper triangular
  !>                      'L' if A is lower triangular
  !>                      anything else if A is a full matrix
  !>                      Please note: This pertains to the original A (as set in the calling program)
  !>                                   whereas the transpose of A is used for calculations
  !>                      If uplo_a is 'U' or 'L', the other triangle is not used at all,
  !>                      i.e. it may contain arbitrary numbers
  !> \param uplo_c        'U' if only the upper diagonal part of C is needed
  !>                      'L' if only the upper diagonal part of C is needed
  !>                      anything else if the full matrix C is needed
  !>                      Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
  !>                                    written to a certain extent, i.e. one shouldn't rely on the content there!
  !> \param na            Number of rows/columns of A, number of rows of B and C
  !> \param ncb           Number of real  of B and C
  !> \param   a           double complex matrix a
  !> \param lda           leading dimension of matrix a
  !> \param ldaCols       columns of matrix a
  !> \param b             double real matrix b
  !> \param ldb           leading dimension of matrix b
  !> \param ldbCols       columns of matrix b
  !> \param c             double real  matrix c
  !> \param ldc           leading dimension of matrix c
  !> \param ldcCols       columns of matrix c
526
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
527
  abstract interface
528
    subroutine elpa_hermitian_multiply_d_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
529
                                          c, ldc, ldcCols, error)
530
531
      use iso_c_binding
      import elpa_t
532
      implicit none
533
534
535
536
537
538
539
540
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
      integer(kind=c_int), intent(in) :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double)             :: a(lda,*), b(ldb,*), c(ldc,*)
#else
      real(kind=c_double)             :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
#endif
541
      integer, optional               :: error
542
543
544
    end subroutine
  end interface

545
  !> \brief abstract definition of interface to compute C : = A**T * B
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   uplo_a      'U' if A is upper triangular
  !>                      'L' if A is lower triangular
  !>                      anything else if A is a full matrix
  !>                      Please note: This pertains to the original A (as set in the calling program)
  !>                                   whereas the transpose of A is used for calculations
  !>                      If uplo_a is 'U' or 'L', the other triangle is not used at all,
  !>                      i.e. it may contain arbitrary numbers
  !> \param uplo_c        'U' if only the upper diagonal part of C is needed
  !>                      'L' if only the upper diagonal part of C is needed
  !>                      anything else if the full matrix C is needed
  !>                      Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
  !>                                    written to a certain extent, i.e. one shouldn't rely on the content there!
  !> \param na            Number of rows/columns of A, number of rows of B and C
  !> \param ncb           Number of real  of B and C
  !> \param   a           single complex matrix a
  !> \param lda           leading dimension of matrix a
  !> \param ldaCols       columns of matrix a
  !> \param b             single real matrix b
  !> \param ldb           leading dimension of matrix b
  !> \param ldbCols       columns of matrix b
  !> \param c             single real  matrix c
  !> \param ldc           leading dimension of matrix c
  !> \param ldcCols       columns of matrix c
572
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
573
  abstract interface
574
    subroutine elpa_hermitian_multiply_f_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
575
                                          c, ldc, ldcCols, error)
576
577
      use iso_c_binding
      import elpa_t
578
      implicit none
579
580
581
582
583
584
585
586
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
      integer(kind=c_int), intent(in) :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)              :: a(lda,*), b(ldb,*), c(ldc,*)
#else
      real(kind=c_float)              :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
#endif
587
      integer, optional               :: error
588
589
590
    end subroutine
  end interface

591
  !> \brief abstract definition of interface to compute C : = A**H * B
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   uplo_a      'U' if A is upper triangular
  !>                      'L' if A is lower triangular
  !>                      anything else if A is a full matrix
  !>                      Please note: This pertains to the original A (as set in the calling program)
  !>                                   whereas the transpose of A is used for calculations
  !>                      If uplo_a is 'U' or 'L', the other triangle is not used at all,
  !>                      i.e. it may contain arbitrary numbers
  !> \param uplo_c        'U' if only the upper diagonal part of C is needed
  !>                      'L' if only the upper diagonal part of C is needed
  !>                      anything else if the full matrix C is needed
  !>                      Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
  !>                                    written to a certain extent, i.e. one shouldn't rely on the content there!
  !> \param na            Number of rows/columns of A, number of rows of B and C
  !> \param ncb           Number of columns  of B and C
  !> \param   a           double complex matrix a
  !> \param lda           leading dimension of matrix a
  !> \param ldaCols       columns of matrix a
  !> \param b             double complex matrix b
  !> \param ldb           leading dimension of matrix b
  !> \param ldbCols       columns of matrix b
  !> \param c             double complex  matrix c
  !> \param ldc           leading dimension of matrix c
  !> \param ldcCols       columns of matrix c
618
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
619
  abstract interface
620
    subroutine elpa_hermitian_multiply_dc_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
621
                                          c, ldc, ldcCols, error)
622
623
      use iso_c_binding
      import elpa_t
624
      implicit none
625
626
627
628
629
630
631
632
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
      integer(kind=c_int), intent(in) :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex)  :: a(lda,*), b(ldb,*), c(ldc,*)
#else
      complex(kind=c_double_complex)  :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
#endif
633
      integer, optional               :: error
634
635
636
    end subroutine
  end interface

637
  !> \brief abstract definition of interface to compute C : = A**H * B
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   uplo_a      'U' if A is upper triangular
  !>                      'L' if A is lower triangular
  !>                      anything else if A is a full matrix
  !>                      Please note: This pertains to the original A (as set in the calling program)
  !>                                   whereas the transpose of A is used for calculations
  !>                      If uplo_a is 'U' or 'L', the other triangle is not used at all,
  !>                      i.e. it may contain arbitrary numbers
  !> \param uplo_c        'U' if only the upper diagonal part of C is needed
  !>                      'L' if only the upper diagonal part of C is needed
  !>                      anything else if the full matrix C is needed
  !>                      Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
  !>                                    written to a certain extent, i.e. one shouldn't rely on the content there!
  !> \param na            Number of rows/columns of A, number of rows of B and C
  !> \param ncb           Number of columns  of B and C
  !> \param   a           single complex matrix a
  !> \param lda           leading dimension of matrix a
  !> \param ldaCols       columns of matrix a
  !> \param b             single complex matrix b
  !> \param ldb           leading dimension of matrix b
  !> \param ldbCols       columns of matrix b
  !> \param c             single complex  matrix c
  !> \param ldc           leading dimension of matrix c
  !> \param ldcCols       columns of matrix c
664
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
665
  abstract interface
666
    subroutine elpa_hermitian_multiply_fc_i (self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
667
                                          c, ldc, ldcCols, error)
668
669
      use iso_c_binding
      import elpa_t
670
      implicit none
671
672
673
674
675
676
677
678
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
      integer(kind=c_int), intent(in) :: na, lda, ldaCols, ldb, ldbCols, ldc, ldcCols, ncb
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_float_complex)   :: a(lda,*), b(ldb,*), c(ldc,*)
#else
      complex(kind=c_float_complex)   :: a(lda,ldaCols), b(ldb,ldbCols), c(ldc,ldcCols)
#endif
679
      integer, optional               :: error
680
681
682
    end subroutine
  end interface

683
  !> \brief abstract definition of interface to do a cholesky decomposition of a double real matrix
684
685
686
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
687
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
688
  abstract interface
689
    subroutine elpa_cholesky_d_i (self, a, error)
690
691
      use iso_c_binding
      import elpa_t
692
      implicit none
693
694
695
696
697
698
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double)             :: a(self%local_nrows,*)
#else
      real(kind=c_double)             :: a(self%local_nrows,self%local_ncols)
#endif
699
      integer, optional               :: error
700
701
702
    end subroutine
  end interface

703
  !> \brief abstract definition of interface to do a cholesky decomposition of a single real matrix
704
705
706
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
707
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
708
  abstract interface
709
    subroutine elpa_cholesky_f_i(self, a, error)
710
711
      use iso_c_binding
      import elpa_t
712
      implicit none
713
714
715
716
717
718
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)              :: a(self%local_nrows,*)
#else
      real(kind=c_float)              :: a(self%local_nrows,self%local_ncols)
#endif
719
      integer, optional               :: error
720
721
722
    end subroutine
  end interface

723
  !> \brief abstract definition of interface to do a cholesky decomposition of a double complex matrix
724
725
726
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
727
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
728
  abstract interface
729
    subroutine elpa_cholesky_dc_i (self, a, error)
730
731
      use iso_c_binding
      import elpa_t
732
      implicit none
733
734
735
736
737
738
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex)  :: a(self%local_nrows,*)
#else
      complex(kind=c_double_complex)  :: a(self%local_nrows,self%local_ncols)
#endif
739
      integer, optional               :: error
740
741
742
    end subroutine
  end interface

743
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
744
745
746
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
747
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
748
  abstract interface
749
    subroutine elpa_cholesky_fc_i (self, a, error)
750
751
      use iso_c_binding
      import elpa_t
752
      implicit none
753
754
755
756
757
758
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_float_complex)   :: a(self%local_nrows,*)
#else
      complex(kind=c_float_complex)   :: a(self%local_nrows,self%local_ncols)
#endif
759
      integer, optional               :: error
760
761
762
    end subroutine
  end interface

763
  !> \brief abstract definition of interface to invert a triangular double real matrix
764
765
766
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
767
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
768
  abstract interface
769
    subroutine elpa_invert_trm_d_i (self, a, error)
770
771
      use iso_c_binding
      import elpa_t
772
      implicit none
773
774
775
776
777
778
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double)             :: a(self%local_nrows,*)
#else
      real(kind=c_double)             :: a(self%local_nrows,self%local_ncols)
#endif
779
      integer, optional               :: error
780
781
782
    end subroutine
  end interface

783
  !> \brief abstract definition of interface to invert a triangular single real matrix
784
785
786
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
787
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
788
  abstract interface
789
    subroutine elpa_invert_trm_f_i (self, a, error)
790
791
      use iso_c_binding
      import elpa_t
792
      implicit none
793
794
795
796
797
798
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)              :: a(self%local_nrows,*)
#else
      real(kind=c_float)              :: a(self%local_nrows,self%local_ncols)
#endif
799
      integer, optional               :: error
800
801
802
    end subroutine
  end interface

803
  !> \brief abstract definition of interface to invert a triangular double complex matrix
804
805
806
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
807
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
808
  abstract interface
809
    subroutine elpa_invert_trm_dc_i (self, a, error)
810
811
      use iso_c_binding
      import elpa_t
812
      implicit none
813
814
815
816
817
818
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_double_complex)  :: a(self%local_nrows,*)
#else
      complex(kind=c_double_complex)  :: a(self%local_nrows,self%local_ncols)
#endif
819
      integer, optional               :: error
820
821
822
    end subroutine
  end interface

823
  !> \brief abstract definition of interface to invert a triangular single complex matrix
824
825
826
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
827
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
828
  abstract interface
829
    subroutine elpa_invert_trm_fc_i (self, a, error)
830
831
      use iso_c_binding
      import elpa_t
832
      implicit none
833
834
835
836
837
838
      class(elpa_t)                   :: self
#ifdef USE_ASSUMED_SIZE
      complex(kind=c_float_complex)   :: a(self%local_nrows,*)
#else
      complex(kind=c_float_complex)   :: a(self%local_nrows,self%local_ncols)
#endif
839
      integer, optional               :: error
840
841
842
    end subroutine
  end interface

843
  !> \brief abstract definition of interface to solve the eigenvalue problem for a double-precision real valued tridiangular matrix
844
845
846
847
848
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   d           double real 1d array: the diagonal elements of a matrix defined in setup
  !> \param   e           double real 1d array: the subdiagonal elements of a matrix defined in setup
  !> \param   q           double real matrix: on output contains the eigenvectors
849
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
850
  abstract interface
851
    subroutine elpa_solve_tridi_d_i (self, d, e, q, error)
852
853
      use iso_c_binding
      import elpa_t
854
      implicit none
855
856
857
858
859
860
861
      class(elpa_t)                   :: self
      real(kind=c_double)             :: d(self%na), e(self%na)
#ifdef USE_ASSUMED_SIZE
      real(kind=c_double)             :: q(self%local_nrows,*)
#else
      real(kind=c_double)             :: q(self%local_nrows,self%local_ncols)
#endif
862
      integer, optional               :: error
863
864
865
    end subroutine
  end interface

866
  !> \brief abstract definition of interface to solve the eigenvalue problem for a single-precision real valued tridiangular matrix
867
868
869
870
871
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   d           single real 1d array: the diagonal elements of a matrix defined in setup
  !> \param   e           single real 1d array: the subdiagonal elements of a matrix defined in setup
  !> \param   q           single real matrix: on output contains the eigenvectors
872
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
873
  abstract interface
874
    subroutine elpa_solve_tridi_f_i (self, d, e, q, error)
875
876
      use iso_c_binding
      import elpa_t
877
      implicit none
878
879
880
881
882
883
884
      class(elpa_t)                   :: self
      real(kind=c_float)              :: d(self%na), e(self%na)
#ifdef USE_ASSUMED_SIZE
      real(kind=c_float)              :: q(self%local_nrows,*)
#else
      real(kind=c_float)              :: q(self%local_nrows,self%local_ncols)
#endif
885
      integer, optional               :: error
886
887
888
    end subroutine
  end interface

889
  !> \brief abstract definition of interface to destroy an ELPA object
890
891
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
892
893
894
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
895
      implicit none
896
897
898
899
900
901
      class(elpa_t) :: self
    end subroutine
  end interface

  contains

902
903
    !> \brief function to intialise the ELPA library
    !> Parameters
904
905
    !> \param   api_version integer: api_version that ELPA should use
    !> \result  error       integer: error code, which can be queried with elpa_strerr
906
907
908
    !
    !c> int elpa_init(int api_version);
    function elpa_init(api_version) result(error) bind(C, name="elpa_init")
909
      use elpa_utilities, only : error_unit
910
      integer, intent(in), value :: api_version
911
912
913
914
915
916
917
918
919
920
921
      integer             :: error

      if (earliest_api_version <= api_version .and. api_version <= current_api_version) then
        initDone = .true.
        error = ELPA_OK
      else
        write(error_unit, "(a,i0,a)") "ELPA: Error API version ", api_version," is not supported by this library"
        error = ELPA_ERROR
      endif
    end function

922
923
    !> \brief function to check whether the ELPA library has been correctly initialised
    !> Parameters
924
    !> \result  state      logical: state is either ELPA_OK or ELPA_ERROR, which can be queried with elpa_strerr
925
    function elpa_initialized() result(state)
926
927
928
929
930
931
      integer :: state
      if (initDone) then
        state = ELPA_OK
      else
        state = ELPA_ERROR
      endif
932
933
    end function

934
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
935
936
937
    !
    !c> void elpa_uninit(void);
    subroutine elpa_uninit() bind(C, name="elpa_uninit")
938
939
    end subroutine

940
    !> \brief helper function for error strings
941
    !> Parameters
942
943
    !> \param   elpa_error  integer: error code to querry
    !> \result  string      string:  error string
944
945
946
947
948
949
950
    function elpa_strerr(elpa_error) result(string)
      use elpa_generated_fortran_interfaces
      integer, intent(in) :: elpa_error
      character(kind=C_CHAR, len=elpa_strlen_c(elpa_strerr_c(elpa_error))), pointer :: string
      call c_f_pointer(elpa_strerr_c(elpa_error), string)
    end function

951
    !> \brief helper function for c strings
952
953
954
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
955
956
957
958
959
960
961
    function elpa_c_string(ptr) result(string)
      use, intrinsic :: iso_c_binding
      type(c_ptr), intent(in) :: ptr
      character(kind=c_char, len=elpa_strlen_c(ptr)), pointer :: string
      call c_f_pointer(ptr, string)
    end function

962
    !> \brief function to convert an integer in its string representation
963
    !> Parameters
964
965
966
967
    !> \param   name        string: the key
    !> \param   value       integer: the value correponding to the key
    !> \param   error       integer, optional: error code, which can be queried with elpa_strerr()
    !> \result  string      string: the string representation
968
969
970
971
972
973
974
    function elpa_int_value_to_string(name, value, error) result(string)
      use elpa_utilities, only : error_unit
      use elpa_generated_fortran_interfaces
      implicit none
      character(kind=c_char, len=*), intent(in) :: name
      integer(kind=c_int), intent(in) :: value
      integer(kind=c_int), intent(out), optional :: error
975
      character(kind=c_char, len=elpa_int_value_to_strlen_c(name // C_NULL_CHAR, value)), pointer :: string
976
977
978
      integer(kind=c_int) :: actual_error
      type(c_ptr) :: ptr

979
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
980
981
982
983
984
985
986
987
988
989
990
991
992
993
      if (c_associated(ptr)) then
        call c_f_pointer(ptr, string)
      else
        nullify(string)
      endif

      if (present(error)) then
        error = actual_error
      else if (actual_error /= ELPA_OK) then
        write(error_unit,'(a,i0,a)') "ELPA: Error converting value '", value, "' to a string for option '" // &
                name // "' and you did not check for errors: " // elpa_strerr(actual_error)
      endif
    end function

994
    !> \brief function to convert a string in its integer representation:
995
    !> Parameters
996
997
998
999
    !> \param   name        string: the key
    !> \param   string      string: the string whose integer representation should be associated with the key
    !> \param   error       integer, optional: error code, which can be queried with elpa_strerr()
    !> \result  value       integer: the integer representation of the string
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
    function elpa_int_string_to_value(name, string, error) result(value)
      use elpa_generated_fortran_interfaces
      use elpa_utilities, only : error_unit
      implicit none
      character(kind=c_char, len=*), intent(in) :: name
      character(kind=c_char, len=*), intent(in), target :: string
      integer(kind=c_int), intent(out), optional :: error
      integer(kind=c_int) :: actual_error

      integer(kind=c_int) :: value
      integer(kind=c_int)   :: i
      type(c_ptr) :: repr

1013
      actual_error = elpa_int_string_to_value_c(name // C_NULL_CHAR, string // C_NULL_CHAR, value)
1014
1015
1016
1017
1018
1019
1020
1021
1022

      if (present(error)) then
        error = actual_error
      else if (actual_error /= ELPA_OK) then
        write(error_unit,'(a)') "ELPA: Error converting string '" // string // "' to value for option '" // &
                name // "' and you did not check for errors: " // elpa_strerr(actual_error)
      endif
    end function

1023
    !> \brief function to get the number of possible choices for an option
1024
    !> Parameters
1025
1026
    !> \param   option_name string:   the option
    !> \result  number      integer:  the total number of possible values to be chosen
1027
1028
1029
1030
1031
1032
1033
    function elpa_option_cardinality(option_name) result(number)
      use elpa_generated_fortran_interfaces
      character(kind=c_char, len=*), intent(in) :: option_name
      integer                                   :: number
      number = elpa_option_cardinality_c(option_name // C_NULL_CHAR)
    end function

1034
    !> \brief function to enumerate an option
1035
    !> Parameters
1036
    !> \param   option_name string: the option
1037
1038
    !> \param   i           integer
    !> \result  option      integer
1039
1040
1041
1042
1043
1044
1045
1046
    function elpa_option_enumerate(option_name, i) result(option)
      use elpa_generated_fortran_interfaces
      character(kind=c_char, len=*), intent(in) :: option_name
      integer, intent(in)                       :: i
      integer                                   :: option
      option = elpa_option_enumerate_c(option_name // C_NULL_CHAR, i)
    end function

1047
end module