elpa_api.F90 50 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
358
359
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
  !>  blocksize, the number of eigenvectors
  !>  to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
  !>  It is possible to change the behaviour of the method by setting tunable parameters with the
  !>  class method "set"
360
361
362
363
364
365
  !> 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
366
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
367
  abstract interface
368
    subroutine elpa_solve_d_i(self, a, ev, q, error)
369
370
      use iso_c_binding
      import elpa_t
371
      implicit none
372
373
374
375
376
377
378
379
      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)

380
      integer, optional   :: error
381
382
383
    end subroutine
  end interface

384
  !> \brief abstract definition of interface to solve single real eigenvalue problem
385
386
387
388
389
390
391
392
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
  !>  blocksize, the number of eigenvectors
  !>  to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
  !>  It is possible to change the behaviour of the method by setting tunable parameters with the
  !>  class method "set"
393
394
395
396
397
398
  !> 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
399
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
400
  abstract interface
401
    subroutine elpa_solve_f_i(self, a, ev, q, error)
402
403
      use iso_c_binding
      import elpa_t
404
      implicit none
405
406
407
408
409
410
411
412
      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)

413
      integer, optional   :: error
414
415
416
    end subroutine
  end interface

417
  !> \brief abstract definition of interface to solve double complex eigenvalue problem
418
419
420
421
422
423
424
425
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
  !>  blocksize, the number of eigenvectors
  !>  to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
  !>  It is possible to change the behaviour of the method by setting tunable parameters with the
  !>  class method "set"
426
427
428
429
430
431
  !> 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
432
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
433
  abstract interface
434
    subroutine elpa_solve_dc_i(self, a, ev, q, error)
435
436
      use iso_c_binding
      import elpa_t
437
      implicit none
438
439
440
441
442
443
444
445
446
      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)

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

451
  !> \brief abstract definition of interface to solve single complex eigenvalue problem
452
453
454
455
456
457
458
459
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
  !>  blocksize, the number of eigenvectors
  !>  to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
  !>  It is possible to change the behaviour of the method by setting tunable parameters with the
  !>  class method "set"
460
461
462
463
464
465
  !> 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
466
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
467
  abstract interface
468
    subroutine elpa_solve_fc_i(self, a, ev, q, error)
469
470
      use iso_c_binding
      import elpa_t
471
      implicit none
472
473
474
475
476
477
478
479
      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)

480
      integer, optional             :: error
481
482
483
    end subroutine
  end interface

Andreas Marek's avatar
Andreas Marek committed
484
485
486
487
488
489
490
491
  !> \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
492
  !> \details
Andreas Marek's avatar
Andreas Marek committed
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
  !>
  !> \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
519
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
520
521
    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)
522
523
      use iso_c_binding
      import elpa_t
524
      implicit none
525
526
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
527
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
528
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
529
      real(kind=c_double)             :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
530
#else
Andreas Marek's avatar
Andreas Marek committed
531
      real(kind=c_double)             :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
532
#endif
533
      integer, optional               :: error
534
535
536
    end subroutine
  end interface

537
  !> \brief abstract definition of interface to compute C : = A**T * B
Andreas Marek's avatar
Andreas Marek committed
538
539
540
541
542
543
544
  !>         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
545
  !> \details
Andreas Marek's avatar
Andreas Marek committed
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
  !>
  !> \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
572
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
573
574
    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)
575
576
      use iso_c_binding
      import elpa_t
577
      implicit none
578
579
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
580
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
581
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
582
      real(kind=c_float)              :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
583
#else
Andreas Marek's avatar
Andreas Marek committed
584
      real(kind=c_float)              :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
585
#endif
586
      integer, optional               :: error
587
588
589
    end subroutine
  end interface

590
  !> \brief abstract definition of interface to compute C : = A**H * B
Andreas Marek's avatar
Andreas Marek committed
591
592
593
594
595
596
597
  !>         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
598
  !> \details
Andreas Marek's avatar
Andreas Marek committed
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
  !>
  !> \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
625
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
626
627
    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)
628
629
      use iso_c_binding
      import elpa_t
630
      implicit none
631
632
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
633
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
634
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
635
      complex(kind=c_double_complex)  :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
