Line data Source code
1 : c*************************************************c
2 : c compute radial integrals of sph. Hamiltonian c
3 : c in scalar relativistic approximation between c
4 : c products of radial function u and spherical c
5 : c Bessel functions c
6 : c c
7 : c integral = < uj1 | H_sra(lzb) | uj2 > c
8 : c*************************************************c
9 : c J.-P. Hanke, Dec. 2015 c
10 : c*************************************************c
11 : MODULE m_wann_uHu_radintsra4
12 : CONTAINS
13 :
14 0 : SUBROUTINE wann_uHu_radintsra4(jmtd,jri,rmsh,dx,e,vr,vm,
15 0 : > uj1,uj2,duj2,lmaxd,lzb,
16 : > integral)
17 :
18 : USE m_intgr, ONLY : intgr3
19 : USE m_constants
20 : IMPLICIT NONE
21 :
22 : REAL difcub
23 : EXTERNAL difcub
24 :
25 : INTEGER, INTENT(IN) :: jri,jmtd,lmaxd
26 : REAL, INTENT(OUT) :: integral
27 : REAL, INTENT(IN) :: vr(jmtd),vm(jmtd),rmsh(jmtd),dx,e
28 : REAL, INTENT(IN) :: uj1(jmtd,2) ! u(b1)*j
29 : REAL, INTENT(IN) :: uj2(jmtd,2) ! u(b2)*j
30 : REAL, INTENT(IN) :: duj2(jmtd,2) ! d/dr (u(b2)j)
31 : INTEGER, INTENT(IN) :: lzb ! l of zentrifugal barrier
32 :
33 0 : REAL, ALLOCATABLE :: x(:)
34 : REAL :: c,c2,cin2,cin
35 : REAL :: ll,xi,vv,mm,vv2
36 : REAL :: sfp
37 : INTEGER :: i,j
38 :
39 0 : c = c_light(1.)
40 0 : c2 = c*c
41 0 : cin = 1./c
42 0 : cin2 = cin*cin
43 0 : ll = lzb*(lzb+1)
44 : sfp = sqrt(4.0*pimach())
45 :
46 0 : allocate( x(jri) )
47 :
48 : ! compute matrix elements of semi-relativistic
49 : ! Hamiltonian [Eq.(3.54) in PhD thesis of P.Kurz]
50 0 : DO i = 1, jri
51 0 : xi = rmsh(i)
52 0 : vv = vr(i) / xi !* sfp
53 0 : vv2= vm(i) / xi !* sfp
54 0 : mm = 1. + 0.5 * cin2 * ( e - vv2 )
55 : x(i) =
56 : > ! large-H-large
57 : > uj1(i,1) * uj2(i,1) * ( 0.5 / mm * ll / xi / xi + vv )
58 : > ! small-H-small
59 : > + uj1(i,2) * uj2(i,2) * ( -2. * c2 + vv )
60 : > ! large-H-small (not symmetrized)
61 : > - c * uj1(i,1) * (2. * uj2(i,2) / xi + duj2(i,2) )
62 : > ! small-H-large (not symmetrized)
63 0 : > + c * uj1(i,2) * duj2(i,1)
64 : ENDDO
65 0 : call intgr3(x,rmsh,dx,jri,integral)
66 :
67 0 : deallocate( x )
68 :
69 0 : END SUBROUTINE wann_uHu_radintsra4
70 : END MODULE m_wann_uHu_radintsra4
|