magma.F90 3.4 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7 8 9 10 11 12 13 14
MODULE m_magma
  use m_juDFT
  INTEGER,PARAMETER :: NGPU_CONST=1
  !**********************************************************
  !     Solve the generalized eigenvalue problem
  !     using the MAGMA library for multiple GPUs
  !**********************************************************
CONTAINS
15
  SUBROUTINE magma_diag(nsize,eig,ne,a_r,b_r,z_r,a_c,b_c,z_c)
16
#ifdef CPP_MAGMA
17
    use magma
18
#endif    
19 20 21 22 23 24 25 26 27
#include"cpp_double.h"
    IMPLICIT NONE

    ! ... Arguments ...

    INTEGER, INTENT (IN) :: nsize
  
    REAL,    INTENT(OUT) :: eig(:)
    INTEGER, INTENT(INOUT) :: ne
28 29 30 31 32

    REAL, OPTIONAL,ALLOCATABLE, INTENT (INOUT) :: a_r(:),b_r(:)
    REAL, OPTIONAL,ALLOCATABLE, INTENT (INOUT) :: z_r(:,:)
    COMPLEX, OPTIONAL,ALLOCATABLE, INTENT (INOUT) :: a_c(:),b_c(:)
    COMPLEX, OPTIONAL,ALLOCATABLE, INTENT (INOUT) :: z_c(:,:)
33

34
#ifdef CPP_MAGMA
35 36 37 38 39 40 41 42

    ! ... Local Variables ..
    INTEGER iind,ind1,ind2,info,lwork,liwork,lrwork,err,i,mout(1)
    REAL eigTemp(nsize)
    LOGICAL:: initialized=.false.

    REAL,    ALLOCATABLE :: rwork(:)
    INTEGER, ALLOCATABLE :: iwork(:)
43 44 45

    REAL, ALLOCATABLE :: largea_r(:,:),largeb_r(:,:)
    COMPLEX, ALLOCATABLE :: largea_c(:,:),largeb_c(:,:)
46
    COMPLEX,ALLOCATABLE :: work(:)
47 48 49
    
    LOGICAL :: l_real
    l_real=present(a_r)
50 51 52 53 54 55 56 57 58 59 60 61

    print *,"MAGMA start"
    IF (.NOT.initialized) THEN
       initialized=.true.
       call magmaf_init()
       print *,"MAGMA init"
    ENDIF

    !**********************************
    !expand from packed to full storage
    !**********************************
    !hamiltonian
62 63 64
    if (l_real) THEN
       call packed_to_full(nsize,a_r,largea_r)
       call packed_to_full(nsize,b_r,largeb_r)
65
       !deallocate(a_r,b_r)
66 67 68
    ELSE
       call packed_to_full(nsize,a_c,largea_c)
       call packed_to_full(nsize,b_c,largeb_c)
69
       !deallocate(a_c,b_c)
70 71 72
    Endif

    if (l_real) call juDFT_error("REAL diagonalization not implemented in magma.F90")
73 74 75

    !Query the workspace size 
    allocate(work(1),rwork(1),iwork(1))
76 77
    print *,"Magma workspace query"
    call flush()
78
    call magmaf_zhegvdx(1,'v','i','l',nsize,largea_c,nsize,largeb_c,nsize,&
79 80 81 82 83 84 85 86 87
         0.0,0.0,1,ne,mout,eigTemp,work,-1,rwork,-1,iwork,-1,err)
    lwork=work(1)
    lrwork=rwork(1)
    liwork=iwork(1)
    print*,"MAGMA:",lwork,lrwork,liwork
    deallocate(work,rwork,iwork)
    allocate(work(lwork),rwork(lrwork),iwork(liwork))
    if (err/=0) call juDFT_error("Failed to allocate workspaces",calledby="magma.F90")
    !Now the diagonalization
88 89
    print *,"Magma diagonalization"
    print *,nsize,shape(largea_c),shape(eigTemp),ne
90
    call magmaf_zhegvdx(1,'v','i','l',nsize,largea_c,nsize,largeb_c,nsize,&
91 92 93 94 95 96 97
         0.0,0.0,1,ne,mout,eigTemp,work,lwork,rwork,lrwork,iwork,liwork,err)
    print*,"MAGMA info:",err
    if (err/=0) call juDFT_error("Magma failed to diagonalize Hamiltonian")
    print *,"MAGMA mout:",mout

    DO i = 1, ne
       eig(i) = eigTemp(i)
98
       z_c(:nsize,i)=largea_c(:nsize,i)
99
    END DO
100
    !call judft_error("Eigenvectors are not calculated in MAGMA")
101
#endif
102 103 104
  END SUBROUTINE magma_diag
END MODULE m_magma