types_potden.f90 8.79 KB
Newer Older
Daniel Wortmann's avatar
Daniel Wortmann committed
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_types_potden

  !> Data type for the density or the potential
   TYPE t_potden
     INTEGER             :: iter  
     INTEGER             :: potdenType
12
     COMPLEX,ALLOCATABLE :: pw(:,:),pw_w(:,:)
Daniel Wortmann's avatar
Daniel Wortmann committed
13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
     REAL,ALLOCATABLE    :: mt(:,:,:,:)
     REAL,ALLOCATABLE    :: vacz(:,:,:)
     COMPLEX,ALLOCATABLE :: vacxy(:,:,:,:)
     !For angles of density/potential in noco case
     REAL,ALLOCATABLE  :: theta_pw(:)
     REAL,ALLOCATABLE  :: phi_pw(:)
     REAL,ALLOCATABLE  :: theta_vacz(:,:)
     REAL,ALLOCATABLE  :: phi_vacz(:,:)
     REAL,ALLOCATABLE  :: theta_vacxy(:,:,:)
     REAL,ALLOCATABLE  :: phi_vacxy(:,:,:)
     

     ! For density matrix and associated potential matrix
     COMPLEX, ALLOCATABLE :: mmpMat(:,:,:,:)

     !this type contains two init routines that should be used to allocate
     !memory. You can either specify the datatypes or give the dimensions as integers
     !See implementation below!
   CONTAINS
     PROCEDURE :: init_potden_types
     PROCEDURE :: init_potden_simple
     PROCEDURE :: resetpotden
     GENERIC   :: init=>init_potden_types,init_potden_simple
36
     PROCEDURE :: copy_both_spin
37
     PROCEDURE :: sum_both_spin
38 39
     procedure :: SpinsToChargeAndMagnetisation
     procedure :: ChargeAndMagnetisationToSpins
40 41
     procedure :: addPotDen
     procedure :: subPotDen
Daniel Wortmann's avatar
Daniel Wortmann committed
42 43 44
  END TYPE t_potden

CONTAINS
45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75

  SUBROUTINE sum_both_spin(this,that)
    IMPLICIT NONE
    CLASS(t_potden),INTENT(INOUT)   :: this
    TYPE(t_potden),INTENT(INOUT),OPTIONAL :: that

    IF (PRESENT(that)) THEN
       IF (SIZE(this%pw,2)>1) THEN
          that%mt(:,0:,:,1)=this%mt(:,0:,:,1)+this%mt(:,0:,:,2)
          that%pw(:,1)=this%pw(:,1)+this%pw(:,2)
          that%vacz(:,:,1)=this%vacz(:,:,1)+this%vacz(:,:,2)
          that%vacxy(:,:,:,1)=this%vacxy(:,:,:,1)+this%vacxy(:,:,:,2)
          IF (ALLOCATED(that%pw_w).AND.ALLOCATED(this%pw_w)) that%pw_w(:,1)=this%pw_w(:,1)+this%pw_w(:,2)
       ELSE
          that%mt(:,0:,:,1)=this%mt(:,0:,:,1)
          that%pw(:,1)=this%pw(:,1)
          that%vacz(:,:,1)=this%vacz(:,:,1)
          that%vacxy(:,:,:,1)=this%vacxy(:,:,:,1)
          IF (ALLOCATED(that%pw_w).AND.ALLOCATED(this%pw_w)) that%pw_w(:,1)=this%pw_w(:,1)
       ENDIF
    ELSE
       IF (SIZE(this%pw,2)>1) THEN
          this%mt(:,0:,:,1)=this%mt(:,0:,:,1)+this%mt(:,0:,:,2)
          this%pw(:,1)=this%pw(:,1)+this%pw(:,2)
          this%vacz(:,:,1)=this%vacz(:,:,1)+this%vacz(:,:,2)
          this%vacxy(:,:,:,1)=this%vacxy(:,:,:,1)+this%vacxy(:,:,:,2)
          IF (ALLOCATED(this%pw_w)) this%pw_w(:,1)=this%pw_w(:,1)+this%pw_w(:,2)
       ENDIF
    END IF
  END SUBROUTINE sum_both_spin
    
76 77 78 79 80
  SUBROUTINE copy_both_spin(this,that)
    IMPLICIT NONE
    CLASS(t_potden),INTENT(IN)   :: this
    TYPE(t_potden),INTENT(INOUT) :: that

81 82 83 84 85 86
    that%mt(:,0:,:,1)=this%mt(:,0:,:,1)
    that%pw(:,1)=this%pw(:,1)
    that%vacz(:,:,1)=this%vacz(:,:,1)
    that%vacxy(:,:,:,1)=this%vacxy(:,:,:,1)
    IF (ALLOCATED(that%pw_w).AND.ALLOCATED(this%pw_w)) that%pw_w(:,1)=this%pw_w(:,1)
    
