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-2022 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

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

  [gpu]= -- optional -- Flag (0 or 1): use GPU version

SIDE EFFECTS

  elpa_hdl(type<elpa_hdl_t>)= ELPA handle

SOURCE

202 subroutine elpa_func_allocate(elpa_hdl,gpu)
203 
204 !Arguments ------------------------------------
205  integer,intent(in),optional :: gpu
206  type(elpa_hdl_t),intent(inout) :: elpa_hdl
207 
208 !Local variables-------------------------------
209  integer :: err
210 
211 ! *********************************************************************
212 
213  err=0
214 
215 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
216  elpa_hdl%elpa => elpa_allocate()
217  if (err==ELPA_OK.and.present(gpu)) call elpa_hdl%elpa%set("gpu",gpu,err)
218 #else
219  if (err==0.and.present(gpu)) elpa_hdl%gpu=gpu
220 #endif
221 
222  call elpa_func_error_handler(err_code=err,err_varname="gpu")
223 
224  elpa_hdl%is_allocated=.true.
225 
226 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.
  elpa_hdl(type<elpa_hdl_t>)=handler for ELPA object

SOURCE

751 subroutine elpa_func_cholesky_complex(elpa_hdl,aa)
752 
753 !Arguments ------------------------------------
754 !scalars
755  type(elpa_hdl_t),intent(inout) :: elpa_hdl
756 !arrays
757  complex(dpc),intent(inout) :: aa(:,:)
758 
759 !Local variables-------------------------------
760  integer :: err
761  logical :: success
762 
763 ! *********************************************************************
764 
765  success=.true. ; err=0
766 
767  if (.not.elpa_hdl%is_allocated) then
768    ABI_BUG('ELPA handle not allocated!')
769  end if
770  if (.not.elpa_hdl%matrix_is_set) then
771    ABI_BUG('Matrix not set in ELPA handle!')
772  end if
773 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
774  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
775 #else
776  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
777 #endif
778 
779 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
780  call elpa_hdl%elpa%cholesky(aa,err)
781  success=(err==ELPA_OK)
782 #elif (defined HAVE_LINALG_ELPA_2016)
783  success = elpa_cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
784 &                                elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
785 #elif (defined HAVE_LINALG_ELPA_2015_11)
786  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
787 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
788 #elif (defined HAVE_LINALG_ELPA_2015_02)
789  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
790 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
791 #elif (defined HAVE_LINALG_ELPA_2014)
792  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
793 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
794 #elif (defined HAVE_LINALG_ELPA_2013)
795  call cholesky_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
796 &                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
797 #else
798 !ELPA-LEGACY-2017
799  success = elpa_cholesky_complex_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
800 &                                       elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
801 #endif
802 
803  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_cholesky_complex!')
804 
805 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

SOURCE

674 subroutine elpa_func_cholesky_real(elpa_hdl,aa)
675 
676 !Arguments ------------------------------------
677 !scalars
678  type(elpa_hdl_t),intent(inout) :: elpa_hdl
679 !arrays
680  real(dp),intent(inout) :: aa(:,:)
681 
682 !Local variables-------------------------------
683  integer :: err
684  logical :: success
685 
686 ! *********************************************************************
687 
688  success=.true. ; err=0
689 
690  if (.not.elpa_hdl%is_allocated) then
691    ABI_BUG('ELPA handle not allocated!')
692  end if
693  if (.not.elpa_hdl%matrix_is_set) then
694    ABI_BUG('Matrix not set in ELPA handle!')
695  end if
696 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
697  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
698 #else
699  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
700 #endif
701 
702 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
703  call elpa_hdl%elpa%cholesky(aa,err)
704  success=(err==ELPA_OK)
705 #elif (defined HAVE_LINALG_ELPA_2016)
706  success = elpa_cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
707 &                             elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
708 #elif (defined HAVE_LINALG_ELPA_2015_11)
709  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
710                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
711 #elif (defined HAVE_LINALG_ELPA_2015_02)
712  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
713                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
714 #elif (defined HAVE_LINALG_ELPA_2014)
715  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
716 &                   elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
717 #elif (defined HAVE_LINALG_ELPA_2013)
718  call cholesky_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
719 &                   elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
720 #else
721 !ELPA-LEGACY-2017
722  success = elpa_cholesky_real_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
723 &                                    elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
724 #endif
725 
726  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_cholesky_real!')
727 
728 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

