vmts.F90 5.18 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21
MODULE m_vmts
  !     *******************************************************************
  !     this subroutine calculates the lattice harmonics expansion coeffi-*
  !     cients of the coulomb potential for all atom types                *
  !                                c.l.fu, r.podloucky                    *
  !     *******************************************************************
CONTAINS
  SUBROUTINE vmts(mpi,stars,sphhar,atoms,&
       &                sym,cell,oneD,&
       &                vpw,rho,&
       &                vr)

#include"cpp_double.h"
    USE m_constants
    USE m_types
    USE m_intgr, ONLY : intgr2
    USE m_phasy1
    USE m_sphbes
    USE m_od_phasy

    IMPLICIT NONE
22

23 24 25 26 27 28 29 30
    !     .. Scalar Arguments ...
    TYPE(t_mpi),INTENT(IN)     :: mpi
    TYPE(t_stars),INTENT(IN)   :: stars
    TYPE(t_sphhar),INTENT(IN)  :: sphhar
    TYPE(t_atoms),INTENT(IN)   :: atoms
    TYPE(t_sym),INTENT(IN)     :: sym
    TYPE(t_cell),INTENT(IN)    :: cell
    TYPE(t_oneD),INTENT(IN)    :: oneD
31

32
    !     .. Array Arguments ..
33 34 35
    COMPLEX, INTENT (IN) :: vpw(:)!(stars%ng3,input%jspins)
    REAL,    INTENT (IN) :: rho(:,0:,:)!(atoms%jmtd,0:sphhar%nlhd,atoms%ntype)
    REAL,    INTENT (OUT):: vr(:,0:,:)!(atoms%jmtd,0:sphhar%nlhd,atoms%ntype)
36

37 38 39 40
    !     .. Local Scalars ..
    COMPLEX cp,sm
    REAL rmt2l,rmtl,ror,rr,rrlr,fpl21
    INTEGER i,jm,k,l,l21,lh ,n,nd ,lm,n1,nat,m
41

42
    !     .. Local Arrays ..
Daniel Wortmann's avatar
Daniel Wortmann committed
43
    COMPLEX vtl(0:sphhar%nlhd,atoms%ntype)
44
    !$ COMPLEX vtl_loc(0:sphhar%nlhd,atoms%ntype)
Daniel Wortmann's avatar
Daniel Wortmann committed
45
    COMPLEX pylm( (atoms%lmaxd+1)**2 ,atoms%ntype )
46 47 48 49 50 51
    REAL    f1r(atoms%jmtd),f2r(atoms%jmtd),x1r(atoms%jmtd),x2r(atoms%jmtd)
    REAL    sbf(0:atoms%lmaxd),rrl(atoms%jmtd),rrl1(atoms%jmtd)
#ifdef CPP_MPI
    INCLUDE 'mpif.h'
    INTEGER ierr(3)
    COMPLEX, ALLOCATABLE :: c_b(:)
52

53 54 55
    ! ..  External Subroutines
    EXTERNAL MPI_REDUCE
#endif
56
    integer :: OMP_GET_NUM_THREADS, OMP_GET_THREAD_NUM
57 58 59 60

    ! calculate lattice harmonics expansion coefficients of the
    ! interstitial coulomb potential on the sphere boundaries

61
    vtl(:,:) = CMPLX(0.0,0.0)
62 63 64

#ifdef CPP_MPI
    CALL MPI_BCAST(vpw,SIZE(vpw),CPP_MPI_COMPLEX,0,mpi,ierr)
65 66
#endif

67
    ! g=0 component
68 69
    IF (mpi%irank == 0) THEN
       DO n = 1,atoms%ntype
70
          vtl(0,n) = sfp_const*vpw(1)
71 72
       ENDDO
    ENDIF
73 74

    ! g.ne.0 components
75 76 77 78 79 80
    !$OMP PARALLEL DEFAULT(NONE) &
    !$OMP& SHARED(mpi,stars,vpw,oneD,atoms,sym,cell,sphhar,vtl) &
    !$OMP& PRIVATE(k,cp,pylm,nat,n,sbf,nd,lh,sm,jm,m,lm,l)&
    !$OMP& PRIVATE(vtl_loc) 
    !$ vtl_loc(:,:) = CMPLX(0.d0,0.d0)
    !$OMP DO
81
    DO k = mpi%irank+2, stars%ng3, mpi%isize
