vgen_xcpot.F90 7.92 KB
Newer Older
1 2 3 4 5
!--------------------------------------------------------------------------------
! 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.
!--------------------------------------------------------------------------------
6
MODULE m_vgen_xcpot
7 8
  USE m_juDFT
CONTAINS
9 10
  SUBROUTINE vgen_xcpot(hybrid,input,xcpot,DIMENSION, atoms,sphhar,stars,&
       vacuum,sym, obsolete,cell,oneD,sliceplot,mpi,noco,den,denRot,vTot,vx,results)
11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
    !     ***********************************************************
    !     FLAPW potential generator                           *
    !     ***********************************************************
    !     calculates the density-potential integrals needed for the
    !     total energy
    !     TE_VCOUL  :   charge density-coulomb potential integral
    !     TE_VEFF:   charge density-effective potential integral
    !     TE_EXC :   charge density-ex-corr.energy density integral
    !     ***********************************************************
    USE m_constants
    USE m_intnv
    USE m_vmtxcg
    USE m_vmtxc
    USE m_vvacxc
    USE m_vvacxcg
    USE m_visxc
    USE m_visxcg
    USE m_checkdopall
29
    USE m_cdn_io 
30
    USE m_types
Daniel Wortmann's avatar
Daniel Wortmann committed
31
   
32
    IMPLICIT NONE
33
    CLASS(t_xcpot),INTENT(IN)       :: xcpot
34 35 36 37 38 39
    TYPE(t_hybrid),INTENT(IN)       :: hybrid
    TYPE(t_mpi),INTENT(IN)          :: mpi
    TYPE(t_dimension),INTENT(IN)    :: dimension
    TYPE(t_oneD),INTENT(IN)         :: oneD
    TYPE(t_obsolete),INTENT(IN)     :: obsolete
    TYPE(t_sliceplot),INTENT(IN)    :: sliceplot
40
    TYPE(t_input),INTENT(IN)        :: input  
41 42 43 44 45 46
    TYPE(t_vacuum),INTENT(IN)       :: vacuum
    TYPE(t_noco),INTENT(IN)         :: noco
    TYPE(t_sym),INTENT(IN)          :: sym
    TYPE(t_stars),INTENT(IN)        :: stars
    TYPE(t_cell),INTENT(IN)         :: cell
    TYPE(t_sphhar),INTENT(IN)       :: sphhar
47
    TYPE(t_atoms),INTENT(IN)        :: atoms 
48
    TYPE(t_potden), INTENT(IN)      :: den, denRot
49 50
    TYPE(t_potden),INTENT(INOUT)    :: vTot,vx
    TYPE(t_results),INTENT(INOUT),OPTIONAL   :: results
51 52 53 54
    !     ..


    !     .. Local type instances ..
55
    TYPE(t_potden)    :: workDen,exc,veff
56 57
    !     .. Local Scalars ..
    INTEGER i,i3,irec2,irec3,ivac,j,js,k,k3,lh,n,nzst1
58 59
    INTEGER ifftd,ifftd2, ifftxc3d
    INTEGER jsp,l
60 61 62 63 64
#ifdef CPP_MPI
    include 'mpif.h'
    integer:: ierr
#endif

65
    CALL exc%init_potden_types(stars,atoms,sphhar,vacuum,1,.false.,1) !one spin only
Daniel Wortmann's avatar
Daniel Wortmann committed
66 67 68
    ALLOCATE(exc%pw_w(stars%ng3,1))
    IF (PRESENT(results)) THEN
       CALL veff%init(stars,atoms,sphhar,vacuum,input%jspins,.FALSE.,1)
69
       ALLOCATE(veff%pw_w,mold=veff%pw)
Daniel Wortmann's avatar
Daniel Wortmann committed
70
    ENDIF
71

72
    !     ******** exchange correlation potential******************
73 74


75
    !     ---> vacuum region
76 77 78 79 80 81 82 83 84 85 86 87
    IF (mpi%irank == 0) THEN
       IF (input%film) THEN

          CALL timestart("Vxc in vacuum")

          ifftd2 = 9*stars%mx1*stars%mx2
          IF (oneD%odi%d1) ifftd2 = 9*stars%mx3*oneD%odi%M

          IF (.NOT.xcpot%is_gga()) THEN  ! LDA

             IF (.NOT.oneD%odi%d1) THEN

88 89
                CALL vvacxc(ifftd2,stars,vacuum,xcpot,input,noco,Den,&
                     vTot,exc)
90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107

             ELSE
                CALL judft_error("OneD broken")
                !           CALL vvacxc(&
                !     &                 stars,oneD%M,vacuum,odi%n2d,dimension,ifftd2,&
                !     &                 xcpot,input,odi%nq2,&
                !     &                 odi%nst2,den,noco,&
                !     &                 odi%kimax2%igf,odl%pgf,&
                !     &                 vTot%vacxy,vTot%vacz,&
                !     &                 excxy,excz)

             ENDIF

          ELSE      ! GGA

             IF (oneD%odi%d1) THEN
                CALL judft_error("OneD broken")

108 109 110
                !CALL vvacxcg(ifftd2,stars,vacuum,noco,oneD,&
                !     cell,xcpot,input,obsolete,workDen, ichsmrg,&
                !     vTot%vacxy,vTot%vacz,rhmn, exc%vacxy,exc%vacz)
