1! 2! Copyright (C) 2010 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!--------------------------------------------------------------- 9subroutine start_scheq(lam,e,b,grid,ze2,solution) 10 !--------------------------------------------------------------- 11 ! 12 ! determines the wave-function in the first two points by 13 ! series developement. It receives as input: 14 ! lam the angular momentum 15 ! e the energy in Ry 16 ! b(0:3) the coefficients of a polynomial that interpolates the 17 ! potential in the first three points 18 ! grid the mesh 19 ! ze2 the zed of the mesh 20 ! in output solution(1:2) contains the solution in the first two points 21 ! 22USE kinds, ONLY : DP 23USE radial_grids, ONLY : radial_grid_type 24IMPLICIT NONE 25TYPE(radial_grid_type), INTENT(IN) :: grid 26INTEGER, INTENT(IN) :: lam 27REAL(DP), INTENT(IN) :: b(0:3), e 28REAL(DP) :: ze2, xl1, x4l6, x6l12, x8l20, b0e, c1, c2, c3, c4, rr1, rr2 29REAL(DP) :: solution(1:2) 30INTEGER :: l1 31! 32! set up constants and initialize 33! 34l1=lam+1 35xl1=lam+1.0_DP 36x4l6=4.0_dp*lam+6.0_dp 37x6l12=6.0_dp*lam+12.0_dp 38x8l20=8.0_dp*lam+20.0_dp 39! 40! 41b0e=b(0)-e 42c1=0.5_dp*ze2/xl1 43c2=(c1*ze2+b0e)/x4l6 44c3=(c2*ze2+c1*b0e+b(1))/x6l12 45c4=(c3*ze2+c2*b0e+c1*b(1)+b(2))/x8l20 46rr1=(1.0_dp+grid%r(1)*(c1+grid%r(1)*(c2+grid%r(1)*(c3+grid%r(1)*c4))))*grid%r(1)**l1 47rr2=(1.0_dp+grid%r(2)*(c1+grid%r(2)*(c2+grid%r(2)*(c3+grid%r(2)*c4))))*grid%r(2)**l1 48solution(1)=rr1/grid%sqr(1) 49solution(2)=rr2/grid%sqr(2) 50 51return 52end subroutine start_scheq 53