elpa_api.F90 56.3 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 80
    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
      ! general
81 82
      procedure(elpa_setup_i),   deferred, public :: setup          !< method to setup an ELPA object
      procedure(elpa_destroy_i), deferred, public :: destroy        !< method to destroy an ELPA object
83

84
      ! key/value store
85
      generic, public :: set => &                                   !< export a method to set integer/double key/values
86 87
          elpa_set_integer, &
          elpa_set_double
88 89 90 91

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

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

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

100
      ! Actual math routines
Andreas Marek's avatar
Andreas Marek committed
101 102 103
      generic, public :: eigenvectors => &                          !< method eigenvectors for solving the full eigenvalue problem
          elpa_eigenvectors_d, &                                    !< the eigenvalues and (parts of) the eigenvectors are computed
          elpa_eigenvectors_f, &                                    !< for symmetric real valued / hermitian complex valued matrices
104 105
          elpa_eigenvectors_dc, &
          elpa_eigenvectors_fc
106

Andreas Marek's avatar
Andreas Marek committed
107 108 109 110 111 112
      generic, public :: eigenvalues => &                           !< method eigenvalues for solving the eigenvalue problem
          elpa_eigenvalues_d, &                                     !< only the eigenvalues are computed
          elpa_eigenvalues_f, &                                     !< for symmetric real valued / hermitian complex valued matrices
          elpa_eigenvalues_dc, &
          elpa_eigenvalues_fc

113
      generic, public :: hermitian_multiply => &                    !< method for a "hermitian" multiplication of matrices a and b
114
          elpa_hermitian_multiply_d, &                              !< for real valued matrices:   a**T * b
115
          elpa_hermitian_multiply_dc, &                             !< for complex valued matrices a**H * b
116 117
          elpa_hermitian_multiply_f, &
          elpa_hermitian_multiply_fc
118

119
      generic, public :: cholesky => &                              !< method for the cholesky factorisation of matrix a
120 121 122 123
          elpa_cholesky_d, &
          elpa_cholesky_f, &
          elpa_cholesky_dc, &
          elpa_cholesky_fc
124

Andreas Marek's avatar
Andreas Marek committed
125
      generic, public :: invert_triangular => &                     !< method to invert a upper triangular matrix a
126 127 128 129
          elpa_invert_trm_d, &
          elpa_invert_trm_f, &
          elpa_invert_trm_dc, &
          elpa_invert_trm_fc
130

131
      generic, private :: solve_tridi => &                          !< method to solve the eigenvalue problem for a tridiagonal
132
          elpa_solve_tridi_d, &                                     !< matrix
133
          elpa_solve_tridi_f
134 135


136
      !> \brief private methods of elpa_t type. NOT accessible for the user
137 138 139 140
      ! privates
      procedure(elpa_set_integer_i), deferred, private :: elpa_set_integer
      procedure(elpa_set_double_i),  deferred, private :: elpa_set_double

141 142 143
      procedure(elpa_get_integer_i), deferred, private :: elpa_get_integer
      procedure(elpa_get_double_i),  deferred, private :: elpa_get_double

144 145 146 147
      procedure(elpa_eigenvectors_d_i),    deferred, private :: elpa_eigenvectors_d
      procedure(elpa_eigenvectors_f_i),    deferred, private :: elpa_eigenvectors_f
      procedure(elpa_eigenvectors_dc_i), deferred, private :: elpa_eigenvectors_dc
      procedure(elpa_eigenvectors_fc_i), deferred, private :: elpa_eigenvectors_fc
148

Andreas Marek's avatar
Andreas Marek committed
149 150 151 152 153
      procedure(elpa_eigenvalues_d_i),    deferred, private :: elpa_eigenvalues_d
      procedure(elpa_eigenvalues_f_i),    deferred, private :: elpa_eigenvalues_f
      procedure(elpa_eigenvalues_dc_i), deferred, private :: elpa_eigenvalues_dc
      procedure(elpa_eigenvalues_fc_i), deferred, private :: elpa_eigenvalues_fc

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

159 160 161 162
      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