SOURCE

245 subroutine elpa_func_deallocate(elpa_hdl)
246 
247 !Arguments ------------------------------------
248  type(elpa_hdl_t),intent(inout) :: elpa_hdl
249 
250 !Local variables-------------------------------
251  integer :: err
252 
253 ! *********************************************************************
254 
255  err=0
256 
257 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
258  call elpa_deallocate(elpa_hdl%elpa)
259 #endif
260 
261  elpa_hdl%matrix_is_set=.false.
262  elpa_hdl%is_allocated=.false.
263 
264 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

SOURCE

286 subroutine elpa_func_error_handler(err_code,err_msg,err_varname)
287 
288 !Arguments ------------------------------------
289  integer,optional :: err_code
290  character(len=*),optional :: err_msg,err_varname
291 
292 !Local variables-------------------------------
293  integer :: err_code_
294  character(len=500) :: msg
295  character(len=100) :: err_strg
296 
297 ! *********************************************************************
298 
299  err_code_=-100;if (present(err_code)) err_code_=err_code
300  if (err_code_==0) return
301 
302  err_strg=''
303 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
304  if (err_code_==ELPA_ERROR) err_strg='ELPA_ERROR'
305  if (err_code_==ELPA_ERROR_ENTRY_READONLY) err_strg='ELPA_ERROR_ENTRY_READONLY'
306  if (err_code_==ELPA_ERROR_ENTRY_NOT_FOUND) err_strg='ELPA_ERROR_ENTRY_NOT_FOUND'
307  if (err_code_==ELPA_ERROR_ENTRY_ALREADY_SET) err_strg='ELPA_ERROR_ENTRY_ALREADY_SET'
308  if (err_code_==ELPA_ERROR_ENTRY_INVALID_VALUE) err_strg='ELPA_ERROR_ENTRY_INVALID_VALUE'
309  if (err_code_==ELPA_ERROR_ENTRY_NO_STRING_REPRESENTATION) err_strg='ELPA_ERROR_NO_STRING_REPRESENTATION'
310 #endif
311 
312  write(msg,'(a)') 'ELPA library error!'
313  if (present(err_msg)) then
314    if (trim(err_msg)/="") write(msg,'(3a)') trim(msg),ch10,trim(err_msg)
315  end if
316  if (present(err_varname)) then
317    if (trim(err_varname)/="") write(msg,'(4a)') trim(msg),ch10,'Variable: ',trim(err_varname)
318  end if
319  if (trim(err_strg)/="") write(msg,'(4a)') trim(msg),ch10,'Error code: ',trim(err_strg)
320    ABI_ERROR(msg)
321 
322 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

SOURCE

