TABLE OF CONTENTS


ABINIT/m_elpa [ Modules ]

[ Top ] [ Modules ]

NAME

 m_elpa

FUNCTION

 This module contains interfaces to ELPA library methods
 See http://elpa.mpcdf.mpg.de

COPYRIGHT

 Copyright (C) 2016-2018 ABINIT group (MT)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

 18 #if defined HAVE_CONFIG_H
 19 #include "config.h"
 20 #endif
 21 
 22 #include "abi_common.h"
 23 
 24 #if (defined HAVE_LINALG_ELPA) && (defined HAVE_LINALG_ELPA_2017)
 25 #define HAVE_ELPA_FORTRAN2008
 26 #else
 27 #undef HAVE_ELPA_FORTRAN2008
 28 #endif
 29 
 30 module m_elpa
 31 
 32  use defs_basis
 33  use m_errors
 34 
 35 #ifdef HAVE_LINALG_ELPA
 36 #ifdef HAVE_ELPA_FORTRAN2008
 37  use elpa
 38 #else
 39  use elpa1
 40 #endif
 41 #endif
 42 
 43  implicit none
 44 
 45  private
 46 
 47 #ifdef HAVE_LINALG_ELPA
 48 
 49 !Public procedures
 50 !Had to choose names different from those provided by elpa
 51  public :: elpa_func_init                ! Init ELPA
 52  public :: elpa_func_uninit              ! End ELPA
 53  public :: elpa_func_allocate            ! Allocate a ELPA handle and set up MPI information
 54  public :: elpa_func_deallocate          ! Deallocate a ELPA handle
 55  public :: elpa_func_error_handler       ! Manage errors (print a readable message)
 56  public :: elpa_func_get_communicators   ! Get rows and cols communicators (not supposed to be called directly)
 57  public :: elpa_func_set_matrix          ! Set matrix specifications in a ELPA handle
 58  public :: elpa_func_solve_evp_1stage    ! Solve the diagonalization problem (use a ELPA handle)
 59  public :: elpa_func_cholesky            ! Apply Cholesky transformation (use a ELPA handle)
 60  public :: elpa_func_invert_triangular   ! Invert triangular matrix (use a ELPA handle)
 61  public :: elpa_func_hermitian_multiply  ! Perform C := A**H * B (use a ELPA handle)
 62 
 63  interface elpa_func_solve_evp_1stage
 64    module procedure elpa_func_solve_evp_1stage_real
 65    module procedure elpa_func_solve_evp_1stage_complex
 66  end interface elpa_func_solve_evp_1stage
 67 
 68  interface elpa_func_cholesky
 69    module procedure elpa_func_cholesky_real
 70    module procedure elpa_func_cholesky_complex
 71  end interface elpa_func_cholesky
 72 
 73  interface elpa_func_invert_triangular
 74    module procedure elpa_func_invert_triangular_real
 75    module procedure elpa_func_invert_triangular_complex
 76  end interface elpa_func_invert_triangular
 77 
 78  interface elpa_func_hermitian_multiply
 79    module procedure elpa_func_hermitian_multiply_real
 80    module procedure elpa_func_hermitian_multiply_complex
 81  end interface elpa_func_hermitian_multiply
 82 
 83 !ELPA gneralized handle
 84  type,public :: elpa_hdl_t
 85    logical :: is_allocated=.false.
 86    logical :: matrix_is_set=.false.
 87 #ifdef HAVE_ELPA_FORTRAN2008
 88    class(elpa_t),pointer :: elpa
 89 #else
 90    integer :: mpi_comm_parent
 91    integer :: elpa_comm_rows,elpa_comm_cols
 92    integer :: process_row,process_col
 93    integer :: local_nrows,local_ncols
 94    integer :: na,nblk
 95    integer :: gpu=0
 96 #endif
 97  end type elpa_hdl_t
 98 
 99 #endif
100 
101 CONTAINS  !==============================================================================

m_elpa/elpa_func_allocate [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_allocate

FUNCTION

  Allocate a ELPA handle and set it up with communicators specification

INPUTS

  mpi_comm_parent=Global communicator for the calculations
  process_row=Row coordinate of the calling process in the process grid
  process_col=Column coordinate of the calling process in the process grid
  [gpu]= -- optional -- Flag (0 or 1): use GPU version

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)= ELPA handle

PARENTS

      m_slk

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

247 subroutine elpa_func_allocate(elpa_hdl,mpi_comm_parent,process_row,process_col,gpu)
248 
249 
250 !This section has been created automatically by the script Abilint (TD).
251 !Do not modify the following lines by hand.
252 #undef ABI_FUNC
253 #define ABI_FUNC 'elpa_func_allocate'
254 !End of the abilint section
255 
256  implicit none
257 
258 !Arguments ------------------------------------
259  integer,intent(in) :: mpi_comm_parent,process_row,process_col
260  integer,intent(in),optional :: gpu
261  type(elpa_hdl_t),intent(inout) :: elpa_hdl
262 
263 !Local variables-------------------------------
264  integer :: err
265 
266 ! *********************************************************************
267 
268  err=0
269 
270 #ifdef HAVE_ELPA_FORTRAN2008
271  elpa_hdl%elpa => elpa_allocate()
272  if (err==ELPA_OK.and.present(gpu)) call elpa_hdl%elpa%set("gpu",gpu,err)
273 #else
274  if (err==0.and.present(gpu)) elpa_hdl%gpu=gpu
275 #endif
276 
277  call elpa_func_error_handler(err_code=err,err_varname="gpu")
278 
279  elpa_hdl%is_allocated=.true.
280 
281  call elpa_func_get_communicators(elpa_hdl,mpi_comm_parent,process_row,process_col)
282 
283 end subroutine elpa_func_allocate