87 88 89 90 91 92 93 94
    IF (SIZE(that%mt,4)==2) THEN
       that%mt(:,0:,:,2)=this%mt(:,0:,:,1)
       that%pw(:,2)=this%pw(:,1)
       that%vacz(:,:,2)=this%vacz(:,:,1)
       that%vacxy(:,:,:,2)=this%vacxy(:,:,:,1)
       IF (ALLOCATED(that%pw_w).AND.ALLOCATED(this%pw_w)) that%pw_w(:,2)=this%pw_w(:,1)
    END IF
  END SUBROUTINE copy_both_spin
95

96
  subroutine SpinsToChargeAndMagnetisation( den )
97
    implicit none
98 99 100 101 102 103 104 105 106 107 108 109 110 111 112
    class(t_potden), intent(inout)    :: den
    !type(t_potden),  intent(inout) :: charge_magn

    type(t_potden) :: copy

    copy = den

    den%mt(:,0:,:,  1) = copy%mt(:,0:,:,  1) + copy%mt(:,0:,:,  2)
    den%mt(:,0:,:,  2) = copy%mt(:,0:,:,  1) - copy%mt(:,0:,:,  2)
    den%pw(:,       1) = copy%pw(:,       1) + copy%pw(:,       2)
    den%pw(:,       2) = copy%pw(:,       1) - copy%pw(:,       2)
    den%vacz(:,:,   1) = copy%vacz(:,:,   1) + copy%vacz(:,:,   2)
    den%vacz(:,:,   2) = copy%vacz(:,:,   1) - copy%vacz(:,:,   2)
    den%vacxy(:,:,:,1) = copy%vacxy(:,:,:,1) + copy%vacxy(:,:,:,2)
    den%vacxy(:,:,:,2) = copy%vacxy(:,:,:,1) - copy%vacxy(:,:,:,2)
113 114 115

  end subroutine

116
  subroutine ChargeAndMagnetisationToSpins( den )
117
    implicit none
118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135
    class(t_potden), intent(inout)    :: den
    !type(t_potden),  intent(inout) :: spins

    type(t_potden) :: copy

    copy = den

    den%mt(:,0:,:,  1) = ( copy%mt(:,0:,:,  1) + copy%mt(:,0:,:,  2) ) / 2
    den%mt(:,0:,:,  2) = ( copy%mt(:,0:,:,  1) - copy%mt(:,0:,:,  2) ) / 2
    den%pw(:,       1) = ( copy%pw(:,       1) + copy%pw(:,       2) ) / 2
    den%pw(:,       2) = ( copy%pw(:,       1) - copy%pw(:,       2) ) / 2
    den%vacz(:,:,   1) = ( copy%vacz(:,:,   1) + copy%vacz(:,:,   2) ) / 2
    den%vacz(:,:,   2) = ( copy%vacz(:,:,   1) - copy%vacz(:,:,   2) ) / 2
    den%vacxy(:,:,:,1) = ( copy%vacxy(:,:,:,1) + copy%vacxy(:,:,:,2) ) / 2
    den%vacxy(:,:,:,2) = ( copy%vacxy(:,:,:,1) - copy%vacxy(:,:,:,2) ) / 2

  end subroutine

136
  subroutine addPotDen( PotDen3, PotDen1, PotDen2 )
137
    implicit none
138 139 140 141 142 143 144 145 146 147 148 149
    class(t_potden), intent(in)    :: PotDen1
    class(t_potden), intent(in)    :: PotDen2
    class(t_potden), intent(inout) :: PotDen3

    PotDen3%iter       = PotDen1%iter
    PotDen3%potdenType = PotDen1%potdenType
    PotDen3%mt         = PotDen1%mt + PotDen2%mt
    PotDen3%pw         = PotDen1%pw + PotDen2%pw
    PotDen3%vacz       = PotDen1%vacz + PotDen2%vacz
    PotDen3%vacxy      = PotDen1%vacxy + PotDen2%vacxy
    if( allocated( PotDen1%pw_w ) .and. allocated( PotDen2%pw_w ) .and. allocated( PotDen3%pw_w ) ) then
      PotDen3%pw_w = PotDen1%pw_w + PotDen2%pw_w
150
    end if
151 152
  
  end subroutine
153

154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169
  subroutine subPotDen( PotDen3, PotDen1, PotDen2 )
    implicit none
    class(t_potden), intent(in)    :: PotDen1
    class(t_potden), intent(in)    :: PotDen2
    class(t_potden), intent(inout) :: PotDen3
 
    PotDen3%iter       = PotDen1%iter
    PotDen3%potdenType = PotDen1%potdenType
    PotDen3%mt         = PotDen1%mt - PotDen2%mt
    PotDen3%pw         = PotDen1%pw - PotDen2%pw
    PotDen3%vacz       = PotDen1%vacz - PotDen2%vacz
    PotDen3%vacxy      = PotDen1%vacxy - PotDen2%vacxy
    if( allocated( PotDen1%pw_w ) .and. allocated( PotDen2%pw_w ) .and. allocated( PotDen3%pw_w ) ) then
      PotDen3%pw_w = PotDen1%pw_w - PotDen2%pw_w
    end if
 
