elpa_api.F90 37.5 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
module elpa_api
  use elpa_constants
  use, intrinsic :: iso_c_binding
53
  use ftimings
54
55
  implicit none

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

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

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

67
  !> \brief Abstract defintion of the elpa_t type
68
69
70
71
72
73
74
75
76
77
  type, abstract :: elpa_t
    private

    ! \todo: it's ugly that these are public
    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()

78
    type(timer_t), public :: timer
79
    contains
80
      !> \brief methods available with the elpa_t type
81
      ! general
82
83
      procedure(elpa_setup_i),   deferred, public :: setup          !< export a setup method
      procedure(elpa_destroy_i), deferred, public :: destroy        !< export a destroy method
84
85

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

92
93
      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
94

95
      ! actual math routines
96
97
98
      generic, public :: solve => &                                 !< method solve for solving the eigenvalue problem
          elpa_solve_real_double, &                                 !< for symmetric real valued / hermitian complex valued
          elpa_solve_real_single, &                                 !< matrices
99
100
101
          elpa_solve_complex_double, &
          elpa_solve_complex_single

102
103
104
      generic, public :: hermitian_multiply => &                    !< method for a "hermitian" multiplication of matrices a and b
          elpa_multiply_at_b_double, &                              !< for real valued matrices:   a**T * b
          elpa_multiply_ah_b_double, &                              !< for complex valued matrices a**H * b
105
106
107
          elpa_multiply_at_b_single, &
          elpa_multiply_ah_b_single

108
      generic, public :: cholesky => &                              !< method for the cholesky factorisation of matrix a
109
110
111
112
113
          elpa_cholesky_double_real, &
          elpa_cholesky_single_real, &
          elpa_cholesky_double_complex, &
          elpa_cholesky_single_complex

Andreas Marek's avatar
Andreas Marek committed
114
      generic, public :: invert_triangular => &                     !< method to invert a upper triangular matrix a
115
116
117
118
119
          elpa_invert_trm_double_real, &
          elpa_invert_trm_single_real, &
          elpa_invert_trm_double_complex, &
          elpa_invert_trm_single_complex

120
121
      generic, public :: solve_tridi => &                           !< method to solve the eigenvalue problem for a tridiagonal
          elpa_solve_tridi_double_real, &                           !< matrix
122
123
124
          elpa_solve_tridi_single_real


125
      !> \brief private methods of elpa_t type. NOT accessible for the user
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
      ! privates
      procedure(elpa_set_integer_i), deferred, private :: elpa_set_integer
      procedure(elpa_set_double_i),  deferred, private :: elpa_set_double

      procedure(elpa_solve_real_double_i),    deferred, private :: elpa_solve_real_double
      procedure(elpa_solve_real_single_i),    deferred, private :: elpa_solve_real_single
      procedure(elpa_solve_complex_double_i), deferred, private :: elpa_solve_complex_double
      procedure(elpa_solve_complex_single_i), deferred, private :: elpa_solve_complex_single

      procedure(elpa_multiply_at_b_double_i), deferred, private :: elpa_multiply_at_b_double
      procedure(elpa_multiply_at_b_single_i), deferred, private :: elpa_multiply_at_b_single
      procedure(elpa_multiply_ah_b_double_i), deferred, private :: elpa_multiply_ah_b_double
      procedure(elpa_multiply_ah_b_single_i), deferred, private :: elpa_multiply_ah_b_single

      procedure(elpa_cholesky_double_real_i),    deferred, private :: elpa_cholesky_double_real
      procedure(elpa_cholesky_single_real_i),    deferred, private :: elpa_cholesky_single_real
      procedure(elpa_cholesky_double_complex_i), deferred, private :: elpa_cholesky_double_complex
      procedure(elpa_cholesky_single_complex_i), deferred, private :: elpa_cholesky_single_complex

      procedure(elpa_invert_trm_double_real_i),    deferred, private :: elpa_invert_trm_double_real
      procedure(elpa_invert_trm_single_real_i),    deferred, private :: elpa_invert_trm_single_real
      procedure(elpa_invert_trm_double_complex_i), deferred, private :: elpa_invert_trm_double_complex
      procedure(elpa_invert_trm_single_complex_i), deferred, private :: elpa_invert_trm_single_complex

      procedure(elpa_solve_tridi_double_real_i), deferred, private :: elpa_solve_tridi_double_real
      procedure(elpa_solve_tridi_single_real_i), deferred, private :: elpa_solve_tridi_single_real
  end type elpa_t


  interface
    pure function elpa_strlen_c(ptr) result(size) bind(c, name="strlen")
      use, intrinsic :: iso_c_binding
      type(c_ptr), intent(in), value :: ptr
      integer(kind=c_size_t) :: size
    end function
  end interface


  abstract interface