344 subroutine elpa_func_get_communicators(elpa_hdl,mpi_comm_parent,process_row,process_col)
345 
346 !Arguments ------------------------------------
347  integer,intent(in)  :: mpi_comm_parent,process_row,process_col
348  type(elpa_hdl_t),intent(inout) :: elpa_hdl
349 
350 !Local variables-------------------------------
351  integer  :: err
352  character(len=20) :: varname
353 
354 ! *********************************************************************
355 
356  err=0 ; varname=''
357 
358  if (.not.elpa_hdl%is_allocated) then
359    ABI_BUG('ELPA handle not allocated!')
360  end if
361 
362 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
363  if (err==ELPA_OK) then
364    varname='mpi_comm_parent'
365    call elpa_hdl%elpa%set(trim(varname),mpi_comm_parent,err)
366  end if
367  if (err==ELPA_OK) then
368    varname='process_row'
369    call elpa_hdl%elpa%set(trim(varname),process_row,err)
370  end if
371  if (err==ELPA_OK) then
372    varname='process_col'
373    call elpa_hdl%elpa%set(trim(varname),process_col,err)
374  end if
375  if (err==ELPA_OK) then
376    varname=''
377    if (elpa_hdl%elpa%setup()/=ELPA_OK) err=ELPA_ERROR
378  endif
379 #else
380  elpa_hdl%mpi_comm_parent=mpi_comm_parent
381  elpa_hdl%process_row=process_row
382  elpa_hdl%process_col=process_col
383 #if (defined HAVE_LINALG_ELPA_2016)
384  err=elpa_get_communicators(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
385 #elif (defined HAVE_LINALG_ELPA_2015_11) || (defined HAVE_LINALG_ELPA_2015_02)
386  err=get_elpa_row_col_comms(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
387 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
388  call get_elpa_row_col_comms(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
389 #else
390 !ELPA-LEGACY-2017
391  err=elpa_get_communicators(mpi_comm_parent,process_row,process_col,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
392 #endif
393 #endif
394 
395  call elpa_func_error_handler(err_code=err,err_msg='Error in elpa_get_communicators',err_varname=varname)
396 
397 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

SOURCE

1102 subroutine elpa_func_hermitian_multiply_complex(elpa_hdl,uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1103 &                                               cc,local_nrows_c,local_ncols_c)
1104 
1105 !Arguments ------------------------------------
1106 !scalars
1107  integer,intent(in)  :: ncb,local_nrows_b,local_nrows_c,local_ncols_b,local_ncols_c
1108  character*1 :: uplo_a,uplo_c
1109  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1110 !arrays
1111  complex(dpc),intent(in) :: aa(:,:),bb(:,:)
1112  complex(dpc),intent(out) :: cc(:,:)
1113 
1114 !Local variables-------------------------------
1115  integer :: err
1116  logical :: success
1117 
1118 ! *********************************************************************
1119 
1120  success=.true. ; err=0
1121 
1122  if (.not.elpa_hdl%is_allocated) then
1123    ABI_BUG('ELPA handle not allocated!')
1124  end if
1125  if (.not.elpa_hdl%matrix_is_set) then
1126    ABI_BUG('Matrix not set in ELPA handle!')
1127  end if
1128 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
1129  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1130 #else
1131  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1132 #endif
1133  ABI_CHECK(size(bb)==local_nrows_b*local_ncols_b,'BUG: matrix B has wrong sizes!')
1134  ABI_CHECK(size(cc)==local_nrows_c*local_ncols_c,'BUG: matrix C has wrong sizes!')
1135 
1136 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
1137  call elpa_hdl%elpa%hermitian_multiply(uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1138 &                                      cc,local_nrows_c,local_ncols_c,err)
1139  success=(err==ELPA_OK)
1140 #elif (defined HAVE_LINALG_ELPA_2016)
1141  success = elpa_mult_ah_b_complex(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,elpa_hdl%local_ncols,&
1142 &          bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,&
1143 &          cc,local_nrows_c,local_ncols_c)
1144 #elif (defined HAVE_LINALG_ELPA_2015_11) || (defined HAVE_LINALG_ELPA_2015_02)
1145  call mult_ah_b_complex(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1146 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1147 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
1148  call mult_ah_b_complex(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1149 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1150 #else
1151 !ELPA-LEGACY-2017
1152  success =  elpa_mult_ah_b_complex_double(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,&
1153 &           elpa_hdl%local_ncols,bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,&
1154 &           elpa_hdl%elpa_comm_cols,cc,local_nrows_c,local_ncols_c)
1155 #endif
1156 
1157  if (.not.success) call elpa_func_error_handler(err_msg='Error in hermitian_multiply_complex!')
1158 
1159 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

SOURCE

1003 subroutine elpa_func_hermitian_multiply_real(elpa_hdl,uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1004 &                                            cc,local_nrows_c,local_ncols_c)
1005 
1006 !Arguments ------------------------------------
1007 !scalars
1008  integer,intent(in)  :: ncb,local_nrows_b,local_nrows_c,local_ncols_b,local_ncols_c
1009  character*1 :: uplo_a,uplo_c
1010  type(elpa_hdl_t),intent(inout) :: elpa_hdl
1011 !arrays
1012  real(dp),intent(in) :: aa(:,:),bb(:,:)
1013  real(dp),intent(out) :: cc(:,:)
1014 
1015 !Local variables-------------------------------
1016  integer :: err
1017  logical :: success
1018 
1019 ! *********************************************************************
1020 
1021  success=.true. ; err=0
1022 
1023  if (.not.elpa_hdl%is_allocated) then
1024    ABI_BUG('ELPA handle not allocated!')
1025  end if
1026  if (.not.elpa_hdl%matrix_is_set) then
1027    ABI_BUG('Matrix not set in ELPA handle!')
1028  end if
1029 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
1030  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
1031 #else
1032  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
1033 #endif
1034  ABI_CHECK(size(bb)==local_nrows_b*local_ncols_b,'BUG: matrix B has wrong sizes!')
1035  ABI_CHECK(size(cc)==local_nrows_c*local_ncols_c,'BUG: matrix C has wrong sizes!')
1036 
1037 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
1038  call elpa_hdl%elpa%hermitian_multiply(uplo_a,uplo_c,ncb,aa,bb,local_nrows_b,local_ncols_b,&
1039 &                                      cc,local_nrows_c,local_ncols_c,err)
1040  success=(err==ELPA_OK)
1041 #elif (defined HAVE_LINALG_ELPA_2016)
1042  success = elpa_mult_at_b_real(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,elpa_hdl%local_ncols,&
1043 &          bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,&
1044 &          cc,local_nrows_c,local_ncols_c)
1045 #elif (defined HAVE_LINALG_ELPA_2015_11) || (defined HAVE_LINALG_ELPA_2015_02)
1046  call mult_at_b_real(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1047 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1048 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
1049  call mult_at_b_real(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,bb,local_nrows_b,&
1050 &               elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,cc,local_nrows_c)
1051 #else
1052 !ELPA-LEGACY-2017
1053  success =  elpa_mult_at_b_real_double(uplo_a,uplo_c,elpa_hdl%na,ncb,aa,elpa_hdl%local_nrows,&
1054 &           elpa_hdl%local_ncols,bb,local_nrows_b,local_ncols_b,elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,&
1055 &           elpa_hdl%elpa_comm_cols,cc,local_nrows_c,local_ncols_c)
1056 #endif
1057 
1058  if (.not.success) call elpa_func_error_handler(err_msg='Error in hermitian_multiply_real!')
1059 
1060 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

SOURCE

120 subroutine elpa_func_init()
121 
122 !Local variables-------------------------------
123 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
124  integer,parameter :: elpa_min_version=20170403
125 #endif
126  logical :: success
127 !
128 ! *********************************************************************
129 
130  success=.true.
131 
132 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
133 
134 #if (defined HAVE_LINALG_ELPA_2016) || (defined HAVE_LINALG_ELPA_2015_11) || (defined HAVE_LINALG_ELPA_2015_02)
135 !This case is not supposed to happen
136  success=.false.
137  ABI_BUG('Wrong ELPA cpp directives!')
138 #elif (defined HAVE_LINALG_ELPA_2014) || (defined HAVE_LINALG_ELPA_2013)
139 !This case is not supposed to happen
140  success=.false.
141  ABI_BUG('Wrong ELPA cpp directives!')
142 #else
143  success=(elpa_init(elpa_min_version)==ELPA_OK)
144 #endif
145 
146 #else
147 !No init function
148 #endif
149 
150  if (.not.success) then
151    ABI_ERROR('Problem with ELPA (elpa_init)!')
152  end if
153 
154 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

SOURCE

906 subroutine elpa_func_invert_triangular_complex(elpa_hdl,aa)
907 
908 !Arguments ------------------------------------
909 !scalars
910  type(elpa_hdl_t),intent(inout) :: elpa_hdl
911 !arrays
912  complex(dpc),intent(inout) :: aa(:,:)
913 
914 !Local variables-------------------------------
915  integer :: err
916  logical :: success
917 
918 ! *********************************************************************
919 
920  success=.true. ; err=0
921 
922  if (.not.elpa_hdl%is_allocated) then
923    ABI_BUG('ELPA handle not allocated!')
924  end if
925  if (.not.elpa_hdl%matrix_is_set) then
926    ABI_BUG('Matrix not set in ELPA handle!')
927  end if
928 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
929  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
930 #else
931  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
932 #endif
933 
934 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
935  call elpa_hdl%elpa%invert_triangular(aa,err)
936  success=(err==ELPA_OK)
937 #elif (defined HAVE_LINALG_ELPA_2016)
938  success = elpa_invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
939 &                                  elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
940 #elif (defined HAVE_LINALG_ELPA_2015_11)
941  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
942                          elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
943 #elif (defined HAVE_LINALG_ELPA_2015_02)
944  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
945                          elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.,success)
946 #elif (defined HAVE_LINALG_ELPA_2014)
947  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
948 &                        elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,success)
949  success=.true. ! Sometimes get unexpected success=false
950 #elif (defined HAVE_LINALG_ELPA_2013)
951  call invert_trm_complex(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
952 &                        elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
953 #else
954 !ELPA-LEGACY-2017
955  success = elpa_invert_trm_complex_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
956 &                                         elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
957 #endif
958 
959  if (.not.success) call elpa_func_error_handler(err_msg='Error in invert_trianguler_complex!')
960 
961 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

SOURCE

828 subroutine elpa_func_invert_triangular_real(elpa_hdl,aa)
829 
830 !Arguments ------------------------------------
831 !scalars
832  type(elpa_hdl_t),intent(inout) :: elpa_hdl
833 !arrays
834  real(dp),intent(inout) :: aa(:,:)
835 
836 !Local variables-------------------------------
837  integer :: err
838  logical :: success
839 
840 ! *********************************************************************
841 
842  success=.true. ; err=0
843 
844  if (.not.elpa_hdl%is_allocated) then
845    ABI_BUG('ELPA handle not allocated!')
846  end if
847  if (.not.elpa_hdl%matrix_is_set) then
848    ABI_BUG('Matrix not set in ELPA handle!')
849  end if
850 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
851  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
852 #else
853  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
854 #endif
855 
856 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
857  call elpa_hdl%elpa%invert_triangular(aa,err)
858  success=(err==ELPA_OK)
859 #elif (defined HAVE_LINALG_ELPA_2016)
860  success = elpa_invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
861 &                               elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
862 #elif (defined HAVE_LINALG_ELPA_2015_11)
863  call invert_trm_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.,success)
865 #elif (defined HAVE_LINALG_ELPA_2015_02)
866  call invert_trm_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 invert_trm_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  success=.true. ! Sometimes get unexpected success=false
872 #elif (defined HAVE_LINALG_ELPA_2013)
873  call invert_trm_real(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,&
874 &                     elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
875 #else
876 !ELPA-LEGACY-2017
877  success = elpa_invert_trm_real_double(elpa_hdl%na,aa,elpa_hdl%local_nrows,elpa_hdl%nblk,elpa_hdl%local_ncols,&
878 &                                      elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,.false.)
879 #endif
880 
881  if (.not.success) call elpa_func_error_handler(err_msg='Error in invert_trianguler_real!')
882 
883 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

SOURCE

421 subroutine elpa_func_set_matrix(elpa_hdl,na,nblk,local_nrows,local_ncols)
422 
423 !Arguments ------------------------------------
424 !scalars
425  integer,intent(in) :: na,nblk,local_nrows,local_ncols
426  type(elpa_hdl_t),intent(inout) :: elpa_hdl
427 !arrays
428 
429 !Local variables-------------------------------
430  integer :: err
431  character(len=15) :: varname
432 
433 ! *********************************************************************
434 
435  err=0 ; varname=''
436 
437  if (.not.elpa_hdl%is_allocated) then
438    ABI_BUG('ELPA handle not allocated!')
439  end if
440 
441 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
442  if (err==ELPA_OK) then
443    varname="na"
444    call elpa_hdl%elpa%set(trim(varname),na,err)
445  end if
446  if (err==ELPA_OK) then
447    varname="nblk"
448    call elpa_hdl%elpa%set(trim(varname),nblk,err)
449  end if
450  if (err==ELPA_OK) then
451    varname="local_nrows"
452    call elpa_hdl%elpa%set(trim(varname),local_nrows,err)
453  end if
454  if (err==ELPA_OK) then
455    varname="local_ncols"
456    call elpa_hdl%elpa%set(trim(varname),local_ncols,err)
457  end if
458 #else
459  elpa_hdl%na=na
460  elpa_hdl%nblk=nblk
461  elpa_hdl%local_nrows=local_nrows
462  elpa_hdl%local_ncols=local_ncols
463 #endif
464 
465  call elpa_func_error_handler(err_code=err,err_msg='Error during matrix initialization',err_varname=varname)
466 
467  elpa_hdl%matrix_is_set=.true.
468 
469 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

SOURCE

590 subroutine elpa_func_solve_evp_1stage_complex(elpa_hdl,aa,qq,ev,nev)
591 
592 !Arguments ------------------------------------
593 !scalars
594  integer,intent(in)  :: nev
595  type(elpa_hdl_t),intent(inout) :: elpa_hdl
596 !arrays
597  complex(dpc),intent(inout) :: aa(:,:)
598  real(dp),intent(out) :: ev(:)
599  complex(dpc),intent(out) :: qq(:,:)
600 
601 !Local variables-------------------------------
602  integer :: err
603  logical  :: success
604 
605 ! *********************************************************************
606 
607  success=.true. ; err=0
608 
609  if (.not.elpa_hdl%is_allocated) then
610    ABI_BUG('ELPA handle not allocated!')
611  end if
612  if (.not.elpa_hdl%matrix_is_set) then
613    ABI_BUG('Matrix not set in ELPA handle!')
614  end if
615 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
616  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
617  ABI_CHECK(size(qq)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix Q has wrong sizes!')
618  ABI_CHECK(size(ev)==elpa_hdl%elpa%na,'BUG: matrix EV has wrong sizes!')
619 #else
620  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
621  ABI_CHECK(size(qq)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix Q has wrong sizes!')
622  ABI_CHECK(size(ev)==elpa_hdl%na,'BUG: matrix EV has wrong sizes!')
623 #endif
624 
625 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
626  if (err==ELPA_OK) call elpa_hdl%elpa%set('nev',nev,err)
627  if (err==ELPA_OK) call elpa_hdl%elpa%set("solver",ELPA_SOLVER_1STAGE,err)
628  if (err==ELPA_OK) call elpa_hdl%elpa%eigenvectors(aa,ev,qq,err)
629  success=(err==ELPA_OK)
630 #elif  (defined HAVE_LINALG_ELPA_2016)
631  success=elpa_solve_evp_complex_1stage(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
632 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
633 #elif (defined HAVE_LINALG_ELPA_2015_11)
634  success=solve_evp_complex(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
635 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
636 #elif (defined HAVE_LINALG_ELPA_2015_02) || (defined HAVE_LINALG_ELPA_2014)
637  success=solve_evp_complex(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
638 &                       elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
639 #elif (defined HAVE_LINALG_ELPA_2013)
640  call solve_evp_complex(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
641 &                    elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
642 #else
643 !ELPA-LEGACY-2017
644  success=elpa_solve_evp_complex_1stage_double(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
645 &        elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,elpa_hdl%mpi_comm_parent,.false.)
646 #endif
647 
648  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_evp_1stage_complex!')
649 
650 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

SOURCE

500 subroutine elpa_func_solve_evp_1stage_real(elpa_hdl,aa,qq,ev,nev)
501 
502 !Arguments ------------------------------------
503 !scalars
504  integer,intent(in)  :: nev
505  type(elpa_hdl_t),intent(inout) :: elpa_hdl
506 !arrays
507  real(dp),intent(inout) :: aa(:,:)
508  real(dp),intent(out) :: ev(:),qq(:,:)
509 
510 !Local variables-------------------------------
511  integer :: err
512  logical  :: success
513 
514 ! *********************************************************************
515 
516  success=.true. ; err=0
517 
518  if (.not.elpa_hdl%is_allocated) then
519    ABI_BUG('ELPA handle not allocated!')
520  end if
521  if (.not.elpa_hdl%matrix_is_set) then
522    ABI_BUG('Matrix not set in ELPA handle!')
523  end if
524 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
525  ABI_CHECK(size(aa)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix A has wrong sizes!')
526  ABI_CHECK(size(qq)==elpa_hdl%elpa%local_nrows*elpa_hdl%elpa%local_ncols,'BUG: matrix Q has wrong sizes!')
527  ABI_CHECK(size(ev)==elpa_hdl%elpa%na,'BUG: matrix EV has wrong sizes!')
528 #else
529  ABI_CHECK(size(aa)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix A has wrong sizes!')
530  ABI_CHECK(size(qq)==elpa_hdl%local_nrows*elpa_hdl%local_ncols,'BUG: matrix Q has wrong sizes!')
531  ABI_CHECK(size(ev)==elpa_hdl%na,'BUG: matrix EV has wrong sizes!')
532 #endif
533 
534 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
535  if (err==ELPA_OK) call elpa_hdl%elpa%set('nev',nev,err)
536  if (err==ELPA_OK) call elpa_hdl%elpa%set("solver",ELPA_SOLVER_1STAGE,err)
537  if (err==ELPA_OK) call elpa_hdl%elpa%eigenvectors(aa,ev,qq,err)
538  success=(err==ELPA_OK)
539 #elif  (defined HAVE_LINALG_ELPA_2016)
540  success=elpa_solve_evp_real_1stage(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
541 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
542 #elif (defined HAVE_LINALG_ELPA_2015_11)
543  success=solve_evp_real(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
544 &                       elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
545 #elif (defined HAVE_LINALG_ELPA_2015_02) || (defined HAVE_LINALG_ELPA_2014)
546  success=solve_evp_real(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
547 &                       elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
548 #elif (defined HAVE_LINALG_ELPA_2013)
549  call solve_evp_real(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
550 &                    elpa_hdl%nblk,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols)
551 #else
552 !ELPA-LEGACY-2017
553  success=elpa_solve_evp_real_1stage_double(elpa_hdl%na,nev,aa,elpa_hdl%local_nrows,ev,qq,elpa_hdl%local_nrows,&
554 &        elpa_hdl%nblk,elpa_hdl%local_ncols,elpa_hdl%elpa_comm_rows,elpa_hdl%elpa_comm_cols,elpa_hdl%mpi_comm_parent,.false.)
555 #endif
556 
557  if (.not.success) call elpa_func_error_handler(err_msg='Error in solve_evp_1stage_real!')
558 
559 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

SOURCE

172 subroutine elpa_func_uninit()
173 
174 ! *********************************************************************
175 
176 #ifdef HAVE_LINALG_ELPA_FORTRAN2008
177  call elpa_uninit()
178 #else
179 !No uninit function
180 #endif
181 
182 end subroutine elpa_func_uninit