LCOV - code coverage report
Current view: top level - cdn - cdntot.f90 (source / functions) Hit Total Coverage
Test: combined.info Lines: 62 96 64.6 %
Date: 2019-09-08 04:53:50 Functions: 2 6 33.3 %

          Line data    Source code
       1             : MODULE m_cdntot
       2             : !     ********************************************************
       3             : !     calculate the total charge density in the interstial.,
       4             : !     vacuum, and mt regions      c.l.fu
       5             : !     ********************************************************
       6             : CONTAINS
       7         367 :    SUBROUTINE integrate_cdn(stars,atoms,sym,vacuum,input,cell,oneD, integrand, &
       8         367 :                                    q, qis, qmt, qvac, qtot, qistot)
       9             :       USE m_intgr, ONLY : intgr3
      10             :       USE m_constants
      11             :       USE m_qsf
      12             :       USE m_pwint
      13             :       USE m_types
      14             :       USE m_juDFT
      15             :       USE m_convol
      16             :       IMPLICIT NONE
      17             :       TYPE(t_stars),INTENT(IN)  :: stars
      18             :       TYPE(t_atoms),INTENT(IN)  :: atoms
      19             :       TYPE(t_sym),INTENT(IN)    :: sym
      20             :       TYPE(t_vacuum),INTENT(IN) :: vacuum
      21             :       TYPE(t_input),INTENT(IN)  :: input
      22             :       TYPE(t_cell),INTENT(IN)   :: cell
      23             :       TYPE(t_oneD),INTENT(IN)   :: oneD
      24             :       TYPE(t_potden),INTENT(IN) :: integrand
      25             :       REAL, INTENT(out)         :: q(input%jspins), qis(input%jspins), qmt(atoms%ntype,input%jspins),&
      26             :                                    qvac(2,input%jspins), qtot, qistot
      27             :       INTEGER                   :: jsp, j, ivac, nz, n
      28         734 :       REAL                      :: q2(vacuum%nmz), w, rht1(vacuum%nmzd,2,input%jspins)
      29         734 :       COMPLEX                   :: x(stars%ng3)
      30             :       
      31         367 :       qtot = 0.0
      32         367 :       qistot = 0.0
      33        1021 :       DO jsp = 1,input%jspins
      34         654 :          q(jsp) = 0.0
      35             : !     -----mt charge
      36         654 :          CALL timestart("MT")
      37        1952 :          DO n = 1,atoms%ntype
      38        1298 :             CALL intgr3(integrand%mt(:,0,n,jsp),atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),w)
      39        1298 :             qmt(n, jsp) = w*sfp_const
      40        1952 :             q(jsp) = q(jsp) + atoms%neq(n)*qmt(n,jsp)
      41             :          ENDDO
      42         654 :          CALL timestop("MT")
      43             : !     -----vacuum region
      44         654 :          IF (input%film) THEN
      45          46 :             DO ivac = 1,vacuum%nvac
      46       11523 :                DO nz = 1,vacuum%nmz
      47        5773 :                   IF (oneD%odi%d1) THEN
      48             :                      rht1(nz,ivac,jsp) = (cell%z1+(nz-1)*vacuum%delz)*&
      49           0 :                                            integrand%vacz(nz,ivac,jsp)
      50             :                   ELSE
      51        5750 :                      rht1(nz,ivac,jsp) =  integrand%vacz(nz,ivac,jsp)
      52             :                   END IF
      53             :                END DO
      54          23 :                CALL qsf(vacuum%delz,rht1(1,ivac,jsp),q2,vacuum%nmz,0)
      55          23 :                qvac(ivac,jsp) = q2(1)*cell%area
      56          46 :                IF (.NOT.oneD%odi%d1) THEN
      57          23 :                   q(jsp) = q(jsp) + qvac(ivac,jsp)*2./real(vacuum%nvac)
      58             :                ELSE
      59           0 :                   q(jsp) = q(jsp) + cell%area*q2(1)
      60             :                END IF
      61             :             ENDDO
      62             :          END IF
      63             : !     -----is region
      64         654 :          qis(jsp) = 0.
      65             : 
      66         654 :          CALL pwint_all(stars,atoms,sym,oneD,cell,1,stars%ng3,x)
      67      788919 :          DO j = 1,stars%ng3
      68      788919 :             qis(jsp) = qis(jsp) + integrand%pw(j,jsp)*x(j)*stars%nstr(j)
      69             :          ENDDO
      70             : 
      71         654 :          qistot = qistot + qis(jsp)
      72         654 :          q(jsp) = q(jsp) + qis(jsp)
      73        1021 :          qtot = qtot + q(jsp)
      74             :       END DO ! loop over spins
      75         367 :    END SUBROUTINE integrate_cdn
      76             : 
      77           0 :    SUBROUTINE integrate_realspace(xcpot, atoms, sym, sphhar, input, &
      78           0 :                                   stars, cell, oneD, vacuum, noco, mt, is, hint)
      79             :       use m_types
      80             :       use m_mt_tofrom_grid
      81             :       use m_pw_tofrom_grid
      82             :       use m_constants
      83             :       implicit none
      84             :       CLASS(t_xcpot), INTENT(inout)   :: xcpot
      85             :       TYPE(t_atoms),INTENT(IN)      :: atoms
      86             :       TYPE(t_sym), INTENT(in)       :: sym
      87             :       TYPE(t_sphhar), INTENT(IN)    :: sphhar
      88             :       TYPE(t_input), INTENT(IN)     :: input
      89             :       TYPE(t_stars), INTENT(IN)     :: stars
      90             :       TYPE(t_cell), INTENT(IN)      :: cell
      91             :       TYPE(t_oneD), INTENT(in)      :: oneD
      92             :       TYPE(t_vacuum), INTENT(in)    :: vacuum
      93             :       TYPE(t_noco), INTENT(in)      :: noco
      94             :       real, intent(inout)           :: mt(:,:,:), is(:,:)
      95             :       character(len=*), intent(in), optional :: hint
      96             :       integer                       :: n_atm, i
      97             : 
      98           0 :       TYPE(t_potden)                :: tmp_potden
      99           0 :       REAL                          :: q(input%jspins), qis(input%jspins), &
     100           0 :                                        qmt(atoms%ntype,input%jspins), qvac(2,input%jspins),&
     101             :                                        qtot, qistot
     102             : 
     103           0 :       call tmp_potden%init(stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_DEN)
     104           0 :       call init_mt_grid(input%jspins, atoms, sphhar, xcpot, sym)
     105           0 :       do n_atm =1,atoms%ntype
     106             :          call mt_from_grid(atoms, sphhar, n_atm, input%jspins, mt(:,:,n_atm), &
     107           0 :                            tmp_potden%mt(:,0:,n_atm,:))
     108             : 
     109           0 :          do i=1,atoms%jri(n_atm)
     110           0 :             tmp_potden%mt(i,:,n_atm,:) = tmp_potden%mt(i,:,n_atm,:) * atoms%rmsh(i,n_atm)**2
     111             :          enddo
     112             :       enddo
     113           0 :       call finish_mt_grid()
     114             : 
     115           0 :       call init_pw_grid(xcpot, stars, sym, cell)
     116           0 :       call pw_from_grid(xcpot, stars, .False., is, tmp_potden%pw)
     117           0 :       call finish_pw_grid()
     118             : 
     119             :       call integrate_cdn(stars,atoms,sym,vacuum,input,cell,oneD, tmp_potden, &
     120           0 :                                    q, qis, qmt, qvac, qtot, qistot)
     121             : 
     122           0 :       call print_cdn_inte(q, qis, qmt, qvac, qtot, qistot, hint)
     123           0 :    END SUBROUTINE integrate_realspace
     124             : 
     125         367 :    SUBROUTINE cdntot(stars,atoms,sym,vacuum,input,cell,oneD,&
     126             :                      den,l_printData,qtot,qistot)
     127             : 
     128             :       USE m_intgr, ONLY : intgr3
     129             :       USE m_constants
     130             :       USE m_qsf
     131             :       USE m_pwint
     132             :       USE m_types
     133             :       USE m_juDFT
     134             :       USE m_convol
     135             :       USE m_xmlOutput
     136             :       IMPLICIT NONE
     137             : 
     138             : !     .. Scalar Arguments ..
     139             :       TYPE(t_stars),INTENT(IN)  :: stars
     140             :       TYPE(t_atoms),INTENT(IN)  :: atoms
     141             :       TYPE(t_sym),INTENT(IN)    :: sym
     142             :       TYPE(t_vacuum),INTENT(IN) :: vacuum
     143             :       TYPE(t_input),INTENT(IN)  :: input
     144             :       TYPE(t_oneD),INTENT(IN)   :: oneD
     145             :       TYPE(t_cell),INTENT(IN)   :: cell
     146             :       TYPE(t_potden),INTENT(IN) :: den
     147             :       LOGICAL,INTENT(IN)        :: l_printData
     148             :       REAL,INTENT(OUT)          :: qtot,qistot
     149             : 
     150             : !     .. Local Scalars ..
     151             :       COMPLEX x(stars%ng3)
     152         734 :       REAL q(input%jspins),qis(input%jspins),w,mtCharge
     153             :       INTEGER i,ivac,j,jsp,n,nz
     154             : !     ..
     155             : !     .. Local Arrays ..
     156         734 :       REAL qmt(atoms%ntype,input%jspins),qvac(2,input%jspins)
     157         367 :       INTEGER, ALLOCATABLE :: lengths(:,:)
     158             :       CHARACTER(LEN=20) :: attributes(6), names(6)
     159             :       
     160         367 :       CALL timestart("cdntot")
     161             :       call integrate_cdn(stars,atoms,sym,vacuum,input,cell,oneD, den, &
     162         367 :                                    q, qis, qmt, qvac, qtot, qistot)
     163             :  
     164         367 :       IF (input%film) THEN
     165          14 :          ALLOCATE(lengths(4+vacuum%nvac,2))
     166             :       ELSE
     167         353 :          ALLOCATE(lengths(4,2))
     168             :       END IF
     169             : 
     170        1021 :       DO jsp = 1,input%jspins
     171         654 :          WRITE (6,FMT=8000) jsp,q(jsp),qis(jsp), (qmt(n,jsp),n=1,atoms%ntype)
     172         654 :          IF (input%film) WRITE (6,FMT=8010) (i,qvac(i,jsp),i=1,vacuum%nvac)
     173         654 :          mtCharge = SUM(qmt(1:atoms%ntype,jsp) * atoms%neq(1:atoms%ntype))
     174         654 :          names(1) = 'spin'         ; WRITE(attributes(1),'(i0)')    jsp      ; lengths(1,1)=4  ; lengths(1,2)=1
     175         654 :          names(2) = 'total'        ; WRITE(attributes(2),'(f14.7)') q(jsp)   ; lengths(2,1)=5  ; lengths(2,2)=14
     176         654 :          names(3) = 'interstitial' ; WRITE(attributes(3),'(f14.7)') qis(jsp) ; lengths(3,1)=12 ; lengths(3,2)=14
     177         654 :          names(4) = 'mtSpheres'    ; WRITE(attributes(4),'(f14.7)') mtCharge ; lengths(4,1)=9  ; lengths(4,2)=14
     178        1021 :          IF(l_printData) THEN
     179         648 :             IF(input%film) THEN
     180          46 :                DO i = 1, vacuum%nvac
     181          23 :                   WRITE(names(4+i),'(a6,i0)') 'vacuum', i
     182          23 :                   WRITE(attributes(4+i),'(f14.7)') qvac(i,jsp)
     183          23 :                   lengths(4+i,1)=7
     184          46 :                   lengths(4+i,2)=14
     185             :                END DO
     186             :                CALL writeXMLElementFormPoly('spinDependentCharge',names(1:4+vacuum%nvac),&
     187          23 :                                             attributes(1:4+vacuum%nvac),lengths)
     188             :             ELSE
     189         625 :                CALL writeXMLElementFormPoly('spinDependentCharge',names(1:4),attributes(1:4),lengths)
     190             :             END IF
     191             :          END IF
     192             :       END DO ! loop over spins
     193         367 :       WRITE (6,FMT=8020) qtot
     194         367 :       IF(l_printData) THEN
     195         364 :          CALL writeXMLElementFormPoly('totalCharge',(/'value'/),(/qtot/),reshape((/5,20/),(/1,2/)))
     196             :       END IF
     197             : 8000  FORMAT (/,10x,'total charge for spin',i3,'=',f12.6,/,10x,&
     198             :                'interst. charge =   ',f12.6,/,&
     199             :                (10x,'mt charge=          ',4f12.6,/))
     200             : 8010  FORMAT (10x,'vacuum ',i2,'  charge=  ',f12.6)
     201             : 8020  FORMAT (/,10x,'total charge  =',f12.6)
     202             : 
     203         367 :       CALL timestop("cdntot")
     204         367 :    END SUBROUTINE cdntot
     205             : 
     206           0 :    SUBROUTINE print_cdn_inte(q, qis, qmt, qvac, qtot, qistot, hint)
     207             :       use  ieee_arithmetic
     208             :       implicit none
     209             :       REAL, INTENT(in)                       :: q(:), qis(:), qmt(:,:), qvac(:,:), qtot, qistot
     210             :       character(len=*), intent(in), optional :: hint
     211             :       integer                                :: n_mt
     212             :       
     213             : 
     214           0 :       if(present(hint)) write (*,*) "DEN of ", hint
     215           0 :       write (*,*) "q   = ", q
     216           0 :       write (*,*) "qis = ", qis
     217             : 
     218           0 :       write (*,*) "qmt"
     219           0 :       do n_mt = 1,size(qmt, dim=1)
     220           0 :          write (*,*) "mt = ", n_mt, qmt(n_mt,:)
     221             :       enddo
     222             : 
     223           0 :       if(.not. any(ieee_is_nan(qvac))) then
     224           0 :          write (*, *) "qvac",    qvac
     225             :       endif
     226           0 :       write (*, *) "qtot",    qtot
     227           0 :       write (*, *) "qis_tot", qistot
     228           0 :       write (*, *) "-------------------------"
     229           0 :    END SUBROUTINE print_cdn_inte
     230           0 : END MODULE m_cdntot

Generated by: LCOV version 1.13