165
    function elpa_setup_i(self) result(error)
166
167
      import elpa_t
      class(elpa_t), intent(inout) :: self
168
      integer :: error
169
170
171
172
173
    end function
  end interface


  abstract interface
174
    subroutine elpa_set_integer_i(self, name, value, error)
175
176
177
178
179
      use iso_c_binding
      import elpa_t
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
180
      integer, optional               :: error
181
182
183
184
185
    end subroutine
  end interface


  abstract interface
186
    function elpa_get_integer_i(self, name, error) result(value)
187
188
189
190
191
      use iso_c_binding
      import elpa_t
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
192
      integer, intent(out), optional :: error
193
194
195
196
197
    end function
  end interface


  abstract interface
198
    function elpa_is_set_i(self, name) result(error)
199
200
201
      import elpa_t
      class(elpa_t)            :: self
      character(*), intent(in) :: name
202
      integer                  :: error
203
204
205
206
207
    end function
  end interface


  abstract interface
208
209
210
211
212
213
214
    function elpa_can_set_i(self, name, value) result(error)
      import elpa_t, c_int
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
      integer                         :: error
    end function
215
216
217
218
  end interface


  abstract interface
219
    subroutine elpa_set_double_i(self, name, value, error)
220
221
222
223
224
      use iso_c_binding
      import elpa_t
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      real(kind=c_double), intent(in) :: value
225
      integer, optional               :: error
226
227
228
229
230
    end subroutine
  end interface


  abstract interface
231
    function elpa_get_double_i(self, name, error) result(value)
232
233
234
235
236
      use iso_c_binding
      import elpa_t
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      real(kind=c_double)            :: value
237
      integer, intent(out), optional :: error
238
239
240
241
242
243
244
245
246
247
248
249
250
251
    end function
  end interface


  abstract interface
    function elpa_associate_int_i(self, name) result(value)
      use iso_c_binding
      import elpa_t
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
    end function
  end interface

252
253
254
255
256
257
258
  !> \brief abstract defintion of interface to solve double real eigenvalue problem
  !> 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
259
  abstract interface
260
    subroutine elpa_solve_real_double_i(self, a, ev, q, error)
261
262
263
264
265
266
267
268
269
270
      use iso_c_binding
      import elpa_t
      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)

271
      integer, optional   :: error
272
273
274
    end subroutine
  end interface

275
276
277
278
279
280
281
  !> \brief abstract defintion of interface to solve single real eigenvalue problem
  !> 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
282
  abstract interface
283
    subroutine elpa_solve_real_single_i(self, a, ev, q, error)
284
285
286
287
288
289
290
291
292
293
      use iso_c_binding
      import elpa_t
      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)

294
      integer, optional   :: error
295
296
297
    end subroutine
  end interface

298
299
300
301
302
303
304
  !> \brief abstract defintion of interface to solve double complex eigenvalue problem
  !> 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
305
  abstract interface
306
    subroutine elpa_solve_complex_double_i(self, a, ev, q, error)
307
308
309
310
311
312
313
314
315
316
317
      use iso_c_binding
      import elpa_t
      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)

318
      integer, optional              :: error
319
320
321
    end subroutine
  end interface

322
323
324
325
326
327
328
  !> \brief abstract defintion of interface to solve single complex eigenvalue problem
  !> 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
329
  abstract interface
330
    subroutine elpa_solve_complex_single_i(self, a, ev, q, error)
331
332
333
334
335
336
337
338
339
340
      use iso_c_binding
      import elpa_t
      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)

341
      integer, optional             :: error
342
343
344
    end subroutine
  end interface

345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
  !> \brief abstract defintion of interface to compute C : = A**T * B
  !> 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
  !> \param   error       integer, optional : error code
373
374
  abstract interface
    subroutine elpa_multiply_at_b_double_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
375
                                          c, ldc, ldcCols, error)
376
377
378
379
380
381
382
383
384
385
      use iso_c_binding
      import elpa_t
      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
386
      integer, optional               :: error
387
388
389
    end subroutine
  end interface

