hybrid.F90 5.13 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
MODULE m_calc_hybrid
8

9
  USE m_judft
10

11
CONTAINS
12

13 14
  SUBROUTINE calc_hybrid(eig_id,hybrid,kpts,atoms,input,DIMENSION,mpi,noco,cell,oneD,&
                         enpara,results,sym,xcpot,v,iter,iterHF)
15

16 17 18 19 20 21 22 23
    USE m_types
    USE m_mixedbasis
    USE m_coulombmatrix
    USE m_hf_init
    USE m_hf_setup
    USE m_hsfock
    USE m_eig66_io
    USE m_io_hybrid
24

25
    IMPLICIT NONE
26 27 28 29 30 31 32 33

    TYPE(t_xcpot_inbuild), INTENT(IN)    :: xcpot
    TYPE(t_mpi),           INTENT(IN)    :: mpi
    TYPE(t_dimension),     INTENT(IN)    :: DIMENSION
    TYPE(t_oneD),          INTENT(IN)    :: oneD
    TYPE(t_hybrid),        INTENT(INOUT) :: hybrid
    TYPE(t_input),         INTENT(IN)    :: input
    TYPE(t_noco),          INTENT(IN)    :: noco
34
    TYPE(t_enpara),        INTENT(IN)    :: enpara
35 36 37 38 39 40 41
    TYPE(t_results),       INTENT(INOUT) :: results
    TYPE(t_sym),           INTENT(IN)    :: sym  
    TYPE(t_cell),          INTENT(IN)    :: cell
    TYPE(t_kpts),          INTENT(IN)    :: kpts
    TYPE(t_atoms),         INTENT(IN)    :: atoms
    TYPE(t_potden),        INTENT(IN)    :: v

42 43
    INTEGER,               INTENT(IN)    :: iter
    INTEGER,               INTENT(INOUT) :: iterHF
44
    INTEGER,               INTENT(IN)    :: eig_id
45 46

    ! local variables
47
    INTEGER           :: jsp,nk,nred
48 49 50 51 52 53 54 55
    TYPE(t_hybdat)    :: hybdat
    type(t_lapw)      :: lapw
    LOGICAL           :: init_vex=.TRUE. !In first call we have to init v_nonlocal
    LOGICAL           :: l_restart=.FALSE.
    LOGICAL           :: l_zref

    REAL              :: bkpt(3)
    REAL, ALLOCATABLE :: eig_irr(:,:)
56 57

    INQUIRE(file="v_x.mat",exist=hybrid%l_addhf)
58
    CALL open_hybrid_io1(DIMENSION,sym%invs)
59 60 61 62 63

    IF (kpts%nkptf==0) THEN
       CALL judft_error("kpoint-set of full BZ not available",&
                        hint="to generate kpts in the full BZ you should specify a k-mesh in inp.xml")
    END IF
64 65
    
    !Check if new non-local potential shall be generated
66 67
    hybrid%l_subvxc = hybrid%l_hybrid.AND.(.NOT.xcpot%is_name("exx"))
    !If this is the first iteration loop we can not calculate a new non-local potential
68 69
    hybrid%l_calhf = (results%last_distance.GE.0.0).AND.(results%last_distance.LT.input%minDistance)
    IF(.NOT.hybrid%l_calhf) THEN
70
       hybrid%l_subvxc = hybrid%l_subvxc.AND.hybrid%l_addhf
71 72
       RETURN
    ENDIF
73

74 75
    results%te_hfex%core = 0

76
    !Check if we are converged well enough to calculate a new potential
77
#if defined(CPP_MPI)&&defined(CPP_NEVER)
78 79 80
    CALL judft_error("Hybrid functionals do not work in parallel version yet")
    CALL MPI_BCAST(results%last_distance .... 
#endif
81 82

    CALL open_hybrid_io1b(DIMENSION,sym%invs)
83
    hybrid%l_addhf = .TRUE.
84

85 86 87
    !In first iteration allocate some memory
    IF (init_vex) THEN
       ALLOCATE(hybrid%ne_eig(kpts%nkpt),hybrid%nbands(kpts%nkpt),hybrid%nobd(kpts%nkptf))
88
       ALLOCATE(hybrid%nbasm(kpts%nkptf))
89
       ALLOCATE(hybrid%div_vv(DIMENSION%neigd,kpts%nkpt,input%jspins))
90 91
       init_vex=.FALSE.
    END IF
92

93
    hybrid%l_subvxc = (hybrid%l_subvxc.AND.hybrid%l_addhf)
94
    IF(.NOT.ALLOCATED(results%w_iks)) ALLOCATE (results%w_iks(DIMENSION%neigd2,kpts%nkpt,input%jspins))
95

96 97
    IF (.NOT.hybrid%l_calhf) RETURN !use existing non-local potential

98 99
    iterHF = iterHF + 1

100 101
    !Delete broyd files
    CALL system("rm broyd*")
102

103 104
    !check if z-reflection trick can be used

105
    l_zref = (sym%zrfs.AND.(SUM(ABS(kpts%bk(3,:kpts%nkpt))).LT.1e-9).AND..NOT.noco%l_noco)
106

107
    CALL timestart("Preparation for Hybrid functionals")
108 109
!    CALL juDFT_WARN ("Hybrid functionals not working in this version")

110 111
    !construct the mixed-basis
    CALL timestart("generation of mixed basis")
112
    CALL mixedbasis(atoms,kpts,dimension,input,cell,sym,xcpot,hybrid,enpara,mpi,v,l_restart)
113 114
    CALL timestop("generation of mixed basis")

115
    CALL open_hybrid_io2(hybrid,DIMENSION,atoms,sym%invs)
116 117 118 119

    CALL timestart("generation of coulomb matrix")
    CALL coulombmatrix(mpi,atoms,kpts,cell,sym,hybrid,xcpot,l_restart)
    CALL timestop("generation of coulomb matrix")
120

121
    CALL hf_init(hybrid,kpts,atoms,input,DIMENSION,hybdat,sym%invs)
122
    CALL timestop("Preparation for Hybrid functionals")
123

124
    CALL timestart("Calculation of non-local HF potential")
125
    DO jsp = 1,input%jspins
126
       CALL HF_setup(hybrid,input,sym,kpts,dimension,atoms,mpi,noco,cell,oneD,results,jsp,enpara,eig_id,&
127
                     hybdat,iterHF,sym%invs,v%mt(:,0,:,:),eig_irr)
128

129
       DO nk = mpi%n_start,kpts%nkpt,mpi%n_stride
130
          CALL lapw%init(input,noco, kpts,atoms,sym,nk,cell,l_zref)
131
  
132
          CALL hsfock(nk,atoms,hybrid,lapw,DIMENSION,kpts,jsp,input,hybdat,eig_irr,sym,cell,&
133
                      noco,results,iterHF,MAXVAL(hybrid%nobd),xcpot,mpi)
134 135
       END DO
    END DO
136
    CALL timestop("Calculation of non-local HF potential")
137
    CALL close_eig(eig_id)
138
  END SUBROUTINE calc_hybrid
139
END MODULE m_calc_hybrid