163

164 165 166 167
      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
168

169 170
      procedure(elpa_solve_tridi_d_i), deferred, private :: elpa_solve_tridi_d
      procedure(elpa_solve_tridi_f_i), deferred, private :: elpa_solve_tridi_f
171 172 173
  end type elpa_t


174 175 176 177 178
  !> \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
179 180 181
  interface
    pure function elpa_strlen_c(ptr) result(size) bind(c, name="strlen")
      use, intrinsic :: iso_c_binding
182
      implicit none
183 184 185 186 187
      type(c_ptr), intent(in), value :: ptr
      integer(kind=c_size_t) :: size
    end function
  end interface

188 189 190 191 192
  !> \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()
193
  abstract interface
194
    function elpa_setup_i(self) result(error)
195
      import elpa_t
196
      implicit none
197
      class(elpa_t), intent(inout) :: self
198
      integer :: error
199 200 201
    end function
  end interface

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

221 222 223 224 225
  !> \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
226
  !> \param   value       integer : the value corresponding to the key
227
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr()
228
  abstract interface
229
    subroutine elpa_get_integer_i(self, name, value, error)
230 231
      use iso_c_binding
      import elpa_t
232
      implicit none
233 234 235
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int)            :: value
236
      integer, intent(out), optional :: error
237
    end subroutine
238 239
  end interface

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

257 258 259 260 261 262 263 264
  !> \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
265
  abstract interface
266
    function elpa_can_set_i(self, name, value) result(state)
267
      import elpa_t, c_int
268
      implicit none
269 270 271
      class(elpa_t)                   :: self
      character(*), intent(in)        :: name
      integer(kind=c_int), intent(in) :: value
272
      integer                         :: state
273
    end function
274 275
  end interface

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

295 296 297 298 299
  !> \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
300
  !> \param   value       double: the value associated with the key
301
  !> \param   error       integer. optional : error code, which can be queried with elpa_strerr
302
  abstract interface
303
    subroutine elpa_get_double_i(self, name, value, error)
304 305
      use iso_c_binding
      import elpa_t
306
      implicit none
307 308 309
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      real(kind=c_double)            :: value
310
      integer, intent(out), optional :: error
311
    end subroutine
312 313
  end interface

314 315 316 317 318 319
  !> \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
320 321 322 323
  abstract interface
    function elpa_associate_int_i(self, name) result(value)
      use iso_c_binding
      import elpa_t
324
      implicit none
325 326 327 328 329 330
      class(elpa_t)                  :: self
      character(*), intent(in)       :: name
      integer(kind=c_int), pointer   :: value
    end function
  end interface

331 332 333

  ! Timer routines

334 335 336 337 338 339
  !> \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
340 341 342
  abstract interface
    function elpa_get_time_i(self, name1, name2, name3, name4, name5, name6) result(s)
      import elpa_t, c_double
343
      implicit none
344 345 346 347 348 349 350
      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

351 352 353 354
  !> \brief abstract definition of print method for timer
  !> Parameters
  !> \details
  !> \param   self        class(elpa_t): the ELPA object
355 356 357
  abstract interface
    subroutine elpa_print_times_i(self)
      import elpa_t
358
      implicit none
359 360 361 362 363
      class(elpa_t), intent(in) :: self
    end subroutine
  end interface


364
  ! Actual math routines
365

366
  !> \brief abstract definition of interface to solve double real eigenvalue problem
367 368 369 370 371 372 373 374
  !>
  !>  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"
375 376 377 378 379 380
  !> 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
381
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
382
  abstract interface
383
    subroutine elpa_eigenvectors_d_i(self, a, ev, q, error)
384 385
      use iso_c_binding
      import elpa_t
386
      implicit none
387 388 389 390 391 392 393 394
      class(elpa_t)       :: self
#ifdef USE_ASSUMED_SIZE
      real(kind=c_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)

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

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

428
      integer, optional   :: error
429 430 431
    end subroutine
  end interface

432
  !> \brief abstract definition of interface to solve double complex eigenvalue problem