390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
  !> \brief abstract defintion of interface to compute C : = A**T * B
  !> 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
  !> \param   error       integer, optional : error code
418
419
  abstract interface
    subroutine elpa_multiply_at_b_single_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
420
                                          c, ldc, ldcCols, error)
421
422
423
424
425
426
427
428
429
430
      use iso_c_binding
      import elpa_t
      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
431
      integer, optional               :: error
432
433
434
    end subroutine
  end interface

435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
  !> \brief abstract defintion of interface to compute C : = A**H * B
  !> 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
  !> \param   error       integer, optional : error code
463
464
  abstract interface
    subroutine elpa_multiply_ah_b_double_i (self,uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
465
                                          c, ldc, ldcCols, error)
466
467
468
469
470
471
472
473
474
475
      use iso_c_binding
      import elpa_t
      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
476
      integer, optional               :: error
477
478
479
    end subroutine
  end interface

480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
  !> \brief abstract defintion of interface to compute C : = A**H * B
  !> 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
  !> \param   error       integer, optional : error code
508
  abstract interface
509
    subroutine elpa_multiply_ah_b_single_i (self, uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, &
510
                                          c, ldc, ldcCols, error)
511
512
513
514
515
516
517
518
519
520
      use iso_c_binding
      import elpa_t
      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
521
      integer, optional               :: error
522
523
524
    end subroutine
  end interface

525
526
527
528
529
  !> \brief abstract defintion of interface to do a cholesky decomposition of a double real matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
  !> \param   error       integer, optional : error code
530
  abstract interface
531
    subroutine elpa_cholesky_double_real_i (self, a, error)
532
533
534
535
536
537
538
539
      use iso_c_binding
      import elpa_t
      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
540
      integer, optional               :: error
541
542
543
    end subroutine
  end interface

544
545
546
547
548
  !> \brief abstract defintion of interface to do a cholesky decomposition of a single real matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
  !> \param   error       integer, optional : error code
549
  abstract interface
550
    subroutine elpa_cholesky_single_real_i(self, a, error)
551
552
553
554
555
556
557
558
      use iso_c_binding
      import elpa_t
      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
559
      integer, optional               :: error
560
561
562
    end subroutine
  end interface

563
564
565
566
567
  !> \brief abstract defintion of interface to do a cholesky decomposition of a double complex matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
  !> \param   error       integer, optional : error code
568
  abstract interface
569
    subroutine elpa_cholesky_double_complex_i (self, a, error)
570
571
572
573
574
575
576
577
      use iso_c_binding
      import elpa_t
      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
578
      integer, optional               :: error
579
580
581
    end subroutine
  end interface

582
583
584
585
586
  !> \brief abstract defintion of interface to do a cholesky decomposition of a single complex matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
  !> \param   error       integer, optional : error code
587
  abstract interface
588
    subroutine elpa_cholesky_single_complex_i (self, a, error)
589
590
591
592
593
594
595
596
      use iso_c_binding
      import elpa_t
      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
597
      integer, optional               :: error
598
599
600
    end subroutine
  end interface

601
602
603
604
605
  !> \brief abstract defintion of interface to invert a triangular double real matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
  !> \param   error       integer, optional : error code
606
  abstract interface
607
    subroutine elpa_invert_trm_double_real_i (self, a, error)
608
609
610
611
612
613
614
615
      use iso_c_binding
      import elpa_t
      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
616
      integer, optional               :: error
617
618
619
    end subroutine
  end interface

620
621
622
623
624
  !> \brief abstract defintion of interface to invert a triangular single real matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
  !> \param   error       integer, optional : error code
625
  abstract interface
626
    subroutine elpa_invert_trm_single_real_i (self, a, error)
627
628
629
630
631
632
633
634
      use iso_c_binding
      import elpa_t
      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
635
      integer, optional               :: error
636
637
638
    end subroutine
  end interface

639
640
641
642
643
  !> \brief abstract defintion of interface to invert a triangular double complex matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
  !> \param   error       integer, optional : error code
644
  abstract interface
645
    subroutine elpa_invert_trm_double_complex_i (self, a, error)
646
647
648
649
650
651
652
653
      use iso_c_binding
      import elpa_t
      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
654
      integer, optional               :: error
655
656
657
    end subroutine
  end interface

658
659
660
661
662
  !> \brief abstract defintion of interface to invert a triangular single complex matrix
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
  !> \param   error       integer, optional : error code
663
  abstract interface
