elpa_api.F90 45.8 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
  type, abstract :: elpa_t
    private

72
    !< these have to be public for proper bounds checking, sadly
73
74
75
76
77
78
79
    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
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
      ! Timer
96
97
98
      procedure(elpa_get_time_i), deferred, public :: get_time
      procedure(elpa_print_times_i), deferred, public :: print_times

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

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

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

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

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


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

134
135
136
137
      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
138

139
140
141
142
      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
143

144
145
146
147
      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
148

149
150
151
152
      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
153

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


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

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

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

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

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

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

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

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

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

316
317
318

  ! Timer routines

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

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


349
  ! Actual math routines
350

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

372
      integer, optional   :: error
373
374
375
    end subroutine
  end interface

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

397
      integer, optional   :: error
398
399
400
    end subroutine
  end interface

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

423
      integer, optional              :: error
424
425
426
    end subroutine
  end interface

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

448
      integer, optional             :: error
449
450
451
    end subroutine
  end interface

Andreas Marek's avatar
Andreas Marek committed
452
453
454
455
456
457
458
459
  !> \brief abstract definition of interface to compute C : = A**T * B for double real matrices
  !>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
  !>                 B is a (na,ncb) matrix
  !>                 C is a (na,ncb) matrix where optionally only the upper or lower
  !>                   triangle may be computed
  !>
  !> the MPI commicators are already known to the type. Thus the class method "setup" must be called
  !> BEFORE this method is used
460
  !> \details
Andreas Marek's avatar
Andreas Marek committed
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
  !>
  !> \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 global matrix A, number of rows of global matrices B and C
  !> \param ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
  !> \param nrows_a               number of rows of local (sub) matrix a
  !> \param ncols_a               number of columns of local (sub) matrix a
  !> \param b                     matrix b
  !> \param nrows_b               number of rows of local (sub) matrix b
  !> \param ncols_b               number of columns of local (sub) matrix b
  !> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
  !> \param c                     matrix c
  !> \param nrows_c               number of rows of local (sub) matrix c
  !> \param ncols_c               number of columns of local (sub) matrix c
  !> \param error                 optional argument, error code which can be queried with elpa_strerr
487
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
488
489
    subroutine elpa_hermitian_multiply_d_i (self,uplo_a, uplo_c, na, ncb, a, nrows_a, ncols_a, b, nrows_b, ncols_b, &
                                          c, nrows_c, ncols_c, error)
490
491
      use iso_c_binding
      import elpa_t
492
      implicit none
493
494
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
495
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
496
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
497
      real(kind=c_double)             :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
498
#else
Andreas Marek's avatar
Andreas Marek committed
499
      real(kind=c_double)             :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
500
#endif
501
      integer, optional               :: error
502
503
504
    end subroutine
  end interface

505
  !> \brief abstract definition of interface to compute C : = A**T * B
Andreas Marek's avatar
Andreas Marek committed
506
507
508
509
510
511
512
  !>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
  !>                 B is a (na,ncb) matrix
  !>                 C is a (na,ncb) matrix where optionally only the upper or lower
  !>                   triangle may be computed
  !>
  !> the MPI commicators are already known to the type. Thus the class method "setup" must be called
  !> BEFORE this method is used
513
  !> \details
Andreas Marek's avatar
Andreas Marek committed
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
  !>
  !> \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 global matrix A, number of rows of global matrices B and C
  !> \param ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
  !> \param nrows_a               number of rows of local (sub) matrix a
  !> \param ncols_a               number of columns of local (sub) matrix a
  !> \param b                     matrix b
  !> \param nrows_b               number of rows of local (sub) matrix b
  !> \param ncols_b               number of columns of local (sub) matrix b
  !> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
  !> \param c                     matrix c
  !> \param nrows_c               number of rows of local (sub) matrix c
  !> \param ncols_c               number of columns of local (sub) matrix c
  !> \param error                 optional argument, error code which can be queried with elpa_strerr
540
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
541
542
    subroutine elpa_hermitian_multiply_f_i (self,uplo_a, uplo_c, na, ncb, a, nrows_a, ncols_a, b, nrows_b, ncols_b, &
                                          c, nrows_c, ncols_c, error)
543
544
      use iso_c_binding
      import elpa_t
545
      implicit none
546
547
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
548
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
549
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
550
      real(kind=c_float)              :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
551
#else
Andreas Marek's avatar
Andreas Marek committed
552
      real(kind=c_float)              :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
553
#endif
554
      integer, optional               :: error
555
556
557
    end subroutine
  end interface

558
  !> \brief abstract definition of interface to compute C : = A**H * B
Andreas Marek's avatar
Andreas Marek committed
559
560
561
562
563
564
565
  !>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
  !>                 B is a (na,ncb) matrix
  !>                 C is a (na,ncb) matrix where optionally only the upper or lower
  !>                   triangle may be computed
  !>
  !> the MPI commicators are already known to the type. Thus the class method "setup" must be called
  !> BEFORE this method is used
566
  !> \details
