elpa_api.F90 50.1 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
!
!    Copyright 2017, L. Hüdepohl and A. Marek, MPCDF
!
!    This file is part of ELPA.
!
!    The ELPA library was originally created by the ELPA consortium,
!    consisting of the following organizations:
!
!    - Max Planck Computing and Data Facility (MPCDF), formerly known as
!      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!    - Bergische Universität Wuppertal, Lehrstuhl für angewandte
!      Informatik,
!    - Technische Universität München, Lehrstuhl für Informatik mit
!      Schwerpunkt Wissenschaftliches Rechnen ,
!    - Fritz-Haber-Institut, Berlin, Abt. Theorie,
!    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
!      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!      and
!    - IBM Deutschland GmbH
!
!    This particular source code file contains additions, changes and
!    enhancements authored by Intel Corporation which is not part of
!    the ELPA consortium.
!
!    More information can be found here:
!    http://elpa.mpcdf.mpg.de/
!
!    ELPA is free software: you can redistribute it and/or modify
!    it under the terms of the version 3 of the license of the
!    GNU Lesser General Public License as published by the Free
!    Software Foundation.
!
!    ELPA is distributed in the hope that it will be useful,
!    but WITHOUT ANY WARRANTY; without even the implied warranty of
!    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!    GNU Lesser General Public License for more details.
!
!    You should have received a copy of the GNU Lesser General Public License
!    along with ELPA.  If not, see <http://www.gnu.org/licenses/>
!
!    ELPA reflects a substantial effort on the part of the original
!    ELPA consortium, and we ask you to respect the spirit of the
!    license that we chose: i.e., please contribute any changes you
!    may have back to the original ELPA library distribution, and keep
!    any derivatives of ELPA under the same license that we chose for
!    the original distribution, the GNU Lesser General Public License.
!
#include "config-f90.h"
49
50
51
!> \brief Fortran module which provides the definition of the ELPA API. Do not use directly! Use the module "elpa"


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

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

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

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

68
  !> \brief Abstract definition of the elpa_t type
69
70
71
  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
91
92

      generic, public :: get => &                                   !< export a method to get integer/double key/values
          elpa_get_integer, &
          elpa_get_double
93

94
95
      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
96

97
      ! Timer
98
99
100
      procedure(elpa_get_time_i), deferred, public :: get_time
      procedure(elpa_print_times_i), deferred, public :: print_times

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

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

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

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

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


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

136
137
138
      procedure(elpa_get_integer_i), deferred, private :: elpa_get_integer
      procedure(elpa_get_double_i),  deferred, private :: elpa_get_double

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

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

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

154
155
156
157
      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
158

159
160
      procedure(elpa_solve_tridi_d_i), deferred, private :: elpa_solve_tridi_d
      procedure(elpa_solve_tridi_f_i), deferred, private :: elpa_solve_tridi_f
161
162
163
  end type elpa_t


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

178
179
180
181
182
  !> \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()
183
  abstract interface
184
    function elpa_setup_i(self) result(error)
185
      import elpa_t
186
      implicit none
187
      class(elpa_t), intent(inout) :: self
188
      integer :: error
189
190
191
    end function
  end interface

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

211
212
213
214
215
  !> \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
216
  !> \param   value       integer : the value corresponding to the key
217
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr()
218
  abstract interface
219
    subroutine elpa_get_integer_i(self, name, value, error)
220
221
      use iso_c_binding
      import elpa_t
222
      implicit none
223
224
225
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
226
      integer, intent(out), optional :: error
227
    end subroutine
228
229
  end interface

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

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

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

285
286
287
288
289
  !> \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
290
  !> \param   value       double: the value associated with the key
291
  !> \param   error       integer. optional : error code, which can be queried with elpa_strerr
292
  abstract interface
293
    subroutine elpa_get_double_i(self, name, value, error)
294
295
      use iso_c_binding
      import elpa_t
296
      implicit none
297
298
299
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      real(kind=c_double)            :: value
300
      integer, intent(out), optional :: error
301
    end subroutine
302
303
  end interface

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

321
322
323

  ! Timer routines

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

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


354
  ! Actual math routines
355

356
  !> \brief abstract definition of interface to solve double real eigenvalue problem
357
358
359
360
361
362
363
364
  !>
  !>  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"
365
366
367
368
369
370
  !> 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
371
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
372
  abstract interface
373
    subroutine elpa_solve_d_i(self, a, ev, q, error)
374
375
      use iso_c_binding
      import elpa_t
376
      implicit none
377
378
379
380
381
382
383
384
      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)

385
      integer, optional   :: error
386
387
388
    end subroutine
  end interface

389
  !> \brief abstract definition of interface to solve single real eigenvalue problem
390
391
392
393
394
395
396
397
  !>
  !>  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"
398
399
400
401
402
403
  !> 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
404
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
405
  abstract interface
406
    subroutine elpa_solve_f_i(self, a, ev, q, error)
407
408
      use iso_c_binding
      import elpa_t
409
      implicit none
