elpa_api.F90 42.7 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
!> \brief Fortran module which provides the definition of the ELPA API
50
51
52
53
54
module elpa_api
  use elpa_constants
  use, intrinsic :: iso_c_binding
  implicit none

55
56
57
58
  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

59
60
61
62
63
64
65
  logical, private :: initDone = .false.

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

66
  !> \brief Abstract definition of the elpa_t type
67
68
69
  type, abstract :: elpa_t
    private

70
    !< these have to be public for proper bounds checking, sadly
71
72
73
74
75
76
77
    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
78
      !> \brief methods available with the elpa_t type
79
      ! general
80
81
      procedure(elpa_setup_i),   deferred, public :: setup          !< export a setup method
      procedure(elpa_destroy_i), deferred, public :: destroy        !< export a destroy method
82

83
      !< key/value store
84
      generic, public :: set => &                                   !< export a method to set integer/double key/values
85
86
          elpa_set_integer, &
          elpa_set_double
87
88
      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
89

90
91
      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
92

93
      !< Timer
94
95
96
      procedure(elpa_get_time_i), deferred, public :: get_time
      procedure(elpa_print_times_i), deferred, public :: print_times

97
      !< Actual math routines
98
      generic, public :: solve => &                                 !< method solve for solving the eigenvalue problem
99
100
101
102
          elpa_solve_d, &                                 !< for symmetric real valued / hermitian complex valued
          elpa_solve_f, &                                 !< matrices
          elpa_solve_dc, &
          elpa_solve_fc
103

104
      generic, public :: hermitian_multiply => &                    !< method for a "hermitian" multiplication of matrices a and b
105
106
107
108
          elpa_hermitian_multiply_d, &                              !< for real valued matrices:   a**T * b
          elpa_hermitian_multiply_dc, &                              !< for complex valued matrices a**H * b
          elpa_hermitian_multiply_f, &
          elpa_hermitian_multiply_fc
109

110
      generic, public :: cholesky => &                              !< method for the cholesky factorisation of matrix a
111
112
113
114
          elpa_cholesky_d, &
          elpa_cholesky_f, &
          elpa_cholesky_dc, &
          elpa_cholesky_fc
115

Andreas Marek's avatar
Andreas Marek committed
116
      generic, public :: invert_triangular => &                     !< method to invert a upper triangular matrix a
117
118
119
120
          elpa_invert_trm_d, &
          elpa_invert_trm_f, &
          elpa_invert_trm_dc, &
          elpa_invert_trm_fc
121

122
      generic, public :: solve_tridi => &                           !< method to solve the eigenvalue problem for a tridiagonal
123
124
          elpa_solve_tridi_d, &                           !< matrix
          elpa_solve_tridi_f
125
126


127
      !> \brief private methods of elpa_t type. NOT accessible for the user
128
129
130
131
      ! privates
      procedure(elpa_set_integer_i), deferred, private :: elpa_set_integer
      procedure(elpa_set_double_i),  deferred, private :: elpa_set_double

132
133
134
135
      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
136

137
138
139
140
      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
141

142
143
144
145
      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
146

147
148
149
150
      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
151

152
153
      procedure(elpa_solve_tridi_d_i), deferred, private :: elpa_solve_tridi_d
      procedure(elpa_solve_tridi_f_i), deferred, private :: elpa_solve_tridi_f
154
155
156
  end type elpa_t


157
158
159
160
161
  !> \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
162
163
164
  interface
    pure function elpa_strlen_c(ptr) result(size) bind(c, name="strlen")
      use, intrinsic :: iso_c_binding
165
      implicit none
166
167
168
169
170
      type(c_ptr), intent(in), value :: ptr
      integer(kind=c_size_t) :: size
    end function
  end interface

171
172
173
174
175
  !> \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()
176
  abstract interface
177
    function elpa_setup_i(self) result(error)
178
      import elpa_t
179
      implicit none
180
      class(elpa_t), intent(inout) :: self
181
      integer :: error
182
183
184
    end function
  end interface

185
186
187
188
189
190
191
  !> \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()
192
  abstract interface
193
    subroutine elpa_set_integer_i(self, name, value, error)
194
195
      use iso_c_binding
      import elpa_t
196
      implicit none
197
198
199
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
200
      integer, optional               :: error
201
202
203
    end subroutine
  end interface