Andreas Marek's avatar
Andreas Marek committed
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
  !>
  !> \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 global matrix A, number of rows of global matrices B and C
  !> \param ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
  !> \param nrows_a               number of rows of local (sub) matrix a
  !> \param ncols_a               number of columns of local (sub) matrix a
  !> \param b                     matrix b
  !> \param nrows_b               number of rows of local (sub) matrix b
  !> \param ncols_b               number of columns of local (sub) matrix b
  !> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
  !> \param c                     matrix c
  !> \param nrows_c               number of rows of local (sub) matrix c
  !> \param ncols_c               number of columns of local (sub) matrix c
  !> \param error                 optional argument, error code which can be queried with elpa_strerr
593
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
594
595
    subroutine elpa_hermitian_multiply_dc_i (self,uplo_a, uplo_c, na, ncb, a, nrows_a, ncols_a, b, nrows_b, ncols_b, &
                                          c, nrows_c, ncols_c, error)
596
597
      use iso_c_binding
      import elpa_t
598
      implicit none
599
600
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
601
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
602
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
603
      complex(kind=c_double_complex)  :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
604
#else
Andreas Marek's avatar
Andreas Marek committed
605
      complex(kind=c_double_complex)  :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
606
#endif
607
      integer, optional               :: error
608
609
610
    end subroutine
  end interface

611
  !> \brief abstract definition of interface to compute C : = A**H * B
Andreas Marek's avatar
Andreas Marek committed
612
613
614
615
616
617
618
  !>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
  !>                 B is a (na,ncb) matrix
  !>                 C is a (na,ncb) matrix where optionally only the upper or lower
  !>                   triangle may be computed
  !>
  !> the MPI commicators are already known to the type. Thus the class method "setup" must be called
  !> BEFORE this method is used
619
  !> \details
Andreas Marek's avatar
Andreas Marek committed
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
  !>
  !> \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 global matrix A, number of rows of global matrices B and C
  !> \param ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
  !> \param nrows_a               number of rows of local (sub) matrix a
  !> \param ncols_a               number of columns of local (sub) matrix a
  !> \param b                     matrix b
  !> \param nrows_b               number of rows of local (sub) matrix b
  !> \param ncols_b               number of columns of local (sub) matrix b
  !> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
  !> \param c                     matrix c
  !> \param nrows_c               number of rows of local (sub) matrix c
  !> \param ncols_c               number of columns of local (sub) matrix c
  !> \param error                 optional argument, error code which can be queried with elpa_strerr
646
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
647
648
    subroutine elpa_hermitian_multiply_fc_i (self, uplo_a, uplo_c, na, ncb, a, nrows_a, ncols_a, b, nrows_b, ncols_b, &
                                          c, nrows_c, ncols_c, error)
649
650
      use iso_c_binding
      import elpa_t
651
      implicit none
652
653
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
654
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
655
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
656
      complex(kind=c_float_complex)   :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
657
#else
Andreas Marek's avatar
Andreas Marek committed
658
      complex(kind=c_float_complex)   :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
659
#endif
660
      integer, optional               :: error
661
662
663
    end subroutine
  end interface

664
  !> \brief abstract definition of interface to do a cholesky decomposition of a double real matrix
665
666
667
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
668
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
669
  abstract interface
670
    subroutine elpa_cholesky_d_i (self, a, error)
671
672
      use iso_c_binding
      import elpa_t
673
      implicit none
674
675
676
677
678
679
      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
680
      integer, optional               :: error
681
682
683
    end subroutine
  end interface

684
  !> \brief abstract definition of interface to do a cholesky decomposition of a single real matrix
685
686
687
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
688
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
689
  abstract interface
690
    subroutine elpa_cholesky_f_i(self, a, error)
691
692
      use iso_c_binding
      import elpa_t
693
      implicit none
694
695
696
697
698
699
      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
700
      integer, optional               :: error
701
702
703
    end subroutine
  end interface

704
  !> \brief abstract definition of interface to do a cholesky decomposition of a double complex matrix
705
706
707
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
708
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
709
  abstract interface
710
    subroutine elpa_cholesky_dc_i (self, a, error)
711
712
      use iso_c_binding
      import elpa_t
713
      implicit none
714
715
716
717
718
719
      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
720
      integer, optional               :: error
721
722
723
    end subroutine
  end interface

724
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
725
726
727
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
728
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
729
  abstract interface
730
    subroutine elpa_cholesky_fc_i (self, a, error)
731
732
      use iso_c_binding
      import elpa_t
733
      implicit none
734
735
736
737
738
739
      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
740
      integer, optional               :: error
741
742
743
    end subroutine
  end interface

744
  !> \brief abstract definition of interface to invert a triangular double real matrix
745
746
747
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
748
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
749
  abstract interface
750
    subroutine elpa_invert_trm_d_i (self, a, error)
751
752
      use iso_c_binding
      import elpa_t
753
      implicit none
754
755
756
757
758
759
      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
760
      integer, optional               :: error
761
762
763
    end subroutine
  end interface

764
  !> \brief abstract definition of interface to invert a triangular single real matrix