636
#else
Andreas Marek's avatar
Andreas Marek committed
637
      complex(kind=c_double_complex)  :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
638
#endif
639
      integer, optional               :: error
640
641
642
    end subroutine
  end interface

643
  !> \brief abstract definition of interface to compute C : = A**H * B
Andreas Marek's avatar
Andreas Marek committed
644
645
646
647
648
649
650
  !>         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
651
  !> \details
Andreas Marek's avatar
Andreas Marek committed
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
  !>
  !> \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
678
  abstract interface
Andreas Marek's avatar
Andreas Marek committed
679
680
    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)
681
682
      use iso_c_binding
      import elpa_t
683
      implicit none
684
685
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
Andreas Marek's avatar
Andreas Marek committed
686
      integer(kind=c_int), intent(in) :: na, nrows_a, ncols_a, nrows_b, ncols_b, nrows_c, ncols_c, ncb
687
#ifdef USE_ASSUMED_SIZE
Andreas Marek's avatar
Andreas Marek committed
688
      complex(kind=c_float_complex)   :: a(nrows_a,*), b(nrows_b,*), c(nrows_c,*)
689
#else
Andreas Marek's avatar
Andreas Marek committed
690
      complex(kind=c_float_complex)   :: a(nrows_a,ncols_a), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
691
#endif
692
      integer, optional               :: error
693
694
695
    end subroutine
  end interface

696
  !> \brief abstract definition of interface to do a cholesky decomposition of a double real matrix
697
698
699
700
701
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !> 
702
703
704
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
705
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
706
  abstract interface
707
    subroutine elpa_cholesky_d_i (self, a, error)
708
709
      use iso_c_binding
      import elpa_t
710
      implicit none
711
712
713
714
715
716
      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
717
      integer, optional               :: error
718
719
720
    end subroutine
  end interface

721
  !> \brief abstract definition of interface to do a cholesky decomposition of a single real matrix
722
723
724
725
726
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !> 
727
728
729
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
730
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
731
  abstract interface
732
    subroutine elpa_cholesky_f_i(self, a, error)
733
734
      use iso_c_binding
      import elpa_t
735
      implicit none
736
737
738
739
740
741
      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
742
      integer, optional               :: error
743
744
745
    end subroutine
  end interface

746
  !> \brief abstract definition of interface to do a cholesky decomposition of a double complex matrix
747
748
749
750
751
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !> 
752
753
754
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
755
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
756
  abstract interface
757
    subroutine elpa_cholesky_dc_i (self, a, error)
758
759
      use iso_c_binding
      import elpa_t
760
      implicit none
761
762
763
764
765
766
      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
767
      integer, optional               :: error
768
769
770
    end subroutine
  end interface

771
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
772
773
774
775
776
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !> 
777
778
779
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
780
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
781
  abstract interface
782
    subroutine elpa_cholesky_fc_i (self, a, error)
783
784
      use iso_c_binding
      import elpa_t
785
      implicit none
786
787
788
789
790
791
      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
792
      integer, optional               :: error
793
794
795
    end subroutine
  end interface

796
  !> \brief abstract definition of interface to invert a triangular double real matrix
797
798
799
800
801
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
802
803
804
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
805
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
806
  abstract interface
807
    subroutine elpa_invert_trm_d_i (self, a, error)
808
809
      use iso_c_binding
      import elpa_t
810
      implicit none
811
812
813
814
815
816
      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
817
      integer, optional               :: error
818
819
820
    end subroutine
  end interface

821
  !> \brief abstract definition of interface to invert a triangular single real matrix
822
  !> Parameters
823
824
825
826
827
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
828
829
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
830
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
831
  abstract interface
832
    subroutine elpa_invert_trm_f_i (self, a, error)
833
834
      use iso_c_binding
      import elpa_t
835
      implicit none
836
837
838
839
840
841
      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
842
      integer, optional               :: error
843
844
845
    end subroutine
  end interface

846
  !> \brief abstract definition of interface to invert a triangular double complex matrix
847
848
849
850
851
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
852
853
854
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
855
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
856
  abstract interface
857
    subroutine elpa_invert_trm_dc_i (self, a, error)
858
859
      use iso_c_binding
      import elpa_t
860
      implicit none
861
862
863
864
865
866
      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
867
      integer, optional               :: error