82
       cp = vpw(k)*stars%nstr(k)
83
       IF (.NOT.oneD%odi%d1) THEN
84
          CALL phasy1(atoms,stars,sym,cell,k,pylm)
85 86
       ELSE
          !-odim
87 88
          CALL od_phasy(atoms%ntype,stars%ng3,atoms%nat,atoms%lmaxd,atoms%ntype,atoms%neq,atoms%lmax,&
                        atoms%taual,cell%bmat,stars%kv3,k,oneD%odi,oneD%ods,pylm)
89 90
          !+odim
       END IF
91
       
92 93 94 95 96 97 98 99 100 101 102 103
       nat = 1
       DO n = 1,atoms%ntype
          CALL sphbes(atoms%lmax(n),stars%sk3(k)*atoms%rmt(n),sbf)
          nd = atoms%ntypsy(nat)
          DO lh = 0,sphhar%nlh(nd)
             l = sphhar%llh(lh,nd)
             sm = (0.,0.)
             DO jm = 1,sphhar%nmem(lh,nd)
                m = sphhar%mlh(jm,lh,nd)
                lm = l*(l+1) + m + 1 
                sm = sm + CONJG(sphhar%clnu(jm,lh,nd))*pylm(lm,n)
             ENDDO
104
             !$ IF (.false.) THEN
105
             vtl(lh,n) = vtl(lh,n) + cp*sbf(l)*sm
106
             !$ ENDIF
107
             !$ vtl_loc(lh,n) = vtl_loc(lh,n) + cp*sbf(l)*sm
108 109 110 111
          ENDDO
          nat = nat + atoms%neq(n)
       ENDDO
    ENDDO
112 113 114 115 116 117
    !$OMP END DO

    !$OMP CRITICAL
    !$ vtl = vtl + vtl_loc
    !$OMP END CRITICAL
    !$OMP END PARALLEL
118
#ifdef CPP_MPI
Daniel Wortmann's avatar
Daniel Wortmann committed
119
    n1 = (sphhar%nlhd+1)*atoms%ntype
120
    ALLOCATE(c_b(n1))
121
    CALL MPI_REDUCE(vtl,c_b,n1,CPP_MPI_COMPLEX,MPI_SUM,0, mpi%mpi_comm,ierr)
122
    IF (mpi%irank.EQ.0) vtl=reshape(c_b,(/sphhar%nlhd+1,atoms%ntype/))
123 124 125 126 127 128 129 130 131 132 133 134 135 136
    DEALLOCATE (c_b)
#endif

    !     ----> solution of the poisson's equation
    nat = 1
    DO  n = 1,atoms%ntype
       nd = atoms%ntypsy(nat)
       DO  lh = 0,sphhar%nlh(nd)
          l = sphhar%llh(lh,nd)
          l21 = 2*l + 1
          fpl21 = fpi_const/l21
          DO i = 1,atoms%jri(n)
             rrl(i) = atoms%rmsh(i,n)**l
             rrl1(i) = 1./( rrl(i) * atoms%rmsh(i,n) )
137 138
             x1r(i) = rrl(i)*rho(i,lh,n)
             x2r(i) = rrl1(i)*rho(i,lh,n)
139 140 141 142 143 144 145 146
          ENDDO
          CALL intgr2(x1r,atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),f1r)
          CALL intgr2(x2r,atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),f2r)
          rmtl = 1./atoms%rmt(n)**l
          rmt2l = 1./atoms%rmt(n)**l21
          DO  i = 1,atoms%jri(n)
             rrlr = rrl(i)*rmt2l
             ror = rrl(i)*rmtl
147
             vr(i,lh,n) = fpl21 * (rrl1(i)*f1r(i)-rrlr*f1r(atoms%jri(n))+&
148
                                   rrl(i) * (f2r(atoms%jri(n))-f2r(i))) + ror*vtl(lh,n)
149 150 151 152 153 154 155
          ENDDO
       ENDDO
       nat = nat + atoms%neq(n)
    ENDDO
    DO  n = 1,atoms%ntype
       DO  i = 1,atoms%jri(n)
          rr = atoms%rmsh(i,n)/atoms%rmt(n)
156
          vr(i,0,n) = vr(i,0,n)-sfp_const*(1.-rr)/atoms%rmsh(i,n)*atoms%zatom(n)
157 158 159 160
       ENDDO
    ENDDO
  END SUBROUTINE vmts
END MODULE m_vmts