od_abvac.f90 4.86 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 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 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 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
MODULE m_od_abvac
CONTAINS
  SUBROUTINE od_abvac(&
       &     cell,vacuum,DIMENSION,stars,&
       &     oneD,qssbti,&
       &     n2d_1,&
       &     wronk,evac,bkpt,MM,vM,&
       &     vz,kvac3,nv2,&
       &     uz,duz,u,udz,dudz,ddnv,ud)
    !**************************************************************
    !      determines the nesessary values and derivatives on the 
    !      vacuum cylindrical boundary for finding a and b coefficients
    !      for the construcing vacuum charge density in vacden.F
    !                          Y.Mokrousov, 7th of october 2002
    !*************************************************************** 
    USE m_vacuz
    USE m_vacudz
    USE m_types
    IMPLICIT NONE
    TYPE(t_dimension),INTENT(IN)   :: DIMENSION
    TYPE(t_oneD),INTENT(IN)        :: oneD
    TYPE(t_vacuum),INTENT(IN)      :: vacuum
    TYPE(t_stars),INTENT(IN)       :: stars
    TYPE(t_cell),INTENT(IN)        :: cell

    !     .. scalar Arguments..
    REAL wronk
    INTEGER, INTENT (in) :: MM ,vM
    INTEGER, INTENT (in) :: n2d_1,nv2
    REAL,    INTENT (in) :: evac
    !     ..array arguments..

    INTEGER, INTENT (in) :: kvac3(DIMENSION%nv2d)
    REAL,    INTENT (in) :: bkpt(3),qssbti 
    REAL,    INTENT (in) :: vz(vacuum%nmzd) 
    REAL,    INTENT (out):: udz(DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: uz(DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: dudz(DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: duz(DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: u(vacuum%nmzd,DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: ud(vacuum%nmzd,DIMENSION%nv2d,-vM:vM)
    REAL,    INTENT (out):: ddnv(DIMENSION%nv2d,-vM:vM)
    !     ..local scalars..
    REAL ev,scale,xv,yv,vzero,v1
    INTEGER i,ik,jk,jspin,jsp1,jsp2 ,l,m
    INTEGER i1,i2,i3,ind1,ind3
    !     .. local arrays..
    REAL wdz(DIMENSION%nv2d,-vM:vM),wz(DIMENSION%nv2d,-vM:vM)
    REAL dwdz(DIMENSION%nv2d,-vM:vM),dwz(DIMENSION%nv2d,-vM:vM)
    REAL v(3),x(vacuum%nmzd)
    REAL  vr0(vacuum%nmzd)
    REAL w(vacuum%nmzd,DIMENSION%nv2d,-vM:vM),wd(vacuum%nmzd,DIMENSION%nv2d,-vM:vM)

    !     wronksian for the schrodinger equation given by an identity

    wronk = 2.0

    DO  ik = 1,nv2
       DO  m = 0,vM
          v(1) = 0.0
          v(2) = 0.0
          v(3) = bkpt(3) + kvac3(ik) + qssbti
          ev = evac - 0.5*DOT_PRODUCT(v,MATMUL(v,cell%bbmat))

          !     constructing of the 'pseudopotential'

          DO  i=1,vacuum%nmz
             v1 = 1./(8.*((cell%z1+(i-1)*vacuum%delz)**2))&
                  &              -(m*m)/(2.*((cell%z1+(i-1)*vacuum%delz)**2))
             vr0(i) = vz(i)-v1
          ENDDO
          vzero = vr0(vacuum%nmz)

          !     obtaining solutions with the 'pseudopotential'

          CALL vacuz(ev,vr0(1),vzero,vacuum%nmz,vacuum%delz,&
                          wz(ik,m),dwz(ik,m),w(1,ik,m))
          CALL vacudz(ev,vr0(1),vzero,vacuum%nmz,vacuum%delz,&
               wdz(ik,m),dwdz(ik,m),ddnv(ik,m), wd(1,ik,m),dwz(ik,m),w(1,ik,m))

          scale = wronk/(wdz(ik,m)*dwz(ik,m)-&
               &           dwdz(ik,m)*wz(ik,m))
          wdz(ik,m) = scale*wdz(ik,m)
          dwdz(ik,m) = scale*dwdz(ik,m)
          ddnv(ik,m) = scale*ddnv(ik,m)
          IF (m.GT.0) THEN
             wdz(ik,-m) = wdz(ik,m)
             dwdz(ik,-m) = dwdz(ik,m)
             ddnv(ik,-m) = ddnv(ik,m)
          END IF
          DO  i = 1,vacuum%nmz
             wd(i,ik,m) = scale*wd(i,ik,m)
             w(i,ik,m) = scale*w(i,ik,m)
             IF (m.GT.0) THEN
                wd(i,ik,-m) = wd(i,ik,m)
                w(i,ik,-m) = w(i,ik,m)
             END IF
          ENDDO
          !     constructing 'real' solutions

          DO  i=1,vacuum%nmz
             u(i,ik,m)=w(i,ik,m)/SQRT(cell%z1+(i-1)*vacuum%delz)
             ud(i,ik,m)=wd(i,ik,m)/SQRT(cell%z1+(i-1)*vacuum%delz)
             IF (m.GT.0) THEN
                u(i,ik,-m) = u(i,ik,m)
                ud(i,ik,-m) = ud(i,ik,m)
             END IF
          ENDDO
          duz(ik,m)=(-dwz(ik,m))/SQRT(cell%z1)-&
               &           wz(ik,m)/(2.0*((cell%z1)**(1.5)))
          uz(ik,m)=wz(ik,m)/SQRT(cell%z1)
          dudz(ik,m)=(-dwdz(ik,m))/SQRT(cell%z1)-&
               &           wdz(ik,m)/(2.0*((cell%z1)**(1.5)))
          udz(ik,m)=wdz(ik,m)/SQRT(cell%z1)
          IF (m.GT.0) THEN
             duz(ik,-m) = duz(ik,m)
             uz(ik,-m) = uz(ik,m)
             dudz(ik,-m) = dudz(ik,m)
             udz(ik,-m) = udz(ik,m)
          END IF

       ENDDO
    ENDDO
    RETURN
  END SUBROUTINE od_abvac
END MODULE m_od_abvac