433 434 435 436 437 438 439 440
  !>
  !>  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"
441 442 443 444 445 446
  !> 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
447
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
448
  abstract interface
449
    subroutine elpa_eigenvectors_dc_i(self, a, ev, q, error)
450 451
      use iso_c_binding
      import elpa_t
452
      implicit none
453 454 455 456 457 458 459 460 461
      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)

462
      integer, optional              :: error
463 464 465
    end subroutine
  end interface

466
  !> \brief abstract definition of interface to solve single complex eigenvalue problem
467 468 469 470 471 472 473 474
  !>
  !>  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"
475 476 477 478 479 480
  !> 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
481
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
482
  abstract interface
483
    subroutine elpa_eigenvectors_fc_i(self, a, ev, q, error)
484 485
      use iso_c_binding
      import elpa_t
486
      implicit none
487 488 489 490 491 492 493 494
      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)

495
      integer, optional             :: error
496 497 498
    end subroutine
  end interface

Andreas Marek's avatar
Andreas Marek committed
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 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630



  !> \brief abstract definition of interface to solve double real eigenvalue problem
  !>
  !>  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"
  !> 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
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
  abstract interface
    subroutine elpa_eigenvalues_d_i(self, a, ev, error)
      use iso_c_binding
      import elpa_t
      implicit none
      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
      real(kind=c_double) :: ev(self%na)

      integer, optional   :: error
    end subroutine
  end interface

  !> \brief abstract definition of interface to solve single real eigenvalue problem
  !>
  !>  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"
  !> 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
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
  abstract interface
    subroutine elpa_eigenvalues_f_i(self, a, ev, error)
      use iso_c_binding
      import elpa_t
      implicit none
      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
      real(kind=c_float)  :: ev(self%na)

      integer, optional   :: error
    end subroutine
  end interface

  !> \brief abstract definition of interface to solve double complex eigenvalue problem
  !>
  !>  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"
  !> 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
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
  abstract interface
    subroutine elpa_eigenvalues_dc_i(self, a, ev, error)
      use iso_c_binding
      import elpa_t
      implicit none
      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
      real(kind=c_double)            :: ev(self%na)

      integer, optional              :: error
    end subroutine
  end interface

  !> \brief abstract definition of interface to solve single complex eigenvalue problem
  !>
  !>  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"
  !> 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
  !> \result  error       integer, optional : error code, which can be queried with elpa_strerr
  abstract interface
    subroutine elpa_eigenvalues_fc_i(self, a, ev, error)
      use iso_c_binding
      import elpa_t
      implicit none
      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
      real(kind=c_float)            :: ev(self%na)

      integer, optional             :: error
    end subroutine
  end interface

Andreas Marek's avatar
Andreas Marek committed
631
  !> \brief abstract definition of interface to compute C : = A**T * B for double real matrices
632 633 634
  !>         where   A is a square matrix (self%a,self%na) which is optionally upper or lower triangular
  !>                 B is a (self%na,ncb) matrix
  !>                 C is a (self%na,ncb) matrix where optionally only the upper or lower
Andreas Marek's avatar
Andreas Marek committed
635 636 637 638
  !>                   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
639
  !> \details
Andreas Marek's avatar
Andreas Marek committed
640
  !>
641
  !> \param   self                class(elpa_t), the ELPA object
Andreas Marek's avatar
Andreas Marek committed
642 643 644 645 646 647 648 649 650 651 652 653 654 655
  !> \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 ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
656 657
  !> \param self%local_nrows      number of rows of local (sub) matrix a, set with method set("local_nrows,value")
  !> \param self%local_ncols      number of columns of local (sub) matrix a, set with method set("local_ncols,value")
Andreas Marek's avatar
Andreas Marek committed
658 659 660 661 662 663 664 665
  !> \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
666
  abstract interface
667
    subroutine elpa_hermitian_multiply_d_i (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
Andreas Marek's avatar
Andreas Marek committed
668
                                          c, nrows_c, ncols_c, error)
669 670
      use iso_c_binding
      import elpa_t
671
      implicit none
672 673
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
674
      integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
