File: genfdufr.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (67 lines) | stat: -rw-r--r-- 1,626 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67

! Copyright (C) 2008  F. Bultmark, F. Cricchio, L. Nordstrom and J. K. Dewhurst.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

!BOP
! !ROUTINE: genfdufr
! !INTERFACE:
subroutine genfdufr
! !USES:
use modmain
use moddftu
! !DESCRIPTION:
!   Generates the radial functions used to calculate the Slater integrals
!   through a Yukawa potential.
!
! !REVISION HISTORY:
!   Created April 2008 from genapwfr (Francesco Cricchio)
!EOP
!BOC
implicit none
! local variables
integer is,ia,ias
integer nr,nri,ir
integer nn,i,l
real(8) t1
! automatic arrays
real(8) vr(nrmtmax),fr(nrmtmax)
real(8) p0(nrmtmax),p1(nrmtmax),q0(nrmtmax),q1(nrmtmax)
! external functions
real(8) splint
external splint
do i=1,ndftu
  is=idftu(1,i)
  l=idftu(2,i)
  nr=nrmt(is)
  nri=nrmti(is)
  do ia=1,natoms(is)
    ias=idxas(ia,is)
    call rfmtlm(1,nr,nri,vsmt(:,ias),vr)
    vr(1:nr)=vr(1:nr)*y00
! integrate the radial Schrodinger equation
    call rschrodint(solsc,l,fdue(l,ias),nr,rsp(:,is),vr,nn,p0,p1,q0,q1)
! normalise radial functions
    fr(1:nr)=p0(1:nr)**2
    t1=splint(nr,rsp(:,is),fr)
    if (t1.lt.1.d-20) then
      write(*,*)
      write(*,'("Error(genfdufr): degenerate APW radial functions")')
      write(*,'(" for species ",I4)') is
      write(*,'(" atom ",I4)') ia
      write(*,'(" and angular momentum ",I4)') l
      write(*,*)
      stop
    end if
    t1=1.d0/sqrt(abs(t1))
    p0(1:nr)=t1*p0(1:nr)
! divide by r and store in global array
    do ir=1,nr
      fdufr(ir,l,ias)=p0(ir)/rsp(ir,is)
    end do
  end do
end do
return
end subroutine
!EOC