m_elpa/elpa_func_cholesky_complex [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_cholesky_complex

FUNCTION

  Wrapper to elpa_cholesky_complex ELPA function

INPUTS

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix which should be factorized.
                     Distribution is like in Scalapack.
                     Only upper triangle is needs to be set.
                     On return, the upper triangle contains the Cholesky factor
                     and the lower triangle is set to 0.#elif (defined HAVE_LINALG_ELPA_2016)
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

911 subroutine elpa_func_cholesky_complex(elpa_hdl,aa)
912 
913 
914 !This section has been created automatically by the script Abilint (TD).
915 !Do not modify the following lines by hand.
916 #undef ABI_FUNC
917 #define ABI_FUNC 'elpa_func_cholesky_complex'
918 !End of the abilint section
919 
920  implicit none
921 
922 !Arguments ------------------------------------
923 !scalars
924  type(elpa_hdl_t),intent(inout) :: elpa_hdl
925 !arrays
926  complex(dpc),intent(inout) :: aa(:,:)
927 
928 !Local variables-------------------------------
929  integer :: err
930  logical :: success
931 
932 ! *********************************************************************
933 
934  success=.true. ; err=0
935 
936  if (.not.elpa_hdl%is_allocated) then
937    MSG_BUG('ELPA handle not allocated!')
938  end if
939  if (.not.elpa_hdl%matrix_is_set) then
940    MSG_BUG('Matrix not set in ELPA handle!')
941  end if
942 #ifdef HAVE_ELPA_FORTRAN2008
943  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
944 #else
945  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
946 #endif
947 
948 #if  (defined HAVE_LINALG_ELPA_2017)
949  call elpa_hdl%elpa%cholesky(aa,err)
950  success=(err==ELPA_OK)
951 #elif (defined HAVE_LINALG_ELPA_2016)
952  success = elpa_cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
953 &                                elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
954 #elif (defined HAVE_LINALG_ELPA_2015)
955  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
956 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
957 #elif (defined HAVE_LINALG_ELPA_2014)
958  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
959 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
960 #elif (defined HAVE_LINALG_ELPA_2013)
961  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
962 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
963 #else
964 !ELPA-LEGACY-2017
965  success = elpa_cholesky_complex_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
966 &                                       elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
967 #endif
968 
969  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_cholesky_complex!')
970 
971 end subroutine elpa_func_cholesky_complex

m_elpa/elpa_func_cholesky_real [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_cholesky_real

FUNCTION

  Wrapper to elpa_cholesky_real ELPA function

INPUTS

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix which should be factorized.
                     Distribution is like in Scalapack.
                     Only upper triangle is needs to be set.
                     On return, the upper triangle contains the Cholesky factor
                     and the lower triangle is set to 0.
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

822 subroutine elpa_func_cholesky_real(elpa_hdl,aa)
823 
824 
825 !This section has been created automatically by the script Abilint (TD).
826 !Do not modify the following lines by hand.
827 #undef ABI_FUNC
828 #define ABI_FUNC 'elpa_func_cholesky_real'
829 !End of the abilint section
830 
831  implicit none
832 
833 !Arguments ------------------------------------
834 !scalars
835  type(elpa_hdl_t),intent(inout) :: elpa_hdl
836 !arrays
837  real(dp),intent(inout) :: aa(:,:)
838 
839 !Local variables-------------------------------
840  integer :: err
841  logical :: success
842 
843 ! *********************************************************************
844 
845  success=.true. ; err=0
846 
847  if (.not.elpa_hdl%is_allocated) then
848    MSG_BUG('ELPA handle not allocated!')
849  end if
850  if (.not.elpa_hdl%matrix_is_set) then
851    MSG_BUG('Matrix not set in ELPA handle!')
852  end if
853 #ifdef HAVE_ELPA_FORTRAN2008
854  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
855 #else
856  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
857 #endif
858 
859 #if  (defined HAVE_LINALG_ELPA_2017)
860  call elpa_hdl%elpa%cholesky(aa,err)
861  success=(err==ELPA_OK)
862 #elif (defined HAVE_LINALG_ELPA_2016)
863  success = elpa_cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
864 &                             elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
865 #elif (defined HAVE_LINALG_ELPA_2015)
866  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
867                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
868 #elif (defined HAVE_LINALG_ELPA_2014)
869  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
870 &                   elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
871 #elif (defined HAVE_LINALG_ELPA_2013)
872  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
873 &                   elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
874 #else
875 !ELPA-LEGACY-2017
876  success = elpa_cholesky_real_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
877 &                                    elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
878 #endif
879 
880  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_cholesky_real!')
881 
882 end subroutine elpa_func_cholesky_real

m_elpa/elpa_func_deallocate [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_deallocate_matrix

FUNCTION

  Deallocate a ELPA handle

INPUTS

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)= ELPA handle

PARENTS

CHILDREN

SOURCE

306 subroutine elpa_func_deallocate(elpa_hdl)
307 
308 
309 !This section has been created automatically by the script Abilint (TD).
310 !Do not modify the following lines by hand.
311 #undef ABI_FUNC
312 #define ABI_FUNC 'elpa_func_deallocate'
313 !End of the abilint section
314 
315  implicit none
316 
317 !Arguments ------------------------------------
318  type(elpa_hdl_t),intent(inout) :: elpa_hdl
319 
320 !Local variables-------------------------------
321  integer :: err
322 
323 ! *********************************************************************
324 
325  err=0
326 
327 #ifdef HAVE_ELPA_FORTRAN2008
328  call elpa_deallocate(elpa_hdl%elpa)
329 #endif
330 
331  elpa_hdl%matrix_is_set=.false.
332  elpa_hdl%is_allocated=.false.
333 
334 end subroutine elpa_func_deallocate

m_elpa/elpa_func_error_handler [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_error_handler

FUNCTION

  Handle ELPA errors

INPUTS

  [err_code]= --optional-- Error code
  [err_msg]= --optional-- Generic error message
  [err_varname]= -- optional-- Name of the ELPA variable related to the error

OUTPUT

  No output, only printing

PARENTS

      m_elpa

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

363 subroutine elpa_func_error_handler(err_code,err_msg,err_varname)
364 
365 
366 !This section has been created automatically by the script Abilint (TD).
367 !Do not modify the following lines by hand.
368 #undef ABI_FUNC
369 #define ABI_FUNC 'elpa_func_error_handler'
370 !End of the abilint section
371 
372  implicit none
373 
374 !Arguments ------------------------------------
375  integer,optional :: err_code
376  character(len=*),optional :: err_msg,err_varname
377 
378 !Local variables-------------------------------
379  integer :: err_code_
380  character(len=500) :: msg
381  character(len=100) :: err_strg
382 
383 ! *********************************************************************
384 
385  err_code_=-100;if (present(err_code)) err_code_=err_code
386  if (err_code_==0) return
387 
388  err_strg=''
389 #ifdef HAVE_ELPA_FORTRAN2008
390  if (err_code_==ELPA_ERROR) err_strg='ELPA_ERROR'
391  if (err_code_==ELPA_ERROR_ENTRY_READONLY) err_strg='ELPA_ERROR_ENTRY_READONLY'
392  if (err_code_==ELPA_ERROR_ENTRY_NOT_FOUND) err_strg='ELPA_ERROR_ENTRY_NOT_FOUND'
393  if (err_code_==ELPA_ERROR_ENTRY_ALREADY_SET) err_strg='ELPA_ERROR_ENTRY_ALREADY_SET'
394  if (err_code_==ELPA_ERROR_ENTRY_INVALID_VALUE) err_strg='ELPA_ERROR_ENTRY_INVALID_VALUE'
395  if (err_code_==ELPA_ERROR_ENTRY_NO_STRING_REPRESENTATION) err_strg='ELPA_ERROR_NO_STRING_REPRESENTATION'
396 #endif
397 
398  write(msg,'(a)') 'ELPA library error!'
399  if (present(err_msg)) then
400    if (trim(err_msg)/="") write(msg,'(3a)') trim(msg),ch10,trim(err_msg)
401  end if
402  if (present(err_varname)) then
403    if (trim(err_varname)/="") write(msg,'(4a)') trim(msg),ch10,'Variable: ',trim(err_varname)
404  end if
405  if (trim(err_strg)/="") write(msg,'(4a)') trim(msg),ch10,'Error code: ',trim(err_strg)
406    MSG_ERROR(msg)
407 
408 end subroutine elpa_func_error_handler

m_elpa/elpa_func_get_communicators [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_get_communicators

FUNCTION

  Wrapper to elpa_get_communicators ELPA function

INPUTS

  mpi_comm_parent=Global communicator for the calculations (in)
  process_row=Row coordinate of the calling process in the process grid (in)
  process_col=Column coordinate of the calling process in the process grid (in)

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)= ELPA handle

PARENTS

      m_elpa

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

437 subroutine elpa_func_get_communicators(elpa_hdl,mpi_comm_parent,process_row,process_col)
438 
439 
440 !This section has been created automatically by the script Abilint (TD).
441 !Do not modify the following lines by hand.
442 #undef ABI_FUNC
443 #define ABI_FUNC 'elpa_func_get_communicators'
444 !End of the abilint section
445 
446  implicit none
447 
448 !Arguments ------------------------------------
449  integer,intent(in)  :: mpi_comm_parent,process_row,process_col
450  type(elpa_hdl_t),intent(inout) :: elpa_hdl
451 
452 !Local variables-------------------------------
453  integer  :: err
454  character(len=20) :: varname
455 
456 ! *********************************************************************
457 
458  err=0 ; varname=''
459 
460  if (.not.elpa_hdl%is_allocated) then
461    MSG_BUG('ELPA handle not allocated!')
462  end if
463 
464 #if (defined HAVE_LINALG_ELPA_2017)
465  if (err==ELPA_OK) then
466    varname='mpi_comm_parent'
467    call elpa_hdl%elpa%set(trim(varname),mpi_comm_parent,err)
468  end if
469  if (err==ELPA_OK) then
470    varname='process_row'
471    call elpa_hdl%elpa%set(trim(varname),process_row,err)
472  end if
473  if (err==ELPA_OK) then
474    varname='process_col'
475    call elpa_hdl%elpa%set(trim(varname),process_col,err)
476  end if
477  if (err==ELPA_OK) then
478    varname=''
479    if (elpa_hdl%elpa%setup()/=ELPA_OK) err=ELPA_ERROR
480  endif
481 #else
482  elpa_hdl%mpi_comm_parent=mpi_comm_parent
483  elpa_hdl%process_row=process_row
484  elpa_hdl%process_col=process_col
485 #if (defined HAVE_LINALG_ELPA_2016)
486  err=elpa_get_communicators(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
487 #elif (defined HAVE_LINALG_ELPA_2015)
488  err=get_elpa_row_col_comms(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
489 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
490  call get_elpa_row_col_comms(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
491 #else
492 !ELPA-LEGACY-2017
493  err=get_elpa_communicators(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
494 #endif
495 #endif
496 
497  call elpa_func_error_handler(err_code=err,err_msg='Error in elpa_get_communicators',err_varname=varname)
498 
499 end subroutine elpa_func_get_communicators

m_elpa/elpa_func_hermitian_multiply_complex [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_hermitian_multiply_complex

FUNCTION

  Wrapper to elpa_hermitian_multiply_complex ELPA function
  Performs C := A**H * B

INPUTS

 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
  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!
  ncb=Number of columns of B and C
  aa(local_nrows,local_ncols)=Matrix A
  bb(ldb,local_ncols_c)=Matrix B
  local_nrows_b=Local rows of matrix B
  local_ncols_b=Local columns of matrix B
  local_nrows_c=Local rows of matrix C
  local_ncols_c=Local columns of matrix C

OUTPUT

  cc(local_nrows_c,local_ncols_c)=Matrix C

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

1310 subroutine elpa_func_hermitian_multiply_complex(elpa_hdl,uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1311 &                                               cc,local_nrows_c,local_ncols_c)
1312 
1313 
1314 !This section has been created automatically by the script Abilint (TD).
1315 !Do not modify the following lines by hand.
1316 #undef ABI_FUNC
1317 #define ABI_FUNC 'elpa_func_hermitian_multiply_complex'
1318 !End of the abilint section
1319 
1320  implicit none
1321 
1322 !Arguments ------------------------------------
1323 !scalars
1324  integer,intent(in)  :: ncb,local_nrows_b,local_nrows_c,local_ncols_b,local_ncols_c
1325  character*1 :: uplo_a,uplo_c
1326  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1327 !arrays
1328  complex(dpc),intent(in) :: aa(:,:),bb(:,:)
1329  complex(dpc),intent(out) :: cc(:,:)
1330 
1331 !Local variables-------------------------------
1332  integer :: err
1333  logical :: success
1334 
1335 ! *********************************************************************
1336 
1337  success=.true. ; err=0
1338 
1339  if (.not.elpa_hdl%is_allocated) then
1340    MSG_BUG('ELPA handle not allocated!')
1341  end if
1342  if (.not.elpa_hdl%matrix_is_set) then
1343    MSG_BUG('Matrix not set in ELPA handle!')
1344  end if
1345 #ifdef HAVE_ELPA_FORTRAN2008
1346  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1347 #else
1348  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1349 #endif
1350  ABI_CHECK(size(bb)==local_nrows_b*local_ncols_b,'BUG: matrix B has wrong sizes!')
1351  ABI_CHECK(size(cc)==local_nrows_c*local_ncols_c,'BUG: matrix C has wrong sizes!')
1352 
1353 #if  (defined HAVE_LINALG_ELPA_2017)
1354  call elpa_hdl%elpa%hermitian_multiply(uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1355 &                                      cc,local_nrows_c,local_ncols_c,err)
1356  success=(err==ELPA_OK)
1357 #elif (defined HAVE_LINALG_ELPA_2016)
1358  success = elpa_mult_ah_b_complex(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,elpa_hdl%local_ncols,&
1359 &          bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,&
1360 &          cc,local_nrows_c,local_ncols_c)
1361 #elif (defined HAVE_LINALG_ELPA_2015) || (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
1362  call mult_ah_b_complex(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1363 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1364 #else
1365 !ELPA-LEGACY-2017
1366  success =  elpa_mult_ah_b_complex_double(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,&
1367 &           elpa_hdl%local_ncols,bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,&
1368 &           elpa_hdl%elpa_comm_cols,cc,local_nrows_c,local_ncols_c)
1369 #endif
1370 
1371  if (.not.success) call elpa_func_error_handler(err_msg='Error in hermitian_multiply_complex!')
1372 
1373 end subroutine elpa_func_hermitian_multiply_complex

m_elpa/elpa_func_hermitian_multiply_real [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_hermitian_multiply_real

FUNCTION

  Wrapper to elpa_hermitian_multiply_real ELPA function
  Performs C := A**T * B

INPUTS

 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
  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!
  ncb=Number of columns of B and C
  aa(local_nrows,local_ncols)=Matrix A
  bb(ldb,local_ncols_c)=Matrix B
  local_nrows_b=Local rows of matrix B
  local_ncols_b=Local columns of matrix B
  local_nrows_c=Local rows of matrix C
  local_ncols_c=Local columns of matrix C

OUTPUT

  cc(local_nrows_c,local_ncols_c)=Matrix C

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

1199 subroutine elpa_func_hermitian_multiply_real(elpa_hdl,uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1200 &                                            cc,local_nrows_c,local_ncols_c)
1201 
1202 
1203 !This section has been created automatically by the script Abilint (TD).
1204 !Do not modify the following lines by hand.
1205 #undef ABI_FUNC
1206 #define ABI_FUNC 'elpa_func_hermitian_multiply_real'
1207 !End of the abilint section
1208 
1209  implicit none
1210 
1211 !Arguments ------------------------------------
1212 !scalars
1213  integer,intent(in)  :: ncb,local_nrows_b,local_nrows_c,local_ncols_b,local_ncols_c
1214  character*1 :: uplo_a,uplo_c
1215  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1216 !arrays
1217  real(dp),intent(in) :: aa(:,:),bb(:,:)
1218  real(dp),intent(out) :: cc(:,:)
1219 
1220 !Local variables-------------------------------
1221  integer :: err
1222  logical :: success
1223 
1224 ! *********************************************************************
1225 
1226  success=.true. ; err=0
1227 
1228  if (.not.elpa_hdl%is_allocated) then
1229    MSG_BUG('ELPA handle not allocated!')
1230  end if
1231  if (.not.elpa_hdl%matrix_is_set) then
1232    MSG_BUG('Matrix not set in ELPA handle!')
1233  end if
1234 #ifdef HAVE_ELPA_FORTRAN2008
1235  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1236 #else
1237  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1238 #endif
1239  ABI_CHECK(size(bb)==local_nrows_b*local_ncols_b,'BUG: matrix B has wrong sizes!')
1240  ABI_CHECK(size(cc)==local_nrows_c*local_ncols_c,'BUG: matrix C has wrong sizes!')
1241 
1242 #if  (defined HAVE_LINALG_ELPA_2017)
1243  call elpa_hdl%elpa%hermitian_multiply(uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1244 &                                      cc,local_nrows_c,local_ncols_c,err)
1245  success=(err==ELPA_OK)
1246 #elif (defined HAVE_LINALG_ELPA_2016)
1247  success = elpa_mult_at_b_real(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,elpa_hdl%local_ncols,&
1248 &          bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,&
1249 &          cc,local_nrows_c,local_ncols_c)
1250 #elif (defined HAVE_LINALG_ELPA_2015) || (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
1251  call mult_at_b_real(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1252 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1253 #else
1254 !ELPA-LEGACY-2017
1255  success =  elpa_mult_at_b_real_double(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,&
1256 &           elpa_hdl%local_ncols,bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,&
1257 &           elpa_hdl%elpa_comm_cols,cc,local_nrows_c,local_ncols_c)
1258 #endif
1259 
1260  if (.not.success) call elpa_func_error_handler(err_msg='Error in hermitian_multiply_real!')
1261 
1262 end subroutine elpa_func_hermitian_multiply_real

m_elpa/elpa_func_init [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_init

FUNCTION

  Wrapper to elpa_init ELPA function

INPUTS

OUTPUT

PARENTS

      m_abi_linalg

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

128 subroutine elpa_func_init()
129 
130 
131 !This section has been created automatically by the script Abilint (TD).
132 !Do not modify the following lines by hand.
133 #undef ABI_FUNC
134 #define ABI_FUNC 'elpa_func_init'
135 !End of the abilint section
136 
137  implicit none
138 
139 !Arguments ------------------------------------
140 
141 !Local variables-------------------------------
142 #if (defined HAVE_LINALG_ELPA_2017)
143  integer,parameter :: elpa_min_version=20170403
144 #endif
145  logical :: success
146 !
147 ! *********************************************************************
148 
149  success=.true.
150 
151 #if   (defined HAVE_LINALG_ELPA_2017)
152  success=(elpa_init(elpa_min_version)==ELPA_OK)
153 #elif (defined HAVE_LINALG_ELPA_2016) || (defined HAVE_LINALG_ELPA_2015)
154 !No init function
155 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
156 !No init function
157 #else
158  success=.true.
159 #endif
160 
161  if (.not.success) then
162    MSG_ERROR('Problem with ELPA (elpa_init)!')
163  end if
164 
165 end subroutine elpa_func_init

m_elpa/elpa_func_invert_triangular_complex [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_invert_triangular_complex

FUNCTION

  Wrapper to elpa_invert_triangular_complex ELPA function

INPUTS

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix which should be factorized.
                     Distribution is like in Scalapack.
                     Only upper triangle is needs to be set.
                     The lower triangle is not referenced.
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

1090 subroutine elpa_func_invert_triangular_complex(elpa_hdl,aa)
1091 
1092 
1093 !This section has been created automatically by the script Abilint (TD).
1094 !Do not modify the following lines by hand.
1095 #undef ABI_FUNC
1096 #define ABI_FUNC 'elpa_func_invert_triangular_complex'
1097 !End of the abilint section
1098 
1099  implicit none
1100 
1101 !Arguments ------------------------------------
1102 !scalars
1103  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1104 !arrays
1105  complex(dpc),intent(inout) :: aa(:,:)
1106 
1107 !Local variables-------------------------------
1108  integer :: err
1109  logical :: success
1110 
1111 ! *********************************************************************
1112 
1113  success=.true. ; err=0
1114 
1115  if (.not.elpa_hdl%is_allocated) then
1116    MSG_BUG('ELPA handle not allocated!')
1117  end if
1118  if (.not.elpa_hdl%matrix_is_set) then
1119    MSG_BUG('Matrix not set in ELPA handle!')
1120  end if
1121 #ifdef HAVE_ELPA_FORTRAN2008
1122  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1123 #else
1124  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1125 #endif
1126 
1127 #if  (defined HAVE_LINALG_ELPA_2017)
1128  call elpa_hdl%elpa%invert_triangular(aa,err)
1129  success=(err==ELPA_OK)
1130 #elif (defined HAVE_LINALG_ELPA_2016)
1131  success = elpa_invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
1132 &                                  elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
1133 #elif (defined HAVE_LINALG_ELPA_2015)
1134  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1135                          elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
1136 #elif (defined HAVE_LINALG_ELPA_2014)
1137  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1138 &                        elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
1139  success=.true. ! Sometimes get unexpected success=false
1140 #elif (defined HAVE_LINALG_ELPA_2013)
1141  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1142 &                        elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
1143 #else
1144 !ELPA-LEGACY-2017
1145  success = elpa_invert_trm_complex_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
1146 &                                         elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
1147 #endif
1148 
1149  if (.not.success) call elpa_func_error_handler(err_msg='Error in invert_trianguler_complex!')
1150 
1151 end subroutine elpa_func_invert_triangular_complex

m_elpa/elpa_func_invert_triangular_real [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_invert_triangular_real

FUNCTION

  Wrapper to elpa_invert_triangular_real ELPA function

INPUTS

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix which should be factorized.
                     Distribution is like in Scalapack.
                     Only upper triangle is needs to be set.
                     The lower triangle is not referenced.
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

1000 subroutine elpa_func_invert_triangular_real(elpa_hdl,aa)
1001 
1002 
1003 !This section has been created automatically by the script Abilint (TD).
1004 !Do not modify the following lines by hand.
1005 #undef ABI_FUNC
1006 #define ABI_FUNC 'elpa_func_invert_triangular_real'
1007 !End of the abilint section
1008 
1009  implicit none
1010 
1011 !Arguments ------------------------------------
1012 !scalars
1013  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1014 !arrays
1015  real(dp),intent(inout) :: aa(:,:)
1016 
1017 !Local variables-------------------------------
1018  integer :: err
1019  logical :: success
1020 
1021 ! *********************************************************************
1022 
1023  success=.true. ; err=0
1024 
1025  if (.not.elpa_hdl%is_allocated) then
1026    MSG_BUG('ELPA handle not allocated!')
1027  end if
1028  if (.not.elpa_hdl%matrix_is_set) then
1029    MSG_BUG('Matrix not set in ELPA handle!')
1030  end if
1031 #ifdef HAVE_ELPA_FORTRAN2008
1032  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1033 #else
1034  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1035 #endif
1036 
1037 #if  (defined HAVE_LINALG_ELPA_2017)
1038  call elpa_hdl%elpa%invert_triangular(aa,err)
1039  success=(err==ELPA_OK)
1040 #elif (defined HAVE_LINALG_ELPA_2016)
1041  success = elpa_invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
1042 &                               elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
1043 #elif (defined HAVE_LINALG_ELPA_2015)
1044  call invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1045                          elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
1046 #elif (defined HAVE_LINALG_ELPA_2014)
1047  call invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1048 &                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
1049  success=.true. ! Sometimes get unexpected success=false
1050 #elif (defined HAVE_LINALG_ELPA_2013)
1051  call invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
1052 &                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
1053 #else
1054 !ELPA-LEGACY-2017
1055  success = elpa_invert_trm_real_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
1056 &                                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
1057 #endif
1058 
1059  if (.not.success) call elpa_func_error_handler(err_msg='Error in invert_trianguler_real!')
1060 
1061 end subroutine elpa_func_invert_triangular_real

m_elpa/elpa_func_set_matrix [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_set_matrix

FUNCTION

  Set parameters decribing a matrix and it's MPI distribution
  in a ELPA handle

INPUTS

  na=Order of matrix A
  nblk=Blocksize of cyclic distribution, must be the same in both directions!
  local_nrows=Leading dimension of A
  local_ncols=Local columns of matrixes A and Q (eigenvectors)

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

      m_slk

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

530 subroutine elpa_func_set_matrix(elpa_hdl,na,nblk,local_nrows,local_ncols)
531 
532 
533 !This section has been created automatically by the script Abilint (TD).
534 !Do not modify the following lines by hand.
535 #undef ABI_FUNC
536 #define ABI_FUNC 'elpa_func_set_matrix'
537 !End of the abilint section
538 
539  implicit none
540 
541 !Arguments ------------------------------------
542 !scalars
543  integer,intent(in) :: na,nblk,local_nrows,local_ncols
544  type(elpa_hdl_t),intent(inout) :: elpa_hdl
545 !arrays
546 
547 !Local variables-------------------------------
548  integer :: err
549  character(len=15) :: varname
550 
551 ! *********************************************************************
552 
553  err=0 ; varname=''
554 
555  if (.not.elpa_hdl%is_allocated) then
556    MSG_BUG('ELPA handle not allocated!')
557  end if
558 
559 #ifdef HAVE_ELPA_FORTRAN2008
560  if (err==ELPA_OK) then
561    varname="na"
562    call elpa_hdl%elpa%set(trim(varname),na,err)
563  end if
564  if (err==ELPA_OK) then
565    varname="nblk"
566    call elpa_hdl%elpa%set(trim(varname),nblk,err)
567  end if
568  if (err==ELPA_OK) then
569    varname="local_nrows"
570    call elpa_hdl%elpa%set(trim(varname),local_nrows,err)
571  end if
572  if (err==ELPA_OK) then
573    varname="local_ncols"
574    call elpa_hdl%elpa%set(trim(varname),local_ncols,err)
575  end if
576 #else
577  elpa_hdl%na=na
578  elpa_hdl%nblk=nblk
579  elpa_hdl%local_nrows=local_nrows
580  elpa_hdl%local_ncols=local_ncols
581 #endif
582 
583  call elpa_func_error_handler(err_code=err,err_msg='Error during matrix initialization',err_varname=varname)
584 
585  elpa_hdl%matrix_is_set=.true.
586 
587 end subroutine elpa_func_set_matrix

m_elpa/elpa_func_solve_evp_1stage_complex [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_solve_evp_1stage_complex

FUNCTION

  Wrapper to elpa_solve_evp_complex_1stage ELPA function

INPUTS

  nev=Number of eigenvalues needed.

OUTPUT

  ev(na)=Eigenvalues of a, every processor gets the complete set
  qq(local_nrows,local_ncols)=Eigenvectors of aa
                     Distribution is like in Scalapack.
                     Must be always dimensioned to the full size (corresponding to (na,na))
                      even if only a part of the eigenvalues is needed.

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix for which eigenvalues are to be computed.
                    Distribution is like in Scalapack.
                    The full matrix must be set (not only one half like in scalapack).
                    Destroyed on exit (upper and lower half).
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

726 subroutine elpa_func_solve_evp_1stage_complex(elpa_hdl,aa,qq,ev,nev)
727 
728 
729 !This section has been created automatically by the script Abilint (TD).
730 !Do not modify the following lines by hand.
731 #undef ABI_FUNC
732 #define ABI_FUNC 'elpa_func_solve_evp_1stage_complex'
733 !End of the abilint section
734 
735  implicit none
736 
737 !Arguments ------------------------------------
738 !scalars
739  integer,intent(in)  :: nev
740  type(elpa_hdl_t),intent(inout) :: elpa_hdl
741 !arrays
742  complex(dpc),intent(inout) :: aa(:,:)
743  real(dp),intent(out) :: ev(:)
744  complex(dpc),intent(out) :: qq(:,:)
745 
746 !Local variables-------------------------------
747  integer :: err
748  logical  :: success
749 
750 ! *********************************************************************
751 
752  success=.true. ; err=0
753 
754  if (.not.elpa_hdl%is_allocated) then
755    MSG_BUG('ELPA handle not allocated!')
756  end if
757  if (.not.elpa_hdl%matrix_is_set) then
758    MSG_BUG('Matrix not set in ELPA handle!')
759  end if
760 #ifdef HAVE_ELPA_FORTRAN2008
761  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
762  ABI_CHECK(size(qq)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix Q has wrong sizes!')
763  ABI_CHECK(size(ev)==elpa_hdl%elpa%na,'BUG: matrix EV has wrong sizes!')
764 #else
765  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
766  ABI_CHECK(size(qq)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix Q has wrong sizes!')
767  ABI_CHECK(size(ev)==elpa_hdl%na,'BUG: matrix EV has wrong sizes!')
768 #endif
769 
770 #if  (defined HAVE_LINALG_ELPA_2017)
771  if (err==ELPA_OK) call elpa_hdl%elpa%set('nev',nev,err)
772  if (err==ELPA_OK) call elpa_hdl%elpa%set("solver",ELPA_SOLVER_1STAGE,err)
773  if (err==ELPA_OK) call elpa_hdl%elpa%eigenvectors(aa,ev,qq,err)
774  success=(err==ELPA_OK)
775 #elif  (defined HAVE_LINALG_ELPA_2016)
776  success=elpa_solve_evp_complex_1stage(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
777 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
778 #elif (defined HAVE_LINALG_ELPA_2015) || (defined HAVE_LINALG_ELPA_2014)
779  success=solve_evp_complex(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
780 &                       elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
781 #elif (defined HAVE_LINALG_ELPA_2013)
782  call solve_evp_complex(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
783 &                    elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
784 #else
785 !ELPA-LEGACY-2017
786  success=elpa_solve_evp_complex_1stage_double(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
787 &        elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,elpa_hdl%mpi_comm_parent,.false.)
788 #endif
789 
790  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_evp_1stage_complex!')
791 
792 end subroutine elpa_func_solve_evp_1stage_complex

m_elpa/elpa_func_solve_evp_1stage_real [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_solve_evp_1stage_real

FUNCTION

  Wrapper to elpa_solve_evp_real_1stage ELPA function

INPUTS

  nev=Number of eigenvalues needed.

OUTPUT

  ev(na)=Eigenvalues of a, every processor gets the complete set
  qq(local_nrows,local_ncols)=Eigenvectors of aa
                     Distribution is like in Scalapack.
                     Must be always dimensioned to the full size (corresponding to (na,na))
                      even if only a part of the eigenvalues is needed.

SIDE EFFECTS

  aa(local_nrows,local_ncols)=Distributed matrix for which eigenvalues are to be computed.
                    Distribution is like in Scalapack.
                    The full matrix must be set (not only one half like in scalapack).
                    Destroyed on exit (upper and lower half).
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

PARENTS

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

624 subroutine elpa_func_solve_evp_1stage_real(elpa_hdl,aa,qq,ev,nev)
625 
626 
627 !This section has been created automatically by the script Abilint (TD).
628 !Do not modify the following lines by hand.
629 #undef ABI_FUNC
630 #define ABI_FUNC 'elpa_func_solve_evp_1stage_real'
631 !End of the abilint section
632 
633  implicit none
634 
635 !Arguments ------------------------------------
636 !scalars
637  integer,intent(in)  :: nev
638  type(elpa_hdl_t),intent(inout) :: elpa_hdl
639 !arrays
640  real(dp),intent(inout) :: aa(:,:)
641  real(dp),intent(out) :: ev(:),qq(:,:)
642 
643 !Local variables-------------------------------
644  integer :: err
645  logical  :: success
646 
647 ! *********************************************************************
648 
649  success=.true. ; err=0
650 
651  if (.not.elpa_hdl%is_allocated) then
652    MSG_BUG('ELPA handle not allocated!')
653  end if
654  if (.not.elpa_hdl%matrix_is_set) then
655    MSG_BUG('Matrix not set in ELPA handle!')
656  end if
657 #ifdef HAVE_ELPA_FORTRAN2008
658  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
659  ABI_CHECK(size(qq)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix Q has wrong sizes!')
660  ABI_CHECK(size(ev)==elpa_hdl%elpa%na,'BUG: matrix EV has wrong sizes!')
661 #else
662  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
663  ABI_CHECK(size(qq)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix Q has wrong sizes!')
664  ABI_CHECK(size(ev)==elpa_hdl%na,'BUG: matrix EV has wrong sizes!')
665 #endif
666 
667 #if  (defined HAVE_LINALG_ELPA_2017)
668  if (err==ELPA_OK) call elpa_hdl%elpa%set('nev',nev,err)
669  if (err==ELPA_OK) call elpa_hdl%elpa%set("solver",ELPA_SOLVER_1STAGE,err)
670  if (err==ELPA_OK) call elpa_hdl%elpa%eigenvectors(aa,ev,qq,err)
671  success=(err==ELPA_OK)
672 #elif  (defined HAVE_LINALG_ELPA_2016)
673  success=elpa_solve_evp_real_1stage(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
674 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
675 #elif (defined HAVE_LINALG_ELPA_2015) || (defined HAVE_LINALG_ELPA_2014)
676  success=solve_evp_real(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
677 &                       elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
678 #elif (defined HAVE_LINALG_ELPA_2013)
679  call solve_evp_real(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
680 &                    elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
681 #else
682 !ELPA-LEGACY-2017
683  success=elpa_solve_evp_real_1stage_double(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
684 &        elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,elpa_hdl%mpi_comm_parent,.false.)
685 #endif
686 
687  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_evp_1stage_real!')
688 
689 end subroutine elpa_func_solve_evp_1stage_real

m_elpa/elpa_func_uninit [ Functions ]

[ Top ] [ m_elpa ] [ Functions ]

NAME

  elpa_func_uninit

FUNCTION

  Wrapper to elpa_uninit ELPA function

INPUTS

OUTPUT

PARENTS

      m_abi_linalg

CHILDREN

      elpa_func_error_handler,elpa_hdl%elpa%hermitian_multiply
      mult_ah_b_complex

SOURCE

190 subroutine elpa_func_uninit()
191 
192 
193 !This section has been created automatically by the script Abilint (TD).
194 !Do not modify the following lines by hand.
195 #undef ABI_FUNC
196 #define ABI_FUNC 'elpa_func_uninit'
197 !End of the abilint section
198 
199  implicit none
200 
201 !Arguments ------------------------------------
202 
203 !Local variables-------------------------------
204 
205 ! *********************************************************************
206 
207 #if   (defined HAVE_LINALG_ELPA_2017)
208  call elpa_uninit()
209 #elif (defined HAVE_LINALG_ELPA_2016) || (defined HAVE_LINALG_ELPA_2015)
210 !No uninit function
211 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
212 !No uninit function
213 #else
214 !No uninit function
215 #endif
216 
217 end subroutine elpa_func_uninit