LCOV - code coverage report
Current view: top level - src/65_paw - m_paw_pwaves_lmn.F90 (source / functions) Coverage Total Hit
Test: coverage.info Lines: 82.1 % 162 133
Test Date: 2026-11-09 22:30:52 Functions: 50.0 % 4 2
Branches: 41.6 % 522 217

             Branch data     Line data    Source code
       1                 :             : !!****m* ABINIT/m_paw_pwaves_lmn
       2                 :             : !! NAME
       3                 :             : !!  m_paw_pwaves_lmn
       4                 :             : !!
       5                 :             : !! FUNCTION
       6                 :             : !!  This module contains the definition of the paw_pwaves_lmn structured datatype,
       7                 :             : !!  as well as related functions and methods.
       8                 :             : !!  paw_pwaves_lmn is used to store the 3D values of the all-electron and of
       9                 :             : !!  the pseudized part of the PAW partial waves on the set of FFT points falling
      10                 :             : !!  inside the spheres around each atom.
      11                 :             : !!
      12                 :             : !! COPYRIGHT
      13                 :             : !! Copyright (C) 2008-2026 ABINIT group (MG,MT)
      14                 :             : !! This file is distributed under the terms of the
      15                 :             : !! GNU General Public License, see ~abinit/COPYING
      16                 :             : !! or http://www.gnu.org/copyleft/gpl.txt .
      17                 :             : !!
      18                 :             : !! NOTES
      19                 :             : !! * Routines tagged with "@type_name" are strongly connected to the definition of the data type.
      20                 :             : !!   Strongly connected means that the proper functioning of the implementation relies on the
      21                 :             : !!   assumption that the tagged procedure is consistent with the type declaration.
      22                 :             : !!   Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure
      23                 :             : !!   that all the strongly connected routines are changed accordingly to accommodate the modification of the data type.
      24                 :             : !!   Typical examples of strongly connected routines are creation, destruction or reset methods.
      25                 :             : !!
      26                 :             : !! SOURCE
      27                 :             : 
      28                 :             : #if defined HAVE_CONFIG_H
      29                 :             : #include "config.h"
      30                 :             : #endif
      31                 :             : 
      32                 :             : #include "abi_common.h"
      33                 :             : 
      34                 :             : MODULE m_paw_pwaves_lmn
      35                 :             : 
      36                 :             :  use defs_basis
      37                 :             :  use m_errors
      38                 :             :  use m_abicore
      39                 :             :  use m_xmpi
      40                 :             :  use m_sort
      41                 :             : 
      42                 :             :  use m_numeric_tools, only : wrap2_zero_one
      43                 :             :  use m_geometry,      only : xcart2xred
      44                 :             :  use m_pawrad,        only : pawrad_type, nderiv_gen, pawrad_deducer0
      45                 :             :  use m_pawtab,        only : pawtab_type
      46                 :             :  use m_pawfgrtab,     only : pawfgrtab_type
      47                 :             :  use m_paw_numeric,   only : paw_spline, paw_splint
      48                 :             :  use m_paw_sphharm,   only : initylmr
      49                 :             :  use m_paral_atom,    only : get_my_atmtab, free_my_atmtab
      50                 :             : 
      51                 :             :  implicit none
      52                 :             : 
      53                 :             :  private
      54                 :             : !!***
      55                 :             : 
      56                 :             : !----------------------------------------------------------------------
      57                 :             : 
      58                 :             : !!****t* m_paw_pwaves_lmn/paw_pwaves_lmn_t
      59                 :             : !! NAME
      60                 :             : !! paw_pwaves_lmn_t
      61                 :             : !!
      62                 :             : !! FUNCTION
      63                 :             : !!  Datatype used to store the 3D values of the all-electron and of the pseudized part of the
      64                 :             : !!  PAW partial waves on the set of FFT points falling inside the spheres around each atom.
      65                 :             : !!  The data is mainly used for plotting the true PAW wavefunctions in real space.
      66                 :             : !!
      67                 :             : !! SOURCE
      68                 :             : 
      69                 :             :  type,public ::  paw_pwaves_lmn_t
      70                 :             : 
      71                 :             :   integer :: nfgd
      72                 :             : 
      73                 :             :   integer :: lmn_size
      74                 :             : 
      75                 :             :   !$integer :: ngfft(18)
      76                 :             : 
      77                 :             :   integer,allocatable :: r0shift(:,:)
      78                 :             :   ! r0shift(3,nfgd)
      79                 :             : 
      80                 :             :   !real(dp),allocatable:: phk_atm(:,:)
      81                 :             :   ! phk_atmt(2,nfgd)
      82                 :             : 
      83                 :             :   real(dp),allocatable :: phi(:,:)
      84                 :             :   ! phi (nfgd,lmn_size)
      85                 :             :   ! \phi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type).
      86                 :             : 
      87                 :             :   real(dp),allocatable :: tphi(:,:)
      88                 :             :   ! tphi (nfgd,lmn_size)
      89                 :             :   ! \tphi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type).
      90                 :             : 
      91                 :             :   real(dp),allocatable :: phi_gr(:,:,:)
      92                 :             :   ! phi_gr (3,nfgd,lmn_size)
      93                 :             :   ! gradient, in cartesian coordinates, of \phi_{nlm}(ifgd) for each point of the FFT mesh
      94                 :             :   ! located inside the PAW sphere (see pawfgrtab_type).
      95                 :             : 
      96                 :             :   real(dp),allocatable :: tphi_gr(:,:,:)
      97                 :             :   ! tphi_gr (3,nfgd,lmn_size)
      98                 :             :   ! gradient, in cartesian coordinates, of \tphi_{nlm}(ifgd) for each point of the FFT mesh
      99                 :             :   ! located inside the PAW sphere (see pawfgrtab_type).
     100                 :             : 
     101                 :             :  end type paw_pwaves_lmn_t
     102                 :             : 
     103                 :             :  public :: paw_pwaves_lmn_init
     104                 :             :  public :: paw_pwaves_lmn_free
     105                 :             : !!***
     106                 :             : 
     107                 :             : CONTAINS !===========================================================
     108                 :             : !!***
     109                 :             : 
     110                 :             : !----------------------------------------------------------------------
     111                 :             : 
     112                 :             : !!****f* m_paw_pwaves_lmn/paw_pwaves_lmn_init
     113                 :             : !! NAME
     114                 :             : !! paw_pwaves_lmn_init
     115                 :             : !!
     116                 :             : !! FUNCTION
     117                 :             : !!
     118                 :             : !! INPUTS
     119                 :             : !!  mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc
     120                 :             : !!  comm_atom=--optional-- MPI communicator over atoms
     121                 :             : !!
     122                 :             : !! OUTPUT
     123                 :             : !!
     124                 :             : !! SOURCE
     125                 :             : 
     126   [ +  +  -  +  :           3 : subroutine paw_pwaves_lmn_init(Paw_onsite,my_natom,natom,ntypat,rprimd,xcart,Pawtab, &
          -  +  -  +  -  
          +  -  +  +  +  
                   +  - ]
     127                 :           2 : &                              Pawrad,local_pawfgrtab,optgrad,&
     128   [ -  +  -  -  :           1 : &                              mpi_atmtab,comm_atom) ! optional arguments (parallelism)
                   -  - ]
     129                 :             : 
     130                 :             : !Arguments ------------------------------------
     131                 :             : !scalars
     132                 :             :  integer,intent(in) :: my_natom,natom,ntypat
     133                 :             :  integer,optional,intent(in) :: optgrad,comm_atom
     134                 :             : !arrays
     135                 :             :  integer,optional,target,intent(in) :: mpi_atmtab(:)
     136                 :             :  real(dp),intent(in) :: xcart(3,natom),rprimd(3,3)
     137                 :             :  type(pawtab_type),target,intent(in) :: Pawtab(ntypat)
     138                 :             :  type(pawrad_type),intent(in) :: Pawrad(ntypat)
     139                 :             :  type(pawfgrtab_type),intent(in) :: local_pawfgrtab(my_natom)
     140                 :             :  type(paw_pwaves_lmn_t),intent(out) :: Paw_onsite(my_natom)
     141                 :             : 
     142                 :             : !Local variables-------------------------------
     143                 :             : !scalars
     144                 :             :  integer :: itypat,ln_size,lmn_size,mesh_size,inl,iatom,iatom1,my_comm_atom,my_optgrad
     145                 :             :  integer :: nfgd,ifgd,ipsang,option_ylmr,normchoice,ii,jlmn,jl,jm,jlm,jln
     146                 :             :  logical :: my_atmtab_allocated,paral_atom
     147                 :             :  real(dp) :: phj,rR,tphj,ybcbeg,ybcend
     148                 :             : !arrays
     149                 :           1 :  integer, allocatable :: iperm(:)
     150                 :           1 :  integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
     151                 :           1 :  integer,pointer :: my_atmtab(:)
     152                 :             :  real(dp) :: yvals(4),red(3),shift(3)
     153                 :           1 :  real(dp),allocatable :: ff(:),nrm(:),nrm_sort(:),phigrd(:,:),tphigrd(:,:),ylm_tmp(:,:),ylm(:,:),ylm_gr(:,:,:)
     154                 :           1 :  real(dp),allocatable :: rsph_red(:,:),rsph_cart(:,:),phigrd_gr(:,:),tphigrd_gr(:,:),gg(:)
     155                 :           1 :  type(paw_pwaves_lmn_t),allocatable :: Paw_lmn_spline(:)
     156                 :             : 
     157                 :             : ! *************************************************************************
     158                 :             : 
     159                 :             : !@paw_pwaves_lmn_t
     160                 :             : 
     161         [ +  - ]:           1 :  if (my_natom==0) return
     162         [ -  + ]:           1 :  my_optgrad = -1; if (present(optgrad)) my_optgrad = optgrad
     163                 :             : 
     164   [ +  +  +  -  :           2 :  ABI_CHECK(all(local_pawfgrtab(:)%rfgd_allocated==1),"R vectors not allocated in pawfgrtab")
                   -  + ]
     165                 :             : 
     166                 :             :  ! Set up parallelism over atoms
     167   [ -  +  -  - ]:           1 :  paral_atom=(present(comm_atom).and.(my_natom/=natom))
     168         [ -  + ]:           1 :  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
     169         [ -  + ]:           1 :  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
     170                 :           1 :  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
     171                 :             : 
     172                 :             :  ! Prepare the spline. Calculate 2nd derivatives of partial waves for each atom type.
     173   [ -  +  +  -  :           4 :  ABI_MALLOC(Paw_lmn_spline,(ntypat))
          +  -  +  -  -  
             +  -  +  +  
                      + ]
     174                 :             : 
     175         [ +  + ]:           2 :  do itypat=1,ntypat
     176                 :           1 :    ln_size  =Pawtab(itypat)%basis_size
     177                 :           1 :    mesh_size=Pawtab(itypat)%mesh_size
     178                 :             : 
     179   [ -  +  +  -  :           4 :    ABI_MALLOC(Paw_lmn_spline(itypat)%phi ,(mesh_size,ln_size))
          -  +  +  -  +  
          -  +  -  -  +  
             -  +  -  + ]
     180   [ -  +  +  -  :           3 :    ABI_MALLOC(Paw_lmn_spline(itypat)%tphi,(mesh_size,ln_size))
          -  +  +  -  -  
             +  -  +  -  
                      + ]
     181                 :             : 
     182         [ +  + ]:           7 :    do inl=1,ln_size ! Calculate 2nd derivatives of %phi and %tphi for each ln component.
     183                 :           5 :      ybcbeg=zero; ybcend=zero
     184                 :             :      call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl), mesh_size,&
     185                 :           5 : &                    ybcbeg,ybcend,Paw_lmn_spline(itypat)%phi(:,inl))
     186                 :             : 
     187                 :             :      ybcbeg=zero; ybcend=zero
     188                 :             :      call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),mesh_size,&
     189                 :           6 : &                    ybcbeg,ybcend,Paw_lmn_spline(itypat)%tphi(:,inl))
     190                 :             :    end do
     191                 :             :  end do
     192                 :             :  !
     193                 :             :  ! === spline for each atom ===
     194                 :             :  ! * FFT points within PAW sphere depend on the atom site.
     195         [ +  + ]:           2 :  do iatom1=1,my_natom
     196         [ -  + ]:           1 :    iatom=iatom1;if (paral_atom) iatom=my_atmtab(iatom1)
     197                 :             : 
     198                 :           1 :    itypat   = local_pawfgrtab(iatom1)%itypat
     199                 :           1 :    ln_size  = Pawtab(itypat)%basis_size
     200                 :           1 :    lmn_size = Pawtab(itypat)%lmn_size
     201                 :           1 :    mesh_size= Pawrad(itypat)%mesh_size
     202                 :           1 :    nfgd     = local_pawfgrtab(iatom1)%nfgd ! no. of points in the fine grid for this PAW sphere
     203                 :           1 :    indlmn   => Pawtab(itypat)%indlmn
     204                 :             : 
     205                 :             :    ! The points in the PAW sphere might belong to a different unit cell, in this case one has to
     206                 :             :    ! reconstruct the contribution to the AE psi_k given by the lattice-symmetric atom in another sphere.
     207                 :             :    ! Here I wrap rr back into the first unit cell keeping trace of the lattice vector L
     208                 :             :    ! needed so that rr = rr_first_ucell + L
     209                 :             :    ! The contribution to the AE u(r) in the first unit cell has to be multiplied by e^{-ikL}.
     210                 :             : 
     211                 :           1 :    Paw_onsite(iatom1)%nfgd     = nfgd
     212                 :           1 :    Paw_onsite(iatom1)%lmn_size = lmn_size
     213                 :             :    !
     214   [ -  +  +  -  :           3 :    ABI_MALLOC(Paw_onsite(iatom1)%r0shift,(3,nfgd))
          +  -  +  -  -  
             +  -  +  -  
                      + ]
     215                 :             : 
     216   [ -  +  +  -  :           3 :    ABI_MALLOC(rsph_red,(3,nfgd))
          +  -  +  -  -  
                +  -  + ]
     217   [ -  +  +  -  :           2 :    ABI_MALLOC(rsph_cart,(3,nfgd))
             -  +  -  + ]
     218         [ +  + ]:        1430 :    do ifgd=1,nfgd
     219         [ +  + ]:        5717 :      rsph_cart(:,ifgd) = local_pawfgrtab(iatom1)%rfgd(:,ifgd) + xcart(:,iatom)
     220                 :             :    end do
     221                 :           1 :    call xcart2xred(nfgd,rprimd,rsph_cart,rsph_red) ! go to reduced coordinates.
     222         [ -  + ]:           1 :    ABI_FREE(rsph_cart)
     223                 :             : 
     224         [ +  + ]:        1430 :    do ifgd=1,nfgd
     225                 :        1429 :      call wrap2_zero_one(rsph_red(1,ifgd),red(1),shift(1)) ! rr = r_cell + shift
     226                 :        1429 :      call wrap2_zero_one(rsph_red(2,ifgd),red(2),shift(2))
     227                 :        1429 :      call wrap2_zero_one(rsph_red(3,ifgd),red(3),shift(3))
     228         [ +  + ]:        5717 :      Paw_onsite(iatom1)%r0shift(:,ifgd) = NINT(shift)
     229                 :             :      !if (ANY( ABS(shift) > tol12)) then
     230                 :             :        !ABI_WARNING("rmR_red is outside the first unit cell.")
     231                 :             :        !write(ab_out,*)rsph_red(:,ifgd),shift
     232                 :             :      !end if
     233                 :             :    end do
     234         [ -  + ]:           1 :    ABI_FREE(rsph_red)
     235                 :             :    !
     236                 :             :    ! * Obtain |r-R| on fine grid, note that rfgd is given in Cartesian coordinates.
     237   [ -  +  +  -  :           3 :    ABI_MALLOC(nrm,(nfgd))
          +  -  +  -  -  
                +  -  + ]
     238         [ +  + ]:        1430 :    do ifgd=1,nfgd
     239         [ +  + ]:        5717 :      nrm(ifgd) = sqrt(dot_product(local_pawfgrtab(iatom1)%rfgd(:,ifgd),local_pawfgrtab(iatom1)%rfgd(:,ifgd)))
     240                 :             :    end do
     241                 :             :    !
     242                 :             :    ! * Compute Ylm for each r-R vector.
     243                 :           1 :    ipsang = 1 + (Pawtab(itypat)%l_size-1)/2 ! recall l_size=2*l_max-1 where l_max is shifted by 1.
     244   [ -  +  +  -  :           4 :    ABI_MALLOC(ylm_tmp,(ipsang**2,nfgd))
          -  +  +  -  +  
          -  +  -  -  +  
                   -  + ]
     245                 :           1 :    normchoice = 1 ! Use computed norms of input vectors.
     246         [ -  + ]:           1 :    if (my_optgrad==1) then
     247                 :           0 :      option_ylmr=2 ! Compute Ylm(r-R) and its gradient
     248   [ #  #  #  #  :           0 :      ABI_MALLOC(ylm_gr,(3,ipsang**2,nfgd))
          #  #  #  #  #  
          #  #  #  #  #  
                   #  # ]
     249                 :             :    else
     250                 :           1 :      option_ylmr= 1 ! To compute Ylm(r-R).
     251         [ -  + ]:           1 :      ABI_MALLOC(ylm_gr,(3,3,0))
     252                 :             :    end if
     253                 :           1 :    call initylmr(ipsang,normchoice,nfgd,nrm,option_ylmr,local_pawfgrtab(iatom1)%rfgd,ylm_tmp,ylm_gr)
     254                 :             :    !
     255                 :             :    !  Exchange dimensions for better memory access.
     256   [ -  +  +  -  :           3 :    ABI_MALLOC(ylm,(nfgd,ipsang**2))
          -  +  +  -  -  
                +  -  + ]
     257         [ +  + ]:          10 :    do ii=1,ipsang**2
     258         [ +  + ]:       12871 :      ylm(:,ii) = ylm_tmp(ii,:)
     259                 :             :    end do
     260         [ -  + ]:           1 :    ABI_FREE(ylm_tmp)
     261                 :             :    !
     262                 :             :    ! In order to do spline fits, the |r-R| data must be sorted
     263                 :             :    ! Here we sort the nrm points, keeping track of which goes where
     264   [ -  +  +  -  :           2 :    ABI_MALLOC(nrm_sort,(nfgd))
             -  +  -  + ]
     265   [ -  +  +  + ]:        1431 :    nrm_sort = nrm
     266                 :             : 
     267   [ -  +  +  -  :           3 :    ABI_MALLOC(iperm,(nfgd))
          +  -  +  -  -  
                +  -  + ]
     268         [ +  + ]:        1430 :    do ifgd=1,nfgd
     269                 :        1430 :      iperm(ifgd)=ifgd
     270                 :             :    end do
     271                 :             : 
     272                 :           1 :    call sort_dp(nfgd,nrm_sort,iperm,tol8)
     273                 :             : 
     274                 :             :    ! Now make spline fits of phi and tphi  onto the fine grid around the atom
     275   [ -  +  +  -  :           4 :    ABI_MALLOC(phigrd,(nfgd,ln_size))
          -  +  +  -  +  
          -  +  -  -  +  
                   -  + ]
     276   [ -  +  +  -  :           3 :    ABI_MALLOC(tphigrd,(nfgd,ln_size))
          -  +  +  -  -  
                +  -  + ]
     277   [ -  +  +  -  :           2 :    ABI_MALLOC(ff,(nfgd))
             -  +  -  + ]
     278         [ -  + ]:           1 :    if (my_optgrad==1) then
     279   [ #  #  #  #  :           0 :      ABI_MALLOC(phigrd_gr,(nfgd,ln_size))
          #  #  #  #  #  
             #  #  #  #  
                      # ]
     280   [ #  #  #  #  :           0 :      ABI_MALLOC(tphigrd_gr,(nfgd,ln_size))
          #  #  #  #  #  
             #  #  #  #  
                      # ]
     281   [ #  #  #  #  :           0 :      ABI_MALLOC(gg,(mesh_size))
          #  #  #  #  #  
             #  #  #  #  
                      # ]
     282                 :             :    end if
     283                 :             : 
     284         [ +  + ]:           6 :    do inl=1,ln_size
     285                 :             :      !
     286                 :             :      ! * splint phi onto points and reorder indices.
     287                 :             :      call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl),&
     288                 :           5 : &                    Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
     289         [ +  + ]:        7150 :      do ifgd=1,nfgd
     290                 :        7145 :        ii=iperm(ifgd)
     291                 :        7150 :        phigrd(ii,inl) = ff(ifgd)
     292                 :             :      end do
     293                 :             :      !
     294                 :             :      ! * compute d phi/dr, interpolate onto points and reorder indices.
     295         [ -  + ]:           5 :      if (my_optgrad==1) then
     296                 :           0 :        ybcbeg=zero; ybcend=zero
     297                 :           0 :        call nderiv_gen(gg,Pawtab(itypat)%phi(:,inl),Pawrad(itypat))
     298                 :             :        call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,&
     299                 :           0 : &                      Paw_lmn_spline(itypat)%phi(:,inl))
     300                 :             :        call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%phi(:,inl),&
     301                 :           0 : &                      Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
     302         [ #  # ]:           0 :        do ifgd=1,nfgd
     303                 :           0 :          ii=iperm(ifgd)
     304                 :           0 :          phigrd_gr(ii,inl)  = ff(ifgd)
     305                 :             :        end do
     306                 :             :      end if
     307                 :             :      !
     308                 :             :      ! * compute d tphi/dr, interpolate onto points and reorder indices.
     309                 :             :      call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),&
     310                 :           5 : &                    Paw_lmn_spline(itypat)%tphi(:,inl),nfgd,nrm_sort,ff)
     311         [ +  + ]:        7150 :      do ifgd=1,nfgd
     312                 :        7145 :        ii=iperm(ifgd)
     313                 :        7150 :        tphigrd(ii,inl) = ff(ifgd)
     314                 :             :      end do
     315         [ -  + ]:           6 :      if (my_optgrad==1) then
     316                 :           0 :        ybcbeg=zero; ybcend=zero
     317                 :           0 :        call nderiv_gen(gg,Pawtab(itypat)%tphi(:,inl),Pawrad(itypat))
     318                 :             :        call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,&
     319                 :           0 : &                      Paw_lmn_spline(itypat)%tphi(:,inl))
     320                 :             :        call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%tphi(:,inl),&
     321                 :           0 : &                      Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff)
     322         [ #  # ]:           0 :        do ifgd=1,nfgd
     323                 :           0 :          ii=iperm(ifgd)
     324                 :           0 :          tphigrd_gr(ii,inl)  = ff(ifgd)
     325                 :             :        end do
     326                 :             :      end if
     327                 :             :    end do !inl
     328                 :             : 
     329         [ -  + ]:           1 :    ABI_FREE(ff)
     330         [ -  + ]:           1 :    if (my_optgrad==1) then
     331         [ #  # ]:           0 :      ABI_FREE(gg)
     332                 :             :    end if
     333                 :             :    !
     334                 :             :    ! === Calculate AE and PS partial waves inside the sphere ===
     335                 :             :    ! * recall that <r|phi>=u(r)*Slm(r^)/r, hence avoid division by zero except for s-waves.
     336   [ -  +  +  -  :           4 :    ABI_MALLOC(Paw_onsite(iatom1)%phi ,(nfgd,lmn_size))
          -  +  +  -  +  
          -  +  -  -  +  
             -  +  -  + ]
     337   [ -  +  +  -  :           3 :    ABI_MALLOC(Paw_onsite(iatom1)%tphi,(nfgd,lmn_size))
          -  +  +  -  -  
             +  -  +  -  
                      + ]
     338                 :             : 
     339         [ -  + ]:           1 :    if (my_optgrad==1) then
     340   [ #  #  #  #  :           0 :      ABI_MALLOC(Paw_onsite(iatom1)%phi_gr ,(3,nfgd,lmn_size))
          #  #  #  #  #  
          #  #  #  #  #  
             #  #  #  # ]
     341   [ #  #  #  #  :           0 :      ABI_MALLOC(Paw_onsite(iatom1)%tphi_gr,(3,nfgd,lmn_size))
          #  #  #  #  #  
             #  #  #  #  
                      # ]
     342                 :             :    end if
     343                 :             : 
     344         [ +  + ]:          14 :    do jlmn=1,lmn_size
     345                 :          13 :      jl  = indlmn(1,jlmn)
     346                 :          13 :      jm  = indlmn(2,jlmn)
     347                 :          13 :      jlm = indlmn(4,jlmn)
     348                 :          13 :      jln = indlmn(5,jlmn)
     349                 :             : 
     350         [ +  + ]:       18591 :      do ifgd=1,nfgd ! loop over fine grid points in current PAW sphere
     351                 :             :        !if (nrm(ifgd)>tol16) then
     352         [ +  + ]:       18590 :        if (nrm(ifgd)>tol10) then ! tol10 to be consistent with initylmr.
     353                 :       18564 :          rR  = nrm(ifgd) ! value of |r-R|
     354                 :             :          !write(ab_out,*) 'rR:',rR,' phigrd:',phigrd(ifgd,jln),' tphigrd:',tphigrd(ifgd,jln),' ylm:',ylm(ifgd,jlm)
     355                 :       18564 :          phj = phigrd (ifgd,jln)*ylm(ifgd,jlm)/rR
     356                 :       18564 :          tphj= tphigrd(ifgd,jln)*ylm(ifgd,jlm)/rR
     357                 :       18564 :          Paw_onsite(iatom1)%phi (ifgd,jlmn) =  phj
     358                 :       18564 :          Paw_onsite(iatom1)%tphi(ifgd,jlmn) = tphj
     359                 :             : 
     360         [ -  + ]:       18564 :          if (my_optgrad==1) then
     361                 :             :            Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = phigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) &
     362         [ #  # ]:           0 : &            + phigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm)
     363                 :             :            Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = tphigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) &
     364         [ #  # ]:           0 : &            + tphigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm)
     365                 :             :          end if
     366                 :             : 
     367                 :             :        else ! Extrapolate if the point is at the origin
     368                 :          13 :          yvals(1) = zero
     369         [ +  + ]:          13 :          if (jl==0) then
     370         [ +  + ]:           8 :            yvals(2:4) = Pawtab(itypat)%phi(2:4,jln)/Pawrad(itypat)%rad(2:4)
     371                 :           2 :            call pawrad_deducer0(yvals,4,Pawrad(itypat))
     372                 :             :          end if
     373                 :          13 :          Paw_onsite(iatom1)%phi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm)
     374                 :          13 :          yvals(1) = zero
     375         [ +  + ]:          13 :          if (jl==0) then
     376         [ +  + ]:           8 :            yvals(2:4) = Pawtab(itypat)%tphi(2:4,jln)/Pawrad(itypat)%rad(2:4)
     377                 :           2 :            call pawrad_deducer0(yvals,4,pawrad(itypat))
     378                 :             :          end if
     379                 :          13 :          Paw_onsite(iatom1)%tphi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm)
     380                 :             :          ! The gradient is expected to go to zero at the origin
     381         [ -  + ]:          13 :          if (my_optgrad==1) then
     382         [ #  # ]:           0 :            Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = zero
     383         [ #  # ]:           0 :            Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = zero
     384                 :             :          end if
     385                 :             :        end if
     386                 :             : 
     387                 :             :      end do !nfgd
     388                 :             :    end do !jlmn
     389                 :             : 
     390         [ -  + ]:           1 :    ABI_FREE(nrm)
     391         [ -  + ]:           1 :    ABI_FREE(nrm_sort)
     392         [ -  + ]:           1 :    ABI_FREE(iperm)
     393         [ -  + ]:           1 :    ABI_FREE(phigrd)
     394         [ -  + ]:           1 :    ABI_FREE(tphigrd)
     395         [ -  + ]:           1 :    ABI_FREE(ylm)
     396         [ -  + ]:           1 :    ABI_FREE(ylm_gr)
     397         [ -  + ]:           2 :    if (my_optgrad==1) then
     398         [ #  # ]:           0 :      ABI_FREE(phigrd_gr)
     399         [ #  # ]:           0 :      ABI_FREE(tphigrd_gr)
     400                 :             :    end if
     401                 :             :  end do !iatom
     402                 :             :  !
     403                 :             :  !* Free 2nd derivates used for spline.
     404                 :           1 :  call paw_pwaves_lmn_free(Paw_lmn_spline)
     405   [ +  -  +  +  :           2 :  ABI_FREE(Paw_lmn_spline)
          -  +  -  +  -  
          +  -  +  -  +  
                   -  + ]
     406                 :             : 
     407                 :             :  ! Destroy atom table used for parallelism
     408                 :           1 :  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
     409                 :             : 
     410   [ -  +  -  +  :           2 : end subroutine paw_pwaves_lmn_init
          -  +  -  -  -  
          -  -  -  -  -  
          -  -  -  -  -  
                      + ]
     411                 :             : !!***
     412                 :             : 
     413                 :             : !----------------------------------------------------------------------
     414                 :             : 
     415                 :             : !!****f* m_paw_pwaves_lmn/paw_pwaves_lmn_free
     416                 :             : !! NAME
     417                 :             : !! paw_pwaves_lmn_free
     418                 :             : !!
     419                 :             : !! FUNCTION
     420                 :             : !!
     421                 :             : !! INPUTS
     422                 :             : !!
     423                 :             : !! OUTPUT
     424                 :             : !!
     425                 :             : !! SOURCE
     426                 :             : 
     427         [ -  + ]:           2 : subroutine paw_pwaves_lmn_free(Paw_onsite)
     428                 :             : 
     429                 :             : !Arguments ------------------------------------
     430                 :             : !arrays
     431                 :             :  type(paw_pwaves_lmn_t),intent(inout) :: Paw_onsite(:)
     432                 :             : 
     433                 :             : !Local variables-------------------------------
     434                 :             : !scalars
     435                 :             :  integer :: iatom
     436                 :             : ! *************************************************************************
     437                 :             : 
     438                 :             : !@paw_pwaves_lmn_t
     439                 :             : 
     440   [ +  -  +  + ]:           6 :  do iatom=LBOUND(Paw_onsite,DIM=1),UBOUND(Paw_onsite,DIM=1)
     441         [ +  - ]:           2 :    ABI_SFREE(Paw_onsite(iatom)%phi)
     442         [ +  - ]:           2 :    ABI_SFREE(Paw_onsite(iatom)%tphi)
     443         [ +  + ]:           2 :    ABI_SFREE(Paw_onsite(iatom)%r0shift)
     444         [ -  + ]:           2 :    ABI_SFREE(Paw_onsite(iatom)%phi_gr)
     445         [ -  + ]:           4 :    ABI_SFREE(Paw_onsite(iatom)%tphi_gr)
     446                 :             :  end do
     447                 :             : 
     448                 :           2 : end subroutine paw_pwaves_lmn_free
     449                 :             : !!***
     450                 :             : 
     451                 :             : !----------------------------------------------------------------------
     452                 :             : 
     453   [ #  #  #  #  :           0 : END MODULE m_paw_pwaves_lmn
          #  #  #  #  #  
          #  #  #  #  #  
          #  #  #  #  #  
                #  #  # ]
     454                 :             : !!***
        

Generated by: LCOV version 2.3-1