111 112 113

             ELSE
                CALL vvacxcg(ifftd2,stars,vacuum,noco,oneD,&
114
                     cell,xcpot,input,obsolete,Den, &
115
                     vTot, exc)
116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136

             END IF

          END IF
          CALL timestop("Vxc in vacuum")
          !+odim
       END IF
       !     ----------------------------------------
       !     ---> interstitial region
       CALL timestart("Vxc in interstitial")

       ifftd=27*stars%mx1*stars%mx2*stars%mx3

       IF ( (.NOT. obsolete%lwb) .OR. ( .not.xcpot%is_gga() ) ) THEN
          ! no White-Bird-trick

          ifftxc3d = stars%kxc1_fft*stars%kxc2_fft*stars%kxc3_fft

          IF ( .NOT.xcpot%is_gga() ) THEN
             ! LDA

137
             CALL visxc(ifftd,stars,noco,xcpot,input,den,vTot,vx,exc)
138 139 140

          ELSE ! GGA

141
             CALL visxcg(ifftd,stars,sym,ifftxc3d,cell,den,xcpot,input,&
142
                  obsolete,noco,vTot,vx,exc)
143 144 145 146 147 148 149

          END IF

       ELSE
          ! White-Bird-trick

          WRITE(6,'(a)') "W+B trick cancelled out. visxcwb uses at present common block cpgft3.",&
150
               "visxcwb needs to be reprogrammed according to visxcg.f"
151
          CALL juDFT_error("visxcwb",calledby ="vgen")
152

153 154 155 156 157

       END IF
       !


158

159 160 161
       CALL timestop("Vxc in interstitial")

    ENDIF !irank==0
162 163 164
    !
    !     ------------------------------------------
    !     ----> muffin tin spheres region
165 166 167 168 169

    IF (mpi%irank == 0) THEN
       CALL timestart ("Vxc in MT")
    END IF
#ifdef CPP_MPI
170
    CALL MPI_BCAST(den%mt,atoms%jmtd*(1+sphhar%nlhd)*atoms%ntype*dimension%jspd,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)
171
#endif
172 173 174 175 176 177
    IF (xcpot%is_gga()) THEN
       CALL vmtxcg(dimension,mpi,sphhar,atoms, den,xcpot,input,sym,&
            obsolete, vTot,vx,exc)
    ELSE
       CALL vmtxc(DIMENSION,sphhar,atoms, den,xcpot,input,sym, vTot,exc,vx)
    ENDIF
178 179


180 181 182
    !
    ! add MT EXX potential to vr
    !
183 184 185 186 187 188 189 190 191

    IF (mpi%irank == 0) THEN

       CALL timestop ("Vxc in MT")
       !     ------------------------------------------
       !     ---> check continuity of total potential

       IF (input%vchk) THEN
          CALL checkDOPAll(input,dimension,sphhar,stars,atoms,sym,vacuum,oneD,&
192
               cell,vTot,1)
193 194 195 196 197
       END IF

       !
       !============TOTAL======================================
       !
198
       IF (PRESENT(results)) THEN
199 200 201 202
          !
          !     CALCULATE THE INTEGRAL OF n1*Veff1 + n2*Veff2
          !     Veff = Vcoulomb + Vxc
          !
203 204
          IF (noco%l_noco) THEN 
             workDen = denRot
205
          ELSE
206 207
             workden=den
          ENDIF
208

209 210 211 212 213
          veff=vTot
          IF( xcpot%is_hybrid() ) THEN
             veff%pw   = vTot%pw   -xcpot%get_exchange_weight()  * vx%pw
             veff%pw_w = vTot%pw_w - xcpot%get_exchange_weight() * vx%pw_w
             veff%mt   = vTot%mt   - xcpot%get_exchange_weight() * vx%mt
214 215 216
          END IF

          results%te_veff = 0.0
217
          DO  js = 1,input%jspins
218 219
             WRITE (6,FMT=8050) js
8050         FORMAT (/,10x,'density-effective potential integrals for spin ',i2,/)
220 221 222 223
             CALL int_nv(js,stars,vacuum,atoms,sphhar, cell,sym,input,oneD,veff,workden,results%te_veff)
          ENDDO
          WRITE (6,FMT=8060) results%te_veff
8060      FORMAT (/,10x,'total density-effective potential integral :', t40,f20.10)
224
          !
225
          !     CALCULATE THE INTEGRAL OF n*exc
226
          !
227 228 229
          !     ---> perform spin summation of charge densities
          !     ---> for the calculation of Exc
          CALL  workden%sum_both_spin()
230

231 232
          WRITE (6,FMT=8070)
8070      FORMAT (/,10x,'charge density-energy density integrals',/)
233

234 235 236 237
          results%te_exc = 0.0
          CALL int_nv(1,stars,vacuum,atoms,sphhar, cell,sym,input,oneD,&
               exc,workDen, results%te_exc)
          WRITE (6,FMT=8080) results%te_exc
238

239
8080      FORMAT (/,10x,'total charge density-energy density integral :', t40,f20.10)
240

241
       END IF
242

243
    ENDIF ! mpi%irank == 0
244

245 246
  END SUBROUTINE vgen_xcpot
END MODULE m_vgen_xcpot