TABLE OF CONTENTS


ABINIT/psp11nl [ Functions ]

[ Top ] [ Functions ]

NAME

 psp11nl

FUNCTION

 Fourier transform the real space UPF projector functions to reciprocal space

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (MJV)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  lmax=maximum ang momentum for which nonlocal form factor is desired.
   Usually lmax=1, sometimes = 0 (e.g. for oxygen); lmax <= 2 allowed.
  mmax=number of radial grid points for atomic grid
  lnmax= maximum index for all l channel projectors, dimension of ffspl
  lmnmax= maximum index for all projectors, dimension of indlmn
  mqgrid=number of grid points for q grid
  n_proj = total number of NL projectors read in
  proj = projector data times r, on a real space grid
  proj_l = ang mom channel for each projector
  proj_np = max number of points used for each projector
  qgrid(mqgrid)=values at which form factors are returned
  r(mmax)=radial grid values
  drdi=derivative of grid point wrt index
  useylm = input to use m dependency of NL part, or only Legendre polynomials

OUTPUT

  ffspl(mqgrid,2,mpsang)=Kleinman-Bylander form factor f_l(q) and
   second derivative from spline fit for each angular momentum
  indlmn = indexing of each projector, for n, l, m, s, ln, lmn (see pspatm.F90)

NOTES

PARENTS

      upf2abinit

CHILDREN

      ctrap,jbessel,spline

SOURCE

 47 #if defined HAVE_CONFIG_H
 48 #include "config.h"
 49 #endif
 50 
 51 #include "abi_common.h"
 52 
 53 subroutine psp11nl(ffspl,indlmn,mmax,lnmax,lmnmax,mqgrid,n_proj,&
 54 &                  proj, proj_l, proj_np, qgrid, r, drdi, useylm)
 55 
 56  use defs_basis
 57  use m_splines
 58  use m_profiling_abi
 59  use m_errors
 60 
 61  use m_numeric_tools, only : ctrap
 62  use m_paw_numeric, only: jbessel=>paw_jbessel
 63 
 64 !This section has been created automatically by the script Abilint (TD).
 65 !Do not modify the following lines by hand.
 66 #undef ABI_FUNC
 67 #define ABI_FUNC 'psp11nl'
 68 !End of the abilint section
 69 
 70  implicit none
 71 
 72 !Arguments ------------------------------------
 73 !scalars
 74  integer,intent(in) :: mmax, lnmax, lmnmax, mqgrid, useylm, n_proj
 75 !arrays
 76  integer, intent(in) :: proj_l(n_proj)
 77  integer, intent(in) :: proj_np(n_proj)
 78  integer, intent(out) :: indlmn(6,lmnmax)
 79  real(dp),intent(in) :: r(mmax)
 80  real(dp),intent(in) :: drdi(mmax)
 81  real(dp),intent(in) :: proj(mmax,n_proj)
 82  real(dp),intent(in) :: qgrid(mqgrid)
 83  real(dp),intent(inout) :: ffspl(mqgrid,2,lnmax) !vz_i
 84 
 85 !Local variables-------------------------------
 86 !scalars
 87  integer :: iproj, np, ll, llold, ipsang, i_indlmn
 88  integer :: iproj_1l, ir, iq, mm
 89  integer :: bessorder
 90  real(dp) :: res, arg, besfact, dummy, dummy2
 91  real(dp), allocatable :: work(:)
 92  character(len=500) :: message
 93 
 94 !*************************************************************************
 95  bessorder = 0 ! never calculate derivatives of bessel functions
 96 
 97  ffspl = zero
 98  indlmn = 0
 99  i_indlmn = 0
100  llold = -1
101  iproj_1l = 1
102 !big loop over all projectors
103  do iproj = 1, n_proj
104 
105    if (iproj > lmnmax) then
106      write(message,'(a,2i0)') ' Too many projectors found. n_proj, lmnmax =  ',n_proj, lmnmax
107      MSG_ERROR(message)
108    end if
109 
110    np = proj_np(iproj)
111    ABI_ALLOCATE(work,(np))
112    ll = proj_l(iproj)
113    if (ll < llold) then
114      message = 'psp11nl : Error: UPF projectors are not in order of increasing ll'
115      MSG_ERROR(message)
116    else if (ll == llold) then
117      iproj_1l = iproj_1l + 1
118    else
119      iproj_1l = 1
120      llold = ll
121    end if
122 !  determine indlmn for this projector (keep in UPF order and enforce that they are in
123 !  increasing ll)
124    do mm = 1, 2*ll*useylm+1
125      i_indlmn = i_indlmn + 1
126      indlmn(1,i_indlmn) = ll
127      indlmn(2,i_indlmn) = mm-ll*useylm-1
128      indlmn(3,i_indlmn) = iproj_1l
129      indlmn(4,i_indlmn) = ll*ll+(1-useylm)*ll+mm
130      indlmn(5,i_indlmn) = iproj
131      indlmn(6,i_indlmn) = 1 !spin? FIXME: to get j for relativistic cases
132    end do
133 
134 !  FT projectors to reciprocal space q
135    do iq = 1, mqgrid
136      arg = two_pi*qgrid(iq)
137 
138 !    FIXME: add semianalytic form for integral from 0 to first point
139      do ir = 1, np
140        call jbessel(besfact, dummy, dummy2, ll, bessorder, arg*r(ir))
141 !      besfact = sin(arg*r(ir))
142        work(ir) = drdi(ir) * besfact * proj(ir, iproj) * r(ir) !* r(ir)
143      end do
144      call ctrap (np, work, one, res)
145 
146      ffspl(iq, 1, iproj) = res
147    end do
148    ABI_DEALLOCATE(work)
149  end do  ! iproj
150 
151 !add derivative of ffspl(:,1,:) for spline interpolation later
152  ABI_ALLOCATE(work,(mqgrid))
153  do ipsang = 1, lnmax
154    call spline(qgrid,ffspl(:,1,ipsang),mqgrid,zero,zero,ffspl(:,2,ipsang))
155  end do
156  ABI_DEALLOCATE(work)
157 
158 end subroutine psp11nl