204
205
206
207
208
209
210
  !> \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
211
  abstract interface
212
    function elpa_get_integer_i(self, name, error) result(value)
213
214
      use iso_c_binding
      import elpa_t
215
      implicit none
216
217
218
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
219
      integer, intent(out), optional :: error
220
221
222
    end function
  end interface

223
224
225
226
227
228
229
  !> \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
230
  abstract interface
231
    function elpa_is_set_i(self, name) result(state)
232
      import elpa_t
233
      implicit none
234
235
      class(elpa_t)            :: self
      character(*), intent(in) :: name
236
      integer                  :: state
237
238
239
    end function
  end interface

240
241
242
243
244
245
246
247
  !> \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
248
  abstract interface
249
    function elpa_can_set_i(self, name, value) result(state)
250
      import elpa_t, c_int
251
      implicit none
252
253
254
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
255
      integer                         :: state
256
    end function
257
258
  end interface

259
260
261
262
263
264
265
  !> \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
266
  abstract interface
267
    subroutine elpa_set_double_i(self, name, value, error)
268
269
      use iso_c_binding
      import elpa_t
270
      implicit none
271
272
273
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      real(kind=c_double), intent(in) :: value
274
      integer, optional               :: error
275
276
277
    end subroutine
  end interface

278
279
280
281
282
283
284
  !> \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
285
  abstract interface
286
    function elpa_get_double_i(self, name, error) result(value)
287
288
      use iso_c_binding
      import elpa_t
289
      implicit none
290
291
292
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      real(kind=c_double)            :: value
293
      integer, intent(out), optional :: error
294
295
296
    end function
  end interface

297
298
299
300
301
302
  !> \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
303
304
305
306
  abstract interface
    function elpa_associate_int_i(self, name) result(value)
      use iso_c_binding
      import elpa_t
307
      implicit none
308
309
310
311
312
313
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
    end function
  end interface

314
315
316

  ! Timer routines

317
318
319
320
321
322
  !> \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
323
324
325
  abstract interface
    function elpa_get_time_i(self, name1, name2, name3, name4, name5, name6) result(s)
      import elpa_t, c_double
326
      implicit none
327
328
329
330
331
332
333
      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

334
335
336
337
  !> \brief abstract definition of print method for timer
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
338
339
340
  abstract interface
    subroutine elpa_print_times_i(self)
      import elpa_t
341
      implicit none
342
343
344
345
346
      class(elpa_t), intent(in) :: self
    end subroutine
  end interface


347
  !< Actual math routines
348

349
  !> \brief abstract definition of interface to solve double real eigenvalue problem
350
351
352
353
354
355
  !> 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
356
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
357
  abstract interface
358
    subroutine elpa_solve_d_i(self, a, ev, q, error)
359
360
      use iso_c_binding
      import elpa_t
361
      implicit none
362
363
364
365
366
367
368
369
      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)

370
      integer, optional   :: error
371
372
373
    end subroutine
  end interface

374
  !> \brief abstract definition of interface to solve single real eigenvalue problem
375
376
377
378
379
380
  !> 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
381
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
382
  abstract interface
383
    subroutine elpa_solve_f_i(self, a, ev, q, error)
384
385
      use iso_c_binding
      import elpa_t
386
      implicit none
387
388
389
390
391
392
393
394
      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)

395
      integer, optional   :: error
396
397
398
    end subroutine
  end interface

399
  !> \brief abstract definition of interface to solve double complex eigenvalue problem
400
401
402
403
404
405
  !> 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
406
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
407
  abstract interface
408
    subroutine elpa_solve_dc_i(self, a, ev, q, error)
409
410
      use iso_c_binding
      import elpa_t
411
      implicit none
412
413
414
415
416
417
418
419
420
      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)

421
      integer, optional              :: error
422
423
424
    end subroutine
  end interface

425
  !> \brief abstract definition of interface to solve single complex eigenvalue problem
426
427
428
429
430
431
  !> 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
432
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
433
  abstract interface
434
    subroutine elpa_solve_fc_i(self, a, ev, q, error)
435
436
      use iso_c_binding
      import elpa_t
437
      implicit none
438
439
440
441
442
443
444
445
      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)

446
      integer, optional             :: error
447
448
449
    end subroutine
  end interface

450
  !> \brief abstract definition of interface to compute C : = A**T * B
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
  !> 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