675
#ifdef USE_ASSUMED_SIZE
676
      real(kind=c_double)             :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
677
#else
678
      real(kind=c_double)             :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
679
#endif
680
      integer, optional               :: error
681 682 683
    end subroutine
  end interface

684
  !> \brief abstract definition of interface to compute C : = A**T * B
685 686 687
  !>         where   A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
  !>                 B is a (self%na,ncb) matrix
  !>                 C is a (self%na,ncb) matrix where optionally only the upper or lower
Andreas Marek's avatar
Andreas Marek committed
688 689 690 691
  !>                   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
692
  !> \details
Andreas Marek's avatar
Andreas Marek committed
693
  !>
694
  !> \param   self                class(elpa_t), the ELPA object
Andreas Marek's avatar
Andreas Marek committed
695 696 697 698 699 700 701 702 703 704 705 706 707 708
  !> \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 ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
709 710
  !> \param self%local_nrows      number of rows of local (sub) matrix a, set with method set("local_nrows",value)
  !> \param self%local_ncols      number of columns of local (sub) matrix a, set with method set("local_nrows",value)
Andreas Marek's avatar
Andreas Marek committed
711 712 713 714 715 716 717 718
  !> \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
719
  abstract interface
720
    subroutine elpa_hermitian_multiply_f_i (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
Andreas Marek's avatar
Andreas Marek committed
721
                                          c, nrows_c, ncols_c, error)
722 723
      use iso_c_binding
      import elpa_t
724
      implicit none
725 726
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
727
      integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
728
#ifdef USE_ASSUMED_SIZE
729
      real(kind=c_float)              :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
730
#else
731
      real(kind=c_float)              :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
732
#endif
733
      integer, optional               :: error
734 735 736
    end subroutine
  end interface

737
  !> \brief abstract definition of interface to compute C : = A**H * B
738 739 740
  !>         where   A is a square matrix (self%na,self%a) which is optionally upper or lower triangular
  !>                 B is a (self%na,ncb) matrix
  !>                 C is a (self%na,ncb) matrix where optionally only the upper or lower
Andreas Marek's avatar
Andreas Marek committed
741 742 743 744
  !>                   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
745
  !> \details
Andreas Marek's avatar
Andreas Marek committed
746
  !>
747
  !> \param   self                class(elpa_t), the ELPA object
Andreas Marek's avatar
Andreas Marek committed
748 749 750 751 752 753 754 755 756 757 758 759 760 761
  !> \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 ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
762 763
  !> \param self%local_nrows      number of rows of local (sub) matrix a, set with the method set("local_nrows",value)
  !> \param self%local_ncols      number of columns of local (sub) matrix a, set with the method set("local_ncols",value)
Andreas Marek's avatar
Andreas Marek committed
764 765 766 767 768 769 770 771
  !> \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
772
  abstract interface
773
    subroutine elpa_hermitian_multiply_dc_i (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
Andreas Marek's avatar
Andreas Marek committed
774
                                          c, nrows_c, ncols_c, error)
775 776
      use iso_c_binding
      import elpa_t
777
      implicit none
778 779
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
780
      integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
781
#ifdef USE_ASSUMED_SIZE
782
      complex(kind=c_double_complex)  :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
783
#else
784
      complex(kind=c_double_complex)  :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
785
#endif
786
      integer, optional               :: error
787 788 789
    end subroutine
  end interface

790
  !> \brief abstract definition of interface to compute C : = A**H * B
791 792 793
  !>         where   A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
  !>                 B is a (self%na,ncb) matrix
  !>                 C is a (self%na,ncb) matrix where optionally only the upper or lower
Andreas Marek's avatar
Andreas Marek committed
794 795 796 797
  !>                   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
798
  !> \details
Andreas Marek's avatar
Andreas Marek committed
799
  !>
800
  !> \param   self                class(elpa_t), the ELPA object
Andreas Marek's avatar
Andreas Marek committed
801 802 803 804 805 806 807 808 809 810 811 812 813 814
  !> \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 ncb                   Number of columns  of global matrices B and C
  !> \param a                     matrix a
