LCOV - code coverage report
Current view: top level - wannier - wann_real.F (source / functions) Hit Total Coverage
Test: combined.info Lines: 0 82 0.0 %
Date: 2019-09-08 04:53:50 Functions: 0 1 0.0 %

          Line data    Source code
       1             : !--------------------------------------------------------------------------------
       2             : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
       3             : ! This file is part of FLEUR and available as free software under the conditions
       4             : ! of the MIT license as expressed in the LICENSE file in more detail.
       5             : !--------------------------------------------------------------------------------
       6             : 
       7             :       MODULE m_wann_real
       8             : c     ********************************************************
       9             : c     calculates the value of the periodic part of the
      10             : c     wavefunction at the given real-grid point p(:)
      11             : c                          Y.Mokrousov 16.8.6
      12             : c     ********************************************************
      13             :       CONTAINS
      14           0 :       SUBROUTINE wann_real(
      15             :      >                  p,n,na,iv,iflag,bkpt,iband,
      16             :      >                  n3d,nmzxyd,n2d,ntypsd,lmaxd,jmtd,
      17           0 :      >                  natd,ntypd,nmzd,nop,nop2,mrot,tau,invtab,
      18             :      >                  nq3,nvac,invs,z1,delz,nmz,nmzxy,nq2,
      19           0 :      >                  lmax,rmsh,jri,pos,ngopr,ntypsy,nvd,
      20           0 :      >                  omtil,amat,bmat,odi,ods,nlod,llod,nlo,llo,
      21           0 :      >                  ff,gg,flo,acof,bcof,ccof,zMat,
      22           0 :      >                  nv,k1,k2,k3,lmd,nbasfcn,l_ss,qss,jspin,addnoco,
      23             :      <                  xdnout)
      24             : 
      25             :       USE m_types
      26             :       USE m_ylm
      27             :       USE m_constants
      28             :       IMPLICIT NONE
      29             : 
      30             :       TYPE(t_mat),INTENT(IN)        :: zMat
      31             : 
      32             : C     .. Scalar Arguments ..
      33             :       INTEGER, INTENT (IN) :: n3d,nmzxyd,n2d,ntypsd,llod,nlod,iband
      34             :       INTEGER, INTENT (IN) :: lmaxd,jmtd,ntypd,natd,nmzd
      35             :       INTEGER, INTENT (IN) :: iflag,n,na,iv,lmd,nv,nvd,nbasfcn
      36             :       INTEGER, INTENT (IN) :: nq3,nvac,nmz,nmzxy,nq2,nop,nop2
      37             :       LOGICAL, INTENT (IN) :: invs,l_ss
      38             :       REAL,    INTENT (IN) :: z1,delz,omtil,bkpt(3),qss(3)
      39             :       INTEGER, INTENT (IN) :: jspin,addnoco
      40             :       COMPLEX, INTENT (OUT):: xdnout
      41             : c-odim
      42             :       TYPE (od_inp), INTENT (IN) :: odi
      43             :       TYPE (od_sym), INTENT (IN) :: ods
      44             : c+odim
      45             : C     ..
      46             : C     .. Array Arguments ..
      47             :       INTEGER, INTENT (IN) :: ngopr(natd),ntypsy(natd),jri(ntypd)
      48             :       INTEGER, INTENT (IN) :: lmax(ntypd),mrot(3,3,nop),invtab(nop)
      49             :       INTEGER, INTENT (IN) :: nlo(ntypd),llo(nlod,ntypd)
      50             :       REAL,    INTENT (IN) :: amat(3,3),bmat(3,3),pos(3,natd)
      51             :       REAL,    INTENT (IN) :: rmsh(jmtd,ntypd),tau(3,nop)
      52             :       INTEGER, INTENT (IN) :: k1(nvd),k2(nvd),k3(nvd) 
      53             :       COMPLEX, INTENT (IN) :: ccof(-llod:llod,nlod,natd)
      54             :       COMPLEX, INTENT (IN) :: acof(0:lmd,natd)
      55             :       COMPLEX, INTENT (IN) :: bcof(0:lmd,natd)
      56             :       REAL,    INTENT (IN) :: ff(ntypd,jmtd,2,0:lmaxd)
      57             :       REAL,    INTENT (IN) :: gg(ntypd,jmtd,2,0:lmaxd)
      58             :       REAL,    INTENT (IN) :: flo(ntypd,jmtd,2,nlod)
      59             :       REAL,    INTENT (INOUT) :: p(3)
      60             : 
      61             : C     ..
      62             : C     .. Local Scalars ..
      63             :       REAL delta,sx,xx1,xx2,rrr,phi,const,arg,tpi,arg1
      64             :       INTEGER i,j,jp3,jr,k,l,nd,nopa,ivac,lm,m,gzi,kk
      65             :       INTEGER kk1,kk2,kk3
      66             :       COMPLEX const2,s,xd1,xd2,const3
      67             : C     ..
      68             : C     .. Local Arrays ..
      69           0 :       COMPLEX sf2(n2d),sf3(n3d),ylm((lmaxd+1)**2)
      70             :       REAL rcc(3),x(3),rcc2(3)
      71             :       REAL bqpt(3)
      72             : 
      73             : 
      74           0 :       tpi   = 2 * pimach()
      75           0 :       const = 1./(sqrt(omtil))
      76             : 
      77             : c..define the factor e^{-ikr}
      78           0 :       rcc2=matmul(bmat,p)/tpi_const
      79             : 
      80           0 :       bqpt = 0.0
      81             : !      if(l_ss.and.jspin.eq.1) then 
      82             : !         bqpt = -qss/2.0
      83             : !      elseif(l_ss.and.jspin.eq.2) then
      84             : !         bqpt = +qss/2.0
      85             : !      endif
      86             : 
      87             :       arg = -tpi*(   (bkpt(1)+bqpt(1))*rcc2(1) 
      88             :      >             + (bkpt(2)+bqpt(2))*rcc2(2)
      89           0 :      >             + (bkpt(3)+bqpt(3))*rcc2(3)  )
      90             : 
      91           0 :       arg1 = tpi*( bkpt(1)*rcc2(1) + bkpt(2)*rcc2(2) + bkpt(3)*rcc2(3) )
      92           0 :       const2 = cmplx(cos(arg),sin(arg))
      93           0 :       const3 = cmplx(cos(arg1),sin(arg1))
      94             : c     write (6,*) 'bkpt,const2,const3=',bkpt(:),const2,const3
      95             : 
      96           0 :       ivac=iv
      97             : 
      98           0 :       IF (iflag.EQ.0) GO TO 20
      99           0 :       IF (iflag.EQ.1) GO TO 40
     100             : c     ---> interstitial part
     101           0 :       rcc=matmul(bmat,p)/tpi_const
     102           0 :       xdnout = cmplx(0.,0.)
     103             : c     write (6,*) 'nv,nvd=',nv,nvd
     104           0 :       IF (zMat%l_real) THEN
     105           0 :          DO k = 1,nv
     106             : c           write (6,*) 'k1,k2,k3=',k1(k),k2(k),k3(k)
     107             : c           write (6,*) 'z(k,iband)=', z(k,iband)
     108           0 :             arg = tpi * ((k1(k))*rcc(1)+(k2(k))*rcc(2)+(k3(k))*rcc(3))
     109             :             xdnout = xdnout + zMat%data_r(k+addnoco,iband)*
     110           0 :      +                        cmplx(cos(arg),sin(arg))*const
     111             :             IF (((abs(p(1)-2.2).le.0.0001).and.(abs(p(2)).le.0.0001))
     112           0 :      &    .or.((abs(p(2)-2.2).le.0.0001).and.(abs(p(1)).le.0.0001)))then
     113             : c              write (6,*) 'p(i)=',p(1:2)
     114             : c              write (6,*) 'G=',k1(k),k2(k),k3(k)
     115             : c              write (6,*) 'z(k,iband)=',z(k,iband)
     116             : c              write (6,*) 'val=',z(k,iband)*cmplx(cos(arg),sin(arg))
     117             :             ENDIF
     118             :          END DO
     119             :       ELSE
     120           0 :          DO k = 1,nv
     121             : c           write (6,*) 'k1,k2,k3=',k1(k),k2(k),k3(k)
     122             : c           write (6,*) 'z(k,iband)=', z(k,iband)
     123           0 :             arg = tpi * ((k1(k))*rcc(1)+(k2(k))*rcc(2)+(k3(k))*rcc(3))
     124             :             xdnout = xdnout + zMat%data_c(k+addnoco,iband)*
     125           0 :      +                        cmplx(cos(arg),sin(arg))*const
     126             :             IF (((abs(p(1)-2.2).le.0.0001).and.(abs(p(2)).le.0.0001))
     127           0 :      &    .or.((abs(p(2)-2.2).le.0.0001).and.(abs(p(1)).le.0.0001)))then
     128             : c              write (6,*) 'p(i)=',p(1:2)
     129             : c              write (6,*) 'G=',k1(k),k2(k),k3(k)
     130             : c              write (6,*) 'z(k,iband)=',z(k,iband)
     131             : c              write (6,*) 'val=',z(k,iband)*cmplx(cos(arg),sin(arg))
     132             :             ENDIF
     133             :          END DO
     134             :       END IF
     135             : c     write (6,*) 'ir:p(i)',p(:)
     136           0 :       RETURN
     137             : c     ---> vacuum part
     138             :    20 CONTINUE
     139           0 :       xdnout = cmplx(0.,0.)
     140           0 :       return
     141             : c     ----> m.t. part
     142             :    40 CONTINUE
     143             :       
     144             : 
     145           0 :       nd = ntypsy(na)
     146           0 :       nopa = ngopr(na)
     147           0 :       nopa=1
     148             :      
     149           0 :       IF (odi%d1) nopa = ods%ngopr(na)
     150           0 :       sx = 0.0
     151           0 :       DO 50 i = 1,3
     152           0 :          x(i) = p(i) - pos(i,na)
     153           0 :          sx = sx + x(i)*x(i)
     154           0 :    50 CONTINUE
     155           0 :       sx = sqrt(sx)
     156           0 :       IF (nopa.NE.1) THEN
     157             : c... switch to internal units
     158           0 :          rcc=matmul(bmat,p)/tpi_const
     159             : c... rotate into representative
     160           0 :          DO 70 i = 1,3
     161           0 :             p(i) = 0.
     162           0 :             DO 60 j = 1,3
     163           0 :               IF (.NOT.odi%d1) THEN
     164           0 :                p(i) = p(i) + mrot(i,j,nopa)*rcc(j)
     165             :               ELSE
     166           0 :                p(i) = p(i) + ods%mrot(i,j,nopa)*rcc(j)
     167             :               END IF
     168           0 :    60       CONTINUE
     169           0 :    70    CONTINUE
     170             : c... switch back to cartesian units
     171           0 :          x=matmul(amat,p)/tpi_const
     172             :       END IF
     173           0 :       DO 80 j = jri(n),2,-1
     174           0 :          IF (sx.GE.rmsh(j,n)) GO TO 90
     175           0 :    80 CONTINUE
     176           0 :    90 jr = j
     177             :       CALL ylm4(
     178             :      >          lmax(n),x,
     179           0 :      <          ylm)
     180           0 :       xd1 = cmplx(0.,0.)
     181           0 :       xd2 = cmplx(0.,0.)
     182           0 :       DO l = 0,lmax(n)
     183             : c        if (p(1).eq.0. .and. p(2).eq.0. .and. p(3).eq.0)then
     184             : c               write (6,*) 'ff(l,300)=',ff(1,300,1,l)
     185             : c               write (6,*) 'ff(l,300)=',ff(1,300,2,l)
     186             : c               write (6,*) 'gg(l,300)=',gg(1,300,1,l)
     187             : c               write (6,*) 'gg(l,300)=',gg(1,300,2,l)
     188             : c        endif
     189           0 :        DO 110 m = -l,l
     190           0 :         lm = l*(l+1)+m
     191           0 :         s = ylm(lm+1)*(ImagUnit)**l
     192             : c       if (p(1).eq.0. .and. p(2).eq.0. .and. p(3).eq.0)then
     193             : c              write (6,*) 'acof=',acof(lm,1)
     194             : c              write (6,*) 'bcof=',bcof(lm,1)
     195             : c       endif
     196             :         xd1 = xd1 + (acof(lm,na)*cmplx(ff(n,jr,1,l),0.)+
     197             :      +               bcof(lm,na)*cmplx(gg(n,jr,1,l),0.))*s/
     198           0 :      /               (rmsh(jr,n)) 
     199             : c    /               (rmsh(jr,n)*rmsh(jr,n))
     200           0 :         IF (jr.EQ.1) GO TO 110
     201             :         xd2 = xd2 + (acof(lm,na)*cmplx(ff(n,jr+1,1,l),0.)+
     202             :      +               bcof(lm,na)*cmplx(gg(n,jr+1,1,l),0.))*s/  
     203           0 :      /               (rmsh(jr+1,n))
     204             : c    /               (rmsh(jr+1,n)*rmsh(jr+1,n))
     205           0 :   110  CONTINUE
     206             :       ENDDO
     207             : c..contributions from the local orbitals
     208           0 :       IF (nlo(n).GE.1) THEN
     209           0 :        DO l = 1,nlo(n)
     210           0 :         DO 111 m = -llo(l,n),llo(l,n)
     211           0 :          lm = llo(l,n)*(llo(l,n)+1)+m
     212           0 :          s = ylm(lm+1)*(ImagUnit)**llo(l,n) 
     213             :          xd1 = xd1 + ccof(m,l,na)*flo(n,jr,1,l)*s/
     214           0 :      /               (rmsh(jr,n))         
     215           0 :          IF (jr.EQ.1) GO TO 111
     216             :          xd2 = xd2 + ccof(m,l,na)*flo(n,jr+1,1,l)*s/
     217           0 :      /               (rmsh(jr+1,n))         
     218           0 :   111   CONTINUE
     219             :        ENDDO
     220             :       ENDIF    
     221           0 :       IF (jr.EQ.1) THEN
     222           0 :          xdnout = xd1
     223             :       ELSE
     224             :          xdnout = xd1 + (xd2-xd1) *
     225           0 :      +                  (sx-rmsh(jr,n)) / (rmsh(jr+1,n)-rmsh(jr,n))
     226             :          
     227             :       END IF
     228           0 :       xdnout = xdnout*const2
     229             : c     write (6,*) 'mt:p(i)',p(:)
     230             :  8000 FORMAT (2f10.6)
     231             : c
     232           0 :       RETURN
     233             :       END SUBROUTINE wann_real
     234             :       END MODULE m_wann_real

Generated by: LCOV version 1.13