477
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
478
  abstract interface
479
    subroutine elpa_hermitian_multiply_d_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
480
                                          c, ldc, ldcCols, error)
481
482
      use iso_c_binding
      import elpa_t
483
      implicit none
484
485
486
487
488
489
490
491
      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
492
      integer, optional               :: error
493
494
495
    end subroutine
  end interface

496
  !> \brief abstract definition of interface to compute C : = A**T * B
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
  !> 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
523
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
524
  abstract interface
525
    subroutine elpa_hermitian_multiply_f_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
526
                                          c, ldc, ldcCols, error)
527
528
      use iso_c_binding
      import elpa_t
529
      implicit none
530
531
532
533
534
535
536
537
      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
538
      integer, optional               :: error
539
540
541
    end subroutine
  end interface

542
  !> \brief abstract definition of interface to compute C : = A**H * B
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
  !> 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
569
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
570
  abstract interface
571
    subroutine elpa_hermitian_multiply_dc_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
572
                                          c, ldc, ldcCols, error)
573
574
      use iso_c_binding
      import elpa_t
575
      implicit none
576
577
578
579
580
581
582
583
      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
584
      integer, optional               :: error
585
586
587
    end subroutine
  end interface

588
  !> \brief abstract definition of interface to compute C : = A**H * B
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
  !> 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
615
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
616
  abstract interface
617
    subroutine elpa_hermitian_multiply_fc_i (self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
618
                                          c, ldc, ldcCols, error)
619
620
      use iso_c_binding
      import elpa_t
621
      implicit none
622
623
624
625
626
627
628
629
      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
630
      integer, optional               :: error
631
632
633
    end subroutine
  end interface

634
  !> \brief abstract definition of interface to do a cholesky decomposition of a double real matrix
635
636
637
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
638
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
639
  abstract interface
640
    subroutine elpa_cholesky_d_i (self, a, error)
641
642
      use iso_c_binding
      import elpa_t
643
      implicit none
644
645
646
647
648
649
      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
650
      integer, optional               :: error
651
652
653
    end subroutine
  end interface

654
  !> \brief abstract definition of interface to do a cholesky decomposition of a single real matrix
655
656
657
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
658
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
659
  abstract interface
660
    subroutine elpa_cholesky_f_i(self, a, error)
661
662
      use iso_c_binding
      import elpa_t
663
      implicit none
664
665
666
667
668
669
      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
670
      integer, optional               :: error
671
672
673
    end subroutine
  end interface

674
  !> \brief abstract definition of interface to do a cholesky decomposition of a double complex matrix
675
676
677
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
678
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
679
  abstract interface
680
    subroutine elpa_cholesky_dc_i (self, a, error)
681
682
      use iso_c_binding
      import elpa_t
683
      implicit none
684
685
686
687
688
689
      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
690
      integer, optional               :: error
691
692
693
    end subroutine
  end interface

694
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
695
696
697
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
698
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
699
  abstract interface
700
    subroutine elpa_cholesky_fc_i (self, a, error)
701
702
      use iso_c_binding
      import elpa_t
703
      implicit none
704
705
706
707
708
709
      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
710
      integer, optional               :: error
711
712
713
    end subroutine
  end interface

714
  !> \brief abstract definition of interface to invert a triangular double real matrix
715
716
717
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
718
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
719
  abstract interface
720
    subroutine elpa_invert_trm_d_i (self, a, error)
721
722
      use iso_c_binding
      import elpa_t
723
      implicit none
724
725
726
727
728
729
      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
730
      integer, optional               :: error
731
732
733
    end subroutine
  end interface

734
  !> \brief abstract definition of interface to invert a triangular single real matrix
735
736
737
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
738
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
739
  abstract interface
740
    subroutine elpa_invert_trm_f_i (self, a, error)
741
742
      use iso_c_binding
      import elpa_t
743
      implicit none
744
745
746
747
748
749
      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
750
      integer, optional               :: error
751
752
753
    end subroutine
  end interface

754
  !> \brief abstract definition of interface to invert a triangular double complex matrix
755
756
757
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
758
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
759
  abstract interface
760
    subroutine elpa_invert_trm_dc_i (self, a, error)
761
762
      use iso_c_binding
      import elpa_t
763
      implicit none
764
765
766
767
768
769
      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
770
      integer, optional               :: error