170 171
  end subroutine

172
  SUBROUTINE init_potden_types(pd,stars,atoms,sphhar,vacuum,jspins,nocoExtraDim,potden_type)
Daniel Wortmann's avatar
Daniel Wortmann committed
173
    USE m_judft
174
    USE m_types_setup
Daniel Wortmann's avatar
Daniel Wortmann committed
175
    IMPLICIT NONE
176
    CLASS(t_potden),INTENT(OUT):: pd 
Daniel Wortmann's avatar
Daniel Wortmann committed
177 178 179 180
    TYPE(t_atoms),INTENT(IN) :: atoms
    TYPE(t_stars),INTENT(IN) :: stars
    TYPE(t_sphhar),INTENT(IN):: sphhar
    TYPE(t_vacuum),INTENT(IN):: vacuum
Gregor Michalicek's avatar
Gregor Michalicek committed
181
    INTEGER,INTENT(IN)       :: jspins, potden_type
Daniel Wortmann's avatar
Daniel Wortmann committed
182 183 184
    LOGICAL,INTENT(IN)       :: nocoExtraDim

    CALL init_potden_simple(pd,stars%ng3,atoms%jmtd,sphhar%nlhd,atoms%ntype,&
185
                            atoms%n_u,jspins,nocoExtraDim,potden_type,&
Daniel Wortmann's avatar
Daniel Wortmann committed
186 187 188
                            vacuum%nmzd,vacuum%nmzxyd,stars%ng2)
  END SUBROUTINE init_potden_types

189
  SUBROUTINE init_potden_simple(pd,ng3,jmtd,nlhd,ntype,n_u,jspins,nocoExtraDim,potden_type,nmzd,nmzxyd,n2d)
Daniel Wortmann's avatar
Daniel Wortmann committed
190 191 192 193
    USE m_constants
    USE m_judft
    IMPLICIT NONE
    CLASS(t_potden),INTENT(OUT) :: pd
Gregor Michalicek's avatar
Gregor Michalicek committed
194
    INTEGER,INTENT(IN)          :: ng3,jmtd,nlhd,ntype,n_u,jspins,potden_type
195
    LOGICAL,INTENT(IN)          :: nocoExtraDim
Daniel Wortmann's avatar
Daniel Wortmann committed
196 197 198 199 200 201 202 203 204 205 206 207
    INTEGER,INTENT(IN)          :: nmzd,nmzxyd,n2d

    INTEGER:: err(4)

    err=0
    pd%iter=0
    pd%potdenType=potden_type
    IF(ALLOCATED(pd%pw)) DEALLOCATE (pd%pw)
    IF(ALLOCATED(pd%mt)) DEALLOCATE (pd%mt)
    IF(ALLOCATED(pd%vacz)) DEALLOCATE (pd%vacz)
    IF(ALLOCATED(pd%vacxy)) DEALLOCATE (pd%vacxy)
    IF(ALLOCATED(pd%mmpMat)) DEALLOCATE (pd%mmpMat)
Gregor Michalicek's avatar
Gregor Michalicek committed
208 209 210 211
    ALLOCATE (pd%pw(ng3,MERGE(3,jspins,nocoExtraDim)),stat=err(1))
    ALLOCATE (pd%mt(jmtd,0:nlhd,ntype,jspins),stat=err(2))
    ALLOCATE (pd%vacz(nmzd,2,MERGE(4,jspins,nocoExtraDim)),stat=err(3))
    ALLOCATE (pd%vacxy(nmzxyd,n2d-1,2,MERGE(3,jspins,nocoExtraDim)),stat=err(4))
Daniel Wortmann's avatar
Daniel Wortmann committed
212

Gregor Michalicek's avatar
Gregor Michalicek committed
213
    ALLOCATE (pd%mmpMat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,MAX(1,n_u),jspins))
Daniel Wortmann's avatar
Daniel Wortmann committed
214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233

    IF (ANY(err>0)) CALL judft_error("Not enough memory allocating potential or density")
    pd%pw=CMPLX(0.0,0.0)
    pd%mt=0.0
    pd%vacz=0.0
    pd%vacxy=CMPLX(0.0,0.0)
    pd%mmpMat = CMPLX(0.0,0.0)
  END SUBROUTINE init_potden_simple

  SUBROUTINE resetPotDen(pd)

    IMPLICIT NONE

    CLASS(t_potden),INTENT(INOUT) :: pd

    pd%pw=CMPLX(0.0,0.0)
    pd%mt=0.0
    pd%vacz=0.0
    pd%vacxy=CMPLX(0.0,0.0)
    pd%mmpMat = CMPLX(0.0,0.0)
234
    IF (ALLOCATED(pd%pw_w)) DEALLOCATE(pd%pw_w)
Daniel Wortmann's avatar
Daniel Wortmann committed
235
  END SUBROUTINE resetPotDen
236

Daniel Wortmann's avatar
Daniel Wortmann committed
237
END MODULE m_types_potden