elpa_20180525_onenode.F90 3.27 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11
!-------------------------------------------------------------------------------
! 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.
!--------------------------------------------------------------------------------

MODULE m_elpa_onenode
CONTAINS
  SUBROUTINE elpa_diag_onenode(hmat,smat,ne,eig,ev)
    !
    !----------------------------------------------------
12 13
    ! Sigensystem solver 
    !  Uses the ELPA (1 node version)
14 15 16 17 18 19 20 21 22 23
    !
    !
    ! hmat ..... Hamiltonian matrix
    ! smat ..... overlap matrix
    ! ne ....... number of ev's searched (and found) on this node
    !            On input, overall number of ev's searched,
    !            On output, local number of ev's found
    ! eig ...... eigenvalues, output
    ! ev ....... eigenvectors, output
    !
24
    ! U.Alekseeva     Nov. 2018
25 26 27
    !----------------------------------------------------
    USE m_juDFT
    USE m_types_mat
28
#ifdef CPP_ELPA_ONENODE
29
    USE elpa
Uliana Alekseeva's avatar
Uliana Alekseeva committed
30 31 32
#endif
#ifdef CPP_GPU
    USE nvtx
33
#endif
34 35 36 37
    IMPLICIT NONE

    CLASS(t_mat),INTENT(INOUT)    :: hmat,smat
    CLASS(t_mat),ALLOCATABLE,INTENT(OUT)::ev
38
    REAL,INTENT(OUT)              :: eig(:)
39 40
    INTEGER,INTENT(INOUT)         :: ne
    
41
#ifdef CPP_ELPA_ONENODE
42 43 44
    !...  Local variables
    !
    INTEGER           :: err
45 46 47
    REAL,ALLOCATABLE  :: eig2(:)
    TYPE(t_mat)       :: ev_dist
    INTEGER           :: kernel
48 49
    CLASS(elpa_t),pointer :: elpa_obj

Uliana Alekseeva's avatar
Uliana Alekseeva committed
50 51 52
#ifdef CPP_GPU
    call nvtxStartRange("ELPA",5)    
#endif
53 54
    err = elpa_init(20180525)
    elpa_obj => elpa_allocate()
55
       
56 57
    ALLOCATE ( eig2(hmat%matsize1), stat=err ) ! The eigenvalue array
    IF (err.NE.0) CALL juDFT_error('Failed to allocated "eig2"', calledby ='elpa')
58

59 60
    CALL ev_dist%init(hmat)! Eigenvectors
    IF (err.NE.0) CALL juDFT_error('Failed to allocated "ev_dist"',calledby ='elpa')
61
       
62 63 64 65 66 67
    CALL elpa_obj%set("na", hmat%matsize1, err)
    CALL elpa_obj%set("nev", ne, err)
    CALL elpa_obj%set("local_nrows", hmat%matsize1, err)
    CALL elpa_obj%set("local_ncols", hmat%matsize2, err)
    CALL elpa_obj%set("nblk",hmat%matsize1, err)
    CALL elpa_obj%set("blacs_context", -1, err)
68 69 70
#ifdef CPP_GPU
    CALL elpa_obj%set("gpu",1,err)
#endif
71
    err = elpa_obj%setup()
Uliana Alekseeva's avatar
Uliana Alekseeva committed
72 73
    call elpa_obj%get("solver", kernel)
    print *, "elpa uses " // elpa_int_value_to_string("solver", kernel) // " solver"
74

75 76
    CALL hmat%add_transpose(hmat)       
    CALL smat%add_transpose(smat)       
77

Uliana Alekseeva's avatar
Uliana Alekseeva committed
78 79 80
#ifdef CPP_GPU
    call nvtxStartRange("EigVec",7)    
#endif
81 82 83 84 85
    IF (hmat%l_real) THEN
        CALL elpa_obj%generalized_eigenvectors(hmat%data_r,smat%data_r,eig2, ev_dist%data_r, .FALSE.,err)
    ELSE
        CALL elpa_obj%generalized_eigenvectors(hmat%data_c,smat%data_c,eig2, ev_dist%data_c, .FALSE., err)
    ENDIF
Uliana Alekseeva's avatar
Uliana Alekseeva committed
86 87 88
#ifdef CPP_GPU
    call nvtxEndRange!("EigVec",8)    
#endif
89
       
90 91 92
    CALL elpa_deallocate(elpa_obj)
    CALL elpa_uninit()
    ! END of ELPA stuff
93

94 95 96
    eig(1:ne) = eig2(1:ne)
    DEALLOCATE(eig2)
       
97 98
    ALLOCATE(t_mat::ev)
    CALL ev%alloc(hmat%l_real,hmat%matsize1,ne)
99
    CALL ev%copy(ev_dist,1,1)
100

Uliana Alekseeva's avatar
Uliana Alekseeva committed
101 102 103
#ifdef CPP_GPU
    call nvtxEndRange!("ELPA",7)    
#endif
104
#endif
105 106 107
 
END SUBROUTINE elpa_diag_onenode
END MODULE m_elpa_onenode