815 816
  !> \param self%local_nrows      number of rows of local (sub) matrix a, set with class method set("local_nrows",value)
  !> \param self%local_ncols      number of columns of local (sub) matrix a, set with class method set("local_ncols",value)
Andreas Marek's avatar
Andreas Marek committed
817 818 819 820 821 822 823 824
  !> \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
825
  abstract interface
826
    subroutine elpa_hermitian_multiply_fc_i (self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
Andreas Marek's avatar
Andreas Marek committed
827
                                          c, nrows_c, ncols_c, error)
828 829
      use iso_c_binding
      import elpa_t
830
      implicit none
831 832
      class(elpa_t)                   :: self
      character*1                     :: uplo_a, uplo_c
833
      integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
834
#ifdef USE_ASSUMED_SIZE
835
      complex(kind=c_float_complex)   :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
836
#else
837
      complex(kind=c_float_complex)   :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
838
#endif
839
      integer, optional               :: error
840 841 842
    end subroutine
  end interface

843
  !> \brief abstract definition of interface to do a cholesky decomposition of a double real matrix
844 845 846 847
  !>
  !>  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"
848
  !>
849 850 851
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be decomposed
852
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
853
  abstract interface
854
    subroutine elpa_cholesky_d_i (self, a, error)
855 856
      use iso_c_binding
      import elpa_t
857
      implicit none
858 859 860 861 862 863
      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
864
      integer, optional               :: error
865 866 867
    end subroutine
  end interface

868
  !> \brief abstract definition of interface to do a cholesky decomposition of a single real matrix
869 870 871 872 873
  !>
  !>  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"
  !> 
874 875 876
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be decomposed
877
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
878
  abstract interface
879
    subroutine elpa_cholesky_f_i(self, a, error)
880 881
      use iso_c_binding
      import elpa_t
882
      implicit none
883 884 885 886 887 888
      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
889
      integer, optional               :: error
890 891 892
    end subroutine
  end interface

893
  !> \brief abstract definition of interface to do a cholesky decomposition of a double complex matrix
894 895 896 897 898
  !>
  !>  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"
  !> 
899 900 901
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double complex matrix: the matrix to be decomposed
902
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
903
  abstract interface
904
    subroutine elpa_cholesky_dc_i (self, a, error)
905 906
      use iso_c_binding
      import elpa_t
907
      implicit none
908 909 910 911 912 913
      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
914
      integer, optional               :: error
915 916 917
    end subroutine
  end interface

918
  !> \brief abstract definition of interface to do a cholesky decomposition of a single complex matrix
919 920 921 922 923
  !>
  !>  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"
  !> 
924 925 926
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single complex matrix: the matrix to be decomposed
927
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
928
  abstract interface
929
    subroutine elpa_cholesky_fc_i (self, a, error)
930 931
      use iso_c_binding
      import elpa_t
932
      implicit none
933 934 935 936 937 938
      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
939
      integer, optional               :: error
940 941 942
    end subroutine
  end interface

943
  !> \brief abstract definition of interface to invert a triangular double real matrix
944 945 946 947 948
  !>
  !>  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"
  !>
949 950 951
  !> Parameters
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           double real matrix: the matrix to be inverted
952
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
953
  abstract interface
954
    subroutine elpa_invert_trm_d_i (self, a, error)
955 956
      use iso_c_binding
      import elpa_t
957
      implicit none
958 959 960 961 962 963
      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
964
      integer, optional               :: error
965 966 967
    end subroutine
  end interface

968
  !> \brief abstract definition of interface to invert a triangular single real matrix
969
  !> Parameters
970 971 972 973 974
  !>
  !>  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"
  !>
975 976
  !> \param   self        class(elpa_t), the ELPA object
  !> \param   a           single real matrix: the matrix to be inverted
977
  !> \param   error       integer, optional : error code, which can be queried with elpa_strerr
978
  abstract interface
979
    subroutine elpa_invert_trm_f_i (self, a, error)
980 981
      use iso_c_binding
      import elpa_t