LCOV - code coverage report
Current view: top level - cdn - qpw_to_nmt.f90 (source / functions) Hit Total Coverage
Test: combined.info Lines: 65 66 98.5 %
Date: 2019-09-08 04:53:50 Functions: 2 2 100.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_qpwtonmt
       8             :   !***************************************************************
       9             :   !     This subroutine calculates the lattice harmonic expansions
      10             :   !     for the plane wave density inside the spheres.      
      11             :   !
      12             :   !             Stefan Bl"ugel  , IFF, Nov. 1997
      13             :   !***************************************************************
      14             : CONTAINS
      15         172 :   SUBROUTINE qpw_to_nmt(&
      16             :        &                      sphhar,atoms,stars,&
      17             :        &                      sym,cell,oneD,mpi,&
      18         172 :        &                      jspin,l_cutoff,qpwc,&
      19         172 :        &                      rho)
      20             :     !
      21             :     USE m_constants
      22             :     USE m_phasy1
      23             :     USE m_sphbes
      24             :     USE m_types
      25             :     USE m_od_phasy
      26             : 
      27             :     IMPLICIT NONE
      28             :     !     ..
      29             :     !     .. Scalar Arguments ..
      30             :     TYPE(t_sphhar),INTENT(IN)   :: sphhar
      31             :     TYPE(t_atoms),INTENT(IN)    :: atoms
      32             :     TYPE(t_stars),INTENT(IN)    :: stars
      33             :     TYPE(t_cell),INTENT(IN)     :: cell
      34             :     TYPE(t_sym),INTENT(IN)      :: sym
      35             :     TYPE(t_oneD),INTENT(IN)     :: oneD
      36             :     TYPE(t_mpi),INTENT(IN)     :: mpi
      37             : 
      38             :     INTEGER, INTENT (IN) :: jspin,l_cutoff    
      39             :     !     ..
      40             :     !     .. Array Arguments ..
      41             :     COMPLEX, INTENT (IN) :: qpwc(stars%ng3)
      42             :     REAL,    INTENT (INOUT) :: rho(:,0:,:,:) !(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins)
      43             :     !-odim
      44             :     !+odim
      45             :     !     ..
      46             :     !     .. Local Scalars ..
      47             :     REAL,    PARAMETER :: zero = 0.0
      48             :     COMPLEX, PARAMETER :: czero = (0.0,0.0)
      49             :     INTEGER in,j,jl,j0,jm,k,l,lh,m,n,nd,nrm,n1,n2,na,lm
      50             :     REAL    d0,gr,r0,rr
      51             :     COMPLEX cp,sm,cprr2
      52             :     !     ..
      53             :     !     .. Local Arrays ..
      54         172 :     COMPLEX pylm( (atoms%lmaxd+1)**2,atoms%ntype ) !,bsl_c(atoms%jmtd,0:lmaxd)
      55         172 :     REAL    bsl(0:atoms%lmaxd),bsl_r(atoms%jmtd,0:atoms%lmaxd),bsl_i(atoms%jmtd,0:atoms%lmaxd)
      56         344 :     INTEGER mr(atoms%ntype),lmx(atoms%ntype),lmxx(atoms%ntype),ntypsy_o(atoms%ntype)
      57             :     !     ..
      58         172 :     !$      REAL,ALLOCATABLE        :: rho_tmp(:,:,:)
      59             : 
      60             : 
      61             :     !----> cut-off l-expansion of non-spherical charge contribution
      62             :     !      from coretails of neighboring atom for l> l_cutoff
      63             :     !
      64         172 :     na = 1
      65         556 :     DO n = 1,atoms%ntype
      66         384 :        lmx(n) = MIN( atoms%lmax(n) , l_cutoff )
      67         384 :        ntypsy_o(n) = atoms%ntypsy(na)
      68         556 :        na = na + atoms%neq(n)
      69             :     END DO
      70             :     !
      71             :     !----> identify atoms with the same radial mesh
      72             :     !
      73             :     j0 = 0
      74             :     r0 = zero
      75             :     d0 = zero
      76             :     nrm= 0
      77         940 :     DO n = 1 , atoms%ntype
      78         384 :        IF (.NOT.(atoms%jri(n).EQ.j0 .AND. atoms%rmsh(1,n).EQ.r0 &
      79             :             &                          .AND. atoms%dx(n).EQ.d0)) THEN   
      80         204 :           j0 = atoms%jri(n)
      81         204 :           r0 = atoms%rmsh(1,n)
      82         204 :           d0 = atoms%dx(n)
      83         204 :           nrm= nrm + 1
      84         204 :           lmxx(nrm) = lmx(n)
      85             :        END IF
      86         384 :        mr(nrm)=n
      87         556 :        lmxx(nrm) = MAX( lmx(n) , lmx(nrm) )
      88             :     END DO
      89             :     !
      90             :     !=====> Loop over the g-vectors
      91             :     !
      92             :     ! ----> g=0 component
      93         172 :     CALL timestart("qpw_to_nmt")
      94         172 :     IF (mpi%irank == 0) THEN
      95          86 :        cp = qpwc(1)*stars%nstr(1)
      96         278 :        DO  n = 1 , atoms%ntype
      97      130362 :           DO j = 1,atoms%jri(n)
      98      130084 :              rr = atoms%rmsh(j,n)*atoms%rmsh(j,n)
      99      130276 :              rho(j,0,n,jspin) = rho(j,0,n,jspin) + rr*sfp_const*cp
     100             :           ENDDO
     101             :        ENDDO
     102             :     ELSE
     103          86 :        rho(:,:,:,jspin) = 0.0
     104             :     ENDIF
     105             : 
     106             :     ! ----> g.ne.0 components
     107             :     !
     108             :     !     g=|=0 terms : \sum_{g =|= 0} 4 \pi i^l \rho_int(g) r_i^{2} \times
     109             :     !                    j_{l} (gr_i) \exp{ig\xi_i} Y^*_{lm} (g)
     110             :     !
     111             :     ! ---->     calculate structure constant for each atom
     112             :     !     (4pi*i**l/nop(3)*sum(R){exp(iRG(taual-taur)*conjg(ylm(RG)) 
     113             :     !
     114             :     !$OMP PARALLEL DEFAULT(SHARED) PRIVATE(k,cp,pylm,n1,in,n2,j)      &
     115             :     !$OMP  PRIVATE(cprr2,gr,bsl,jl,bsl_r,bsl_i,n,nd,lh,l,sm,jm,m,lm)  &
     116         344 :     !$OMP  PRIVATE(rho_tmp)
     117             : 
     118             : 
     119         172 :     !$    ALLOCATE(rho_tmp(size(rho,1),0:size(rho,2)-1,size(rho,3)))
     120         940 :     !$    rho_tmp=0
     121         172 :     !$OMP DO
     122             :     DO k = mpi%irank+2,stars%ng3,mpi%isize
     123      237209 :        cp = qpwc(k)*stars%nstr(k)
     124      237209 :        IF (.NOT.oneD%odi%d1) THEN
     125      237209 :           CALL phasy1(atoms,stars,sym,cell,k,pylm)
     126             :        ELSE
     127             :           !-odim
     128             :           CALL od_phasy(atoms%ntype,stars%ng3,atoms%nat,atoms%lmaxd,atoms%ntype,atoms%neq,atoms%lmax,&
     129           0 :                atoms%taual,cell%bmat,stars%kv3,k,oneD%odi,oneD%ods,pylm) !keep
     130             :           !+odim
     131             :        END IF
     132             :        !
     133             :        n1 = 1
     134      773515 :        DO in = 1 , nrm
     135      268153 :           n2 = mr(in)
     136     5707947 :           bsl_r = 0.0
     137     5707947 :           bsl_i = 0.0
     138   162269170 :           DO j = 1,atoms%jri(n1)
     139   162001017 :              cprr2 = cp*atoms%rmsh(j,n1)*atoms%rmsh(j,n1)
     140   162001017 :              gr = stars%sk3(k)*atoms%rmsh(j,n1)
     141   162001017 :              CALL sphbes(lmxx(in),gr,bsl)
     142  1649488755 :              DO jl=0,lmxx(in)
     143  1487219585 :                 bsl_r(j,jl) = bsl(jl) *  REAL(cprr2)
     144  1649220602 :                 bsl_i(j,jl) = bsl(jl) * AIMAG(cprr2)
     145             :              END DO
     146             :           END DO
     147     1009919 :           DO n = n1,n2 
     148      370883 :              nd = ntypsy_o(n)
     149    16971686 :              DO lh = 0,sphhar%nlh(nd)
     150    16332650 :                 l = sphhar%llh(lh,nd)
     151    16703533 :                 IF ( l.LE.lmx(n) ) THEN
     152    15326911 :                    sm = czero
     153    46777918 :                    DO jm = 1,sphhar%nmem(lh,nd)
     154    31451007 :                       m  = sphhar%mlh(jm,lh,nd)
     155    31451007 :                       lm = l*(l+1) + m + 1
     156    46777918 :                       sm = sm + CONJG(sphhar%clnu(jm,lh,nd)) *pylm(lm,n)
     157             :                    END DO
     158             :                    !$                if (.false.) THEN
     159             :                    rho(:,lh,n,jspin) = rho(:,lh,n,jspin) + bsl_r(:,l) *  REAL(sm)
     160             :                    rho(:,lh,n,jspin) = rho(:,lh,n,jspin) - bsl_i(:,l) * AIMAG(sm)
     161             :                    !$                 endif
     162  9388331894 :                    !$            rho_tmp(:,lh,n)=rho_tmp(:,lh,n)+bsl_r(:,l)*real(sm)
     163    15326911 :                    !$            rho_tmp(:,lh,n)=rho_tmp(:,lh,n)-bsl_i(:,l)*aimag(sm)
     164             :                 END IF
     165             :              END DO
     166             :           END DO
     167      505362 :           n1 = n2 + 1 
     168             :        ENDDO
     169             :     ENDDO
     170             :     !$OMP END DO
     171         344 :     !$OMP CRITICAL
     172         556 :     !$          rho(:,:,:,jspin)=rho(:,:,:,jspin)+rho_tmp
     173             :     !$OMP END CRITICAL
     174         172 :     !$      DEALLOCATE(rho_tmp)
     175             :     !$OMP END PARALLEL
     176             :     !
     177         172 :     CALL timestop("qpw_to_nmt")
     178             :     !
     179         344 :   END SUBROUTINE qpw_to_nmt
     180             : END MODULE m_qpwtonmt

Generated by: LCOV version 1.13