664
    subroutine elpa_invert_trm_single_complex_i (self, a, error)
665
666
667
668
669
670
671
672
      use iso_c_binding
      import elpa_t
      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
673
      integer, optional               :: error
674
675
676
    end subroutine
  end interface

677
678
679
680
681
682
683
  !> \brief abstract defintion of interface to solve the eigenvalue problem for a real valued tridiangular matrix
  !> 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
  !> \param   error       integer, optional : error code
684
  abstract interface
685
    subroutine elpa_solve_tridi_double_real_i (self, d, e, q, error)
686
687
688
689
690
691
692
693
694
      use iso_c_binding
      import elpa_t
      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
695
      integer, optional               :: error
696
697
698
    end subroutine
  end interface

699
700
701
702
703
704
705
  !> \brief abstract defintion of interface to solve the eigenvalue problem for a real valued tridiangular matrix
  !> 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
  !> \param   error       integer, optional : error code
706
  abstract interface
707
    subroutine elpa_solve_tridi_single_real_i (self, d, e, q, error)
708
709
710
711
712
713
714
715
716
      use iso_c_binding
      import elpa_t
      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
717
      integer, optional               :: error
718
719
720
    end subroutine
  end interface

721
722
723
  !> \brief abstract defintion of interface of subroutine to destroy an ELPA object
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
724
725
726
727
728
729
730
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
      class(elpa_t) :: self
    end subroutine
  end interface

731
732
733
734
  !> \brief abstract defintion of interface of function to get an integer option
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \result  value       integer
735
736
737
738
739
740
741
742
743
744
745
  abstract interface
    function elpa_get_int_i(self) result(value)
      use iso_c_binding
      import elpa_t
      class(elpa_t), intent(in) :: self
      integer(kind=C_INT) :: value
    end function
  end interface

  contains

746
747
748
749
    !> \brief function to intialise the ELPA library
    !> Parameters
    !> \param   api_version integer, api_version that ELPA should use
    !> \result  error       integer
750
751
752
753
754
755
756
757
758
759
760
761
762
763
    function elpa_init(api_version) result(error)
      use elpa_utilities, only : error_unit
      integer, intent(in) :: api_version
      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

764
765
766
    !> \brief function to check whether the ELPA library has been correctly initialised
    !> Parameters
    !> \result  state      logical
767
    function elpa_initialized() result(state)
768
769
770
771
772
773
      integer :: state
      if (initDone) then
        state = ELPA_OK
      else
        state = ELPA_ERROR
      endif
774
775
    end function

776
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
777
778
779
    subroutine elpa_uninit()
    end subroutine

780
781
782
783
    !> \brief helper function for error strings: NOT public to the user
    !> Parameters
    !> \param   elpa_error  integer
    !> \result  string      string
784
785
786
787
788
789
790
    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

791
792
793
794
    !> \brief helper function for c strings: NOT public to the user
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
795
796
797
798
799
800
801
    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

802
803
804
805
806
807
    !> \brief function to convert an integer in its string representation: NOT public to the user
    !> Parameters
    !> \param   name        string
    !> \param   value       integer
    !> \param   error       integer, optional
    !> \result  string      string
808
809
810
811
812
813
814
    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
815
      character(kind=c_char, len=elpa_int_value_to_strlen_c(name // C_NULL_CHAR, value)), pointer :: string
816
817
818
      integer(kind=c_int) :: actual_error
      type(c_ptr) :: ptr

819
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
820
821
822
823
824
825
826
827
828
829
830
831
832
833
      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

834
835
836
837
838
839
    !> \brief function to convert a string in its integer representation: NOT public to the user
    !> Parameters
    !> \param   name        string
    !> \param   string      string
    !> \param   error       integer, optional
    !> \result  value       integer
840
841
842
843
844
845
846
847
848
849
850
851
852
    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

853
      actual_error = elpa_int_string_to_value_c(name // C_NULL_CHAR, string // C_NULL_CHAR, value)
854
855
856
857
858
859
860
861
862

      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

863
864
865
866
    !> \brief function to get the cardinality of an option. NOT public to the user
    !> Parameters
    !> \param   option_name string
    !> \result  number      integer
867
868
869
870
871
872
873
    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

874
875
876
877
878
    !> \brief function to enumerate an option. NOT public to the user
    !> Parameters
    !> \param   option_name string
    !> \param   i           integer
    !> \result  option      integer
879
880
881
882
883
884
885
886
    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

887
end module