410
411
412
413
414
415
416
417
      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)

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

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

452
      integer, optional              :: error
453
454
455
    end subroutine
  end interface

456
  !> \brief abstract definition of interface to solve single complex eigenvalue problem
457
458
459
460
461
462
463
464
  !>
  !>  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"
465
466
467
468
469
470
  !> 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
471
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
472
  abstract interface
473
    subroutine elpa_solve_fc_i(self, a, ev, q, error)
474
475
      use iso_c_binding
      import elpa_t
476
      implicit none
477
478
479
480
481
482
483
484
      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)

485
      integer, optional             :: error
486
487
488
    end subroutine
  end interface

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

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

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

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

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

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

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

776
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
777
778
779
780
781
  !>
  !>  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"
  !> 
782
783
784
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
785
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
786
  abstract interface
787
    subroutine elpa_cholesky_fc_i (self, a, error)
788
789
      use iso_c_binding
      import elpa_t
790
      implicit none
791
792
793
794
795
796
      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
797
      integer, optional               :: error
798
799
800
    end subroutine
  end interface

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

826
  !> \brief abstract definition of interface to invert a triangular single real matrix
827
  !> Parameters
828
829
830
831
832
  !>
  !>  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"
  !>
833
834
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
835
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
836
  abstract interface
837
    subroutine elpa_invert_trm_f_i (self, a, error)
838
839
      use iso_c_binding
      import elpa_t
840
      implicit none
841
842
843
844
845
846
      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
847
      integer, optional               :: error
848
849
850
    end subroutine
  end interface

851
  !> \brief abstract definition of interface to invert a triangular double complex matrix
852
853
854
855
856
  !>
  !>  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"
  !>
857
858
859
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be inverted
860
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
861
  abstract interface
862
    subroutine elpa_invert_trm_dc_i (self, a, error)
863
864
      use iso_c_binding
      import elpa_t
865
      implicit none
866
867
868
869
870
871
      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
872
      integer, optional               :: error
873
874
875
    end subroutine
  end interface

876
  !> \brief abstract definition of interface to invert a triangular single complex matrix
877
878
879
880
881
  !>
  !>  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"
  !>
882
883
884
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be inverted
885
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
886
  abstract interface
887
    subroutine elpa_invert_trm_fc_i (self, a, error)
888
889
      use iso_c_binding
      import elpa_t
890
      implicit none
891
892
893
894
895
896
      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
897
      integer, optional               :: error
898
899
900
    end subroutine
  end interface

901
  !> \brief abstract definition of interface to solve the eigenvalue problem for a double-precision real valued tridiangular matrix
902
903
904
905
906
  !>
  !>  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"
  !>
907
908
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
909
910
  !> \param   d           double real 1d array: the diagonal elements of a matrix defined in setup, on output the eigenvalues
  !>                      in ascending order
911
912
  !> \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
913
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
914
  abstract interface
915
    subroutine elpa_solve_tridi_d_i (self, d, e, q, error)
916
917
      use iso_c_binding
      import elpa_t
918
      implicit none
919
920
921
922
923
924
925
      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
926
      integer, optional               :: error
927
928
929
    end subroutine
  end interface

930
  !> \brief abstract definition of interface to solve the eigenvalue problem for a single-precision real valued tridiangular matrix
931
932
933
934
935
  !>
  !>  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"
  !>
936
937
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
938
939
  !> \param   d           single real 1d array: the diagonal elements of a matrix defined in setup, on output the eigenvalues
  !>                      in ascending order
940
941
  !> \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
942
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
943
  abstract interface
944
    subroutine elpa_solve_tridi_f_i (self, d, e, q, error)
945
946
      use iso_c_binding
      import elpa_t
947
      implicit none
948
949
950
951
952
953
954
      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
955
      integer, optional               :: error
956
957
958
    end subroutine
  end interface

959
  !> \brief abstract definition of interface to destroy an ELPA object
960
961
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
962
963
964
  abstract interface
    subroutine elpa_destroy_i(self)
      import elpa_t
965
      implicit none
966
967
968
969
970
971
      class(elpa_t) :: self
    end subroutine
  end interface

  contains

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

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

1004
    !> \brief subroutine to uninit the ELPA library. Does nothing at the moment. Might do sth. later
1005
1006
1007
    !
    !c> void elpa_uninit(void);
    subroutine elpa_uninit() bind(C, name="elpa_uninit")
1008
1009
    end subroutine

1010
    !> \brief helper function for error strings
1011
    !> Parameters
1012
1013
    !> \param   elpa_error  integer: error code to querry
    !> \result  string      string:  error string
1014
1015
1016
1017
1018
1019
1020
    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

1021
    !> \brief helper function for c strings
1022
1023
1024
    !> Parameters
    !> \param   ptr         type(c_ptr)
    !> \result  string      string
1025
1026
1027
1028
1029
1030
1031
    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

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

1049
      actual_error = elpa_int_value_to_string_c(name // C_NULL_CHAR, value, ptr)
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
      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