868
869
870
    end subroutine
  end interface

871
  !> \brief abstract definition of interface to invert a triangular single complex matrix
872
873
874
875
876
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
877
878
879
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
880
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
881
  abstract interface
882
    subroutine elpa_invert_trm_fc_i (self, a, error)
883
884
      use iso_c_binding
      import elpa_t
885
      implicit none
886
887
888
889
890
891
      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
892
      integer, optional               :: error
893
894
895
    end subroutine
  end interface

896
  !> \brief abstract definition of interface to solve the eigenvalue problem for a double-precision real valued tridiangular matrix
897
898
899
900
901
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
902
903
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
904
905
  !> \param   d           double real 1d array: the diagonal elements of a matrix defined in setup, on output the eigenvalues
  !>                      in ascending order
906
907
  !> \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
908
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
909
  abstract interface
910
    subroutine elpa_solve_tridi_d_i (self, d, e, q, error)
911
912
      use iso_c_binding
      import elpa_t
913
      implicit none
914
915
916
917
918
919
920
      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
921
      integer, optional               :: error
922
923
924
    end subroutine
  end interface

925
  !> \brief abstract definition of interface to solve the eigenvalue problem for a single-precision real valued tridiangular matrix
926
927
928
929
930
  !>
  !>  The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
  !>  block size, and the MPI communicators are already known to the object and MUST be set BEFORE
  !>  with the class method "setup"
  !>
931
932
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
933
934
  !> \param   d           single real 1d array: the diagonal elements of a matrix defined in setup, on output the eigenvalues
  !>                      in ascending order
935
936
  !> \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
937
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
938
  abstract interface
939
    subroutine elpa_solve_tridi_f_i (self, d, e, q, error)
940
941
      use iso_c_binding
      import elpa_t
942
      implicit none
943
944
945
946
947
948
949
      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
950
      integer, optional               :: error
951
952
953
    end subroutine
  end interface

954
  !> \brief abstract definition of interface to destroy an ELPA object
955
956
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
957
958
959
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
960
      implicit none
961
962
963
964
965
966
      class(elpa_t) :: self
    end subroutine
  end interface

  contains

967
968
    !> \brief function to intialise the ELPA library
    !> Parameters
969
970
    !> \param   api_version integer: api_version that ELPA should use
    !> \result  error       integer: error code, which can be queried with elpa_strerr
971
972
973
    !
    !c> int elpa_init(int api_version);
    function elpa_init(api_version) result(error) bind(C, name="elpa_init")
974
      use elpa_utilities, only : error_unit
975
      integer, intent(in), value :: api_version
976
977
978
979
980
981
982
983
984
985
986
      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

987
988
    !> \brief function to check whether the ELPA library has been correctly initialised
    !> Parameters
989
    !> \result  state      logical: state is either ELPA_OK or ELPA_ERROR, which can be queried with elpa_strerr
990
    function elpa_initialized() result(state)
991
992
993
994
995
996
      integer :: state
      if (initDone) then
        state = ELPA_OK
      else
        state = ELPA_ERROR
      endif
997
998
    end function

999
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
1000
1001
1002
    !
    !c> void elpa_uninit(void);
    subroutine elpa_uninit() bind(C, name="elpa_uninit")
1003
1004
    end subroutine

1005
    !> \brief helper function for error strings
1006
    !> Parameters
1007
1008
    !> \param   elpa_error  integer: error code to querry
    !> \result  string      string:  error string
1009
1010
1011
1012
1013
1014
1015
    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

1016
    !> \brief helper function for c strings
1017
1018
1019
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
1020
1021
1022
1023
1024
1025
1026
    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

1027
    !> \brief function to convert an integer in its string representation
1028
    !> Parameters
1029
1030
1031
1032
    !> \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
1033
1034
1035
1036
1037
1038
1039
    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
1040
      character(kind=c_char, len=elpa_int_value_to_strlen_c(name // C_NULL_CHAR, value)), pointer :: string
1041
1042
1043
      integer(kind=c_int) :: actual_error
      type(c_ptr) :: ptr

1044
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
      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

1059
    !> \brief function to convert a string in its integer representation:
1060
    !> Parameters
1061
1062
1063
1064
    !> \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
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075
1076
1077
    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),