765
766
767
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
768
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
769
  abstract interface
770
    subroutine elpa_invert_trm_f_i (self, a, error)
771
772
      use iso_c_binding
      import elpa_t
773
      implicit none
774
775
776
777
778
779
      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
780
      integer, optional               :: error
781
782
783
    end subroutine
  end interface

784
  !> \brief abstract definition of interface to invert a triangular double complex matrix
785
786
787
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
788
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
789
  abstract interface
790
    subroutine elpa_invert_trm_dc_i (self, a, error)
791
792
      use iso_c_binding
      import elpa_t
793
      implicit none
794
795
796
797
798
799
      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
800
      integer, optional               :: error
801
802
803
    end subroutine
  end interface

804
  !> \brief abstract definition of interface to invert a triangular single complex matrix
805
806
807
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
808
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
809
  abstract interface
810
    subroutine elpa_invert_trm_fc_i (self, a, error)
811
812
      use iso_c_binding
      import elpa_t
813
      implicit none
814
815
816
817
818
819
      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
820
      integer, optional               :: error
821
822
823
    end subroutine
  end interface

824
  !> \brief abstract definition of interface to solve the eigenvalue problem for a double-precision real valued tridiangular matrix
825
826
827
828
829
  !> 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
830
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
831
  abstract interface
832
    subroutine elpa_solve_tridi_d_i (self, d, e, q, error)
833
834
      use iso_c_binding
      import elpa_t
835
      implicit none
836
837
838
839
840
841
842
      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
843
      integer, optional               :: error
844
845
846
    end subroutine
  end interface

847
  !> \brief abstract definition of interface to solve the eigenvalue problem for a single-precision real valued tridiangular matrix
848
849
850
851
852
  !> 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
853
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
854
  abstract interface
855
    subroutine elpa_solve_tridi_f_i (self, d, e, q, error)
856
857
      use iso_c_binding
      import elpa_t
858
      implicit none
859
860
861
862
863
864
865
      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
866
      integer, optional               :: error
867
868
869
    end subroutine
  end interface

870
  !> \brief abstract definition of interface to destroy an ELPA object
871
872
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
873
874
875
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
876
      implicit none
877
878
879
880
881
882
      class(elpa_t) :: self
    end subroutine
  end interface

  contains

883
884
    !> \brief function to intialise the ELPA library
    !> Parameters
885
886
    !> \param   api_version integer: api_version that ELPA should use
    !> \result  error       integer: error code, which can be queried with elpa_strerr
887
888
889
    !
    !c> int elpa_init(int api_version);
    function elpa_init(api_version) result(error) bind(C, name="elpa_init")
890
      use elpa_utilities, only : error_unit
891
      integer, intent(in), value :: api_version
892
893
894
895
896
897
898
899
900
901
902
      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

903
904
    !> \brief function to check whether the ELPA library has been correctly initialised
    !> Parameters
905
    !> \result  state      logical: state is either ELPA_OK or ELPA_ERROR, which can be queried with elpa_strerr
906
    function elpa_initialized() result(state)
907
908
909
910
911
912
      integer :: state
      if (initDone) then
        state = ELPA_OK
      else
        state = ELPA_ERROR
      endif
913
914
    end function

915
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
916
917
918
    !
    !c> void elpa_uninit(void);
    subroutine elpa_uninit() bind(C, name="elpa_uninit")
919
920
    end subroutine

921
    !> \brief helper function for error strings
922
    !> Parameters
923
924
    !> \param   elpa_error  integer: error code to querry
    !> \result  string      string:  error string
925
926
927
928
929
930
931
    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

932
    !> \brief helper function for c strings
933
934
935
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
936
937
938
939
940
941
942
    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

943
    !> \brief function to convert an integer in its string representation
944
    !> Parameters
945
946
947
948
    !> \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
949
950
951
952
953
954
955
    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
956
      character(kind=c_char, len=elpa_int_value_to_strlen_c(name // C_NULL_CHAR, value)), pointer :: string
957
958
959
      integer(kind=c_int) :: actual_error
      type(c_ptr) :: ptr

960
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
961
962
963
964
965
966
967
968
969
970
971
972
973
974
      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

975
    !> \brief function to convert a string in its integer representation:
976
    !> Parameters
977
978
979
980
    !> \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
981
982
983
984
985
986
987
988
989
990
991
992
993
    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

994
      actual_error = elpa_int_string_to_value_c(name // C_NULL_CHAR, string // C_NULL_CHAR, value)
995
996
997
998
999
1000
1001
1002
1003

      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

1004
    !> \brief function to get the number of possible choices for an option
1005
    !> Parameters
1006
1007
    !> \param   option_name string:   the option
    !> \result  number      integer:  the total number of possible values to be chosen
1008
1009
1010
1011
1012
1013
1014
    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

1015
    !> \brief function to enumerate an option
1016
    !> Parameters
1017
    !> \param   option_name string: the option
1018
1019
    !> \param   i           integer
    !> \result  option      integer
1020
1021
1022
1023
1024
1025
1026
1027
    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

1028
end module