771
772
773
    end subroutine
  end interface

774
  !> \brief abstract definition of interface to invert a triangular single complex matrix
775
776
777
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
778
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
779
  abstract interface
780
    subroutine elpa_invert_trm_fc_i (self, a, error)
781
782
      use iso_c_binding
      import elpa_t
783
      implicit none
784
785
786
787
788
789
      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
790
      integer, optional               :: error
791
792
793
    end subroutine
  end interface

794
  !> \brief abstract definition of interface to solve the eigenvalue problem for a double-precision real valued tridiangular matrix
795
796
797
798
799
  !> 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
800
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
801
  abstract interface
802
    subroutine elpa_solve_tridi_d_i (self, d, e, q, error)
803
804
      use iso_c_binding
      import elpa_t
805
      implicit none
806
807
808
809
810
811
812
      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
813
      integer, optional               :: error
814
815
816
    end subroutine
  end interface

817
  !> \brief abstract definition of interface to solve the eigenvalue problem for a single-precision real valued tridiangular matrix
818
819
820
821
822
  !> 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
823
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
824
  abstract interface
825
    subroutine elpa_solve_tridi_f_i (self, d, e, q, error)
826
827
      use iso_c_binding
      import elpa_t
828
      implicit none
829
830
831
832
833
834
835
      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
836
      integer, optional               :: error
837
838
839
    end subroutine
  end interface

840
  !> \brief abstract definition of interface to destroy an ELPA object
841
842
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
843
844
845
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
846
      implicit none
847
848
849
850
851
852
      class(elpa_t) :: self
    end subroutine
  end interface

  contains

853
854
    !> \brief function to intialise the ELPA library
    !> Parameters
855
856
    !> \param   api_version integer: api_version that ELPA should use
    !> \result  error       integer: error code, which can be queried with elpa_strerr
857
858
859
    !
    !c> int elpa_init(int api_version);
    function elpa_init(api_version) result(error) bind(C, name="elpa_init")
860
      use elpa_utilities, only : error_unit
861
      integer, intent(in), value :: api_version
862
863
864
865
866
867
868
869
870
871
872
      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

873
874
    !> \brief function to check whether the ELPA library has been correctly initialised
    !> Parameters
875
    !> \result  state      logical: state is either ELPA_OK or ELPA_ERROR, which can be queried with elpa_strerr
876
    function elpa_initialized() result(state)
877
878
879
880
881
882
      integer :: state
      if (initDone) then
        state = ELPA_OK
      else
        state = ELPA_ERROR
      endif
883
884
    end function

885
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
886
887
888
    !
    !c> void elpa_uninit(void);
    subroutine elpa_uninit() bind(C, name="elpa_uninit")
889
890
    end subroutine

891
    !> \brief helper function for error strings
892
    !> Parameters
893
894
    !> \param   elpa_error  integer: error code to querry
    !> \result  string      string:  error string
895
896
897
898
899
900
901
    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

902
    !> \brief helper function for c strings
903
904
905
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
906
907
908
909
910
911
912
    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

913
    !> \brief function to convert an integer in its string representation
914
    !> Parameters
915
916
917
918
    !> \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
919
920
921
922
923
924
925
    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
926
      character(kind=c_char, len=elpa_int_value_to_strlen_c(name // C_NULL_CHAR, value)), pointer :: string
927
928
929
      integer(kind=c_int) :: actual_error
      type(c_ptr) :: ptr

930
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
931
932
933
934
935
936
937
938
939
940
941
942
943
944
      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

945
    !> \brief function to convert a string in its integer representation:
946
    !> Parameters
947
948
949
950
    !> \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
951
952
953
954
955
956
957
958
959
960
961
962
963
    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

964
      actual_error = elpa_int_string_to_value_c(name // C_NULL_CHAR, string // C_NULL_CHAR, value)
965
966
967
968
969
970
971
972
973

      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

974
    !> \brief function to get the number of possible choices for an option
975
    !> Parameters
976
977
    !> \param   option_name string:   the option
    !> \result  number      integer:  the total number of possible values to be chosen
978
979
980
981
982
983
984
    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

985
    !> \brief function to enumerate an option
986
    !> Parameters
987
    !> \param   option_name string: the option
988
989
    !> \param   i           integer
    !> \result  option      integer
990
991
992
993
994
995
996
997
    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

998
end module