mirror of https://gitlab.com/QEF/q-e.git
Cleanup of spherical bessel functions.
Added a routine to compute the derivativs of spherical bessel functions (used by calc_btq). C.S. git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@1450 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
parent
5f208f2c3b
commit
36c6f48b0f
|
@ -22,6 +22,7 @@ scnds.o \
|
|||
simpsn.o \
|
||||
sort.o \
|
||||
sph_bes.o \
|
||||
sph_dbes.o \
|
||||
transto.o \
|
||||
date_and_tim.o \
|
||||
sort_gvec.o \
|
||||
|
|
402
flib/sph_bes.f90
402
flib/sph_bes.f90
|
@ -1,247 +1,283 @@
|
|||
!
|
||||
! Copyright (C) 2001 PWSCF group
|
||||
! Copyright (C) 2001-2004 PWSCF-FPMD-CP90 group
|
||||
! This file is distributed under the terms of the
|
||||
! GNU General Public License. See the file `License'
|
||||
! in the root directory of the present distribution,
|
||||
! or http://www.gnu.org/copyleft/gpl.txt .
|
||||
!
|
||||
!
|
||||
!--------------------------------------------------------------------
|
||||
subroutine sph_bes (msh, r, q, l, jl)
|
||||
!--------------------------------------------------------------------
|
||||
!----------------------------------------------------------------------------
|
||||
SUBROUTINE sph_bes( msh, r, q, l, jl )
|
||||
!----------------------------------------------------------------------------
|
||||
!
|
||||
! input:
|
||||
! msh = number of grid points points
|
||||
! r(1:msh)= radial grid
|
||||
! q = q
|
||||
! l = angular momentum (-1 <= l <= 6)
|
||||
! output:
|
||||
! jl(1:msh) = j_l(q*r(i)) (j_l = spherical bessel function)
|
||||
! ... input:
|
||||
! ... msh = number of grid points points
|
||||
! ... r(1:msh)= radial grid
|
||||
! ... q = q
|
||||
! ... l = angular momentum (-1 <= l <= 6)
|
||||
! ... output:
|
||||
! ... jl(1:msh) = j_l(q*r(i)) (j_l = spherical bessel function)
|
||||
!
|
||||
use kinds, only: DP
|
||||
implicit none
|
||||
USE kinds, ONLY : DP
|
||||
USE constants, ONLY : eps8
|
||||
!
|
||||
integer :: msh, l
|
||||
real(kind=DP) :: r (msh), q, jl (msh)
|
||||
IMPLICIT NONE
|
||||
!
|
||||
real(kind=DP), parameter :: eps = 1.d-8
|
||||
integer :: ir, ir0
|
||||
INTEGER :: msh, l
|
||||
REAL(kind=DP) :: r(msh), q, jl(msh)
|
||||
!
|
||||
#if defined __MASS
|
||||
real(kind=DP) :: qr(msh), sin_qr(msh), cos_qr(msh)
|
||||
INTEGER :: ir, ir0
|
||||
!
|
||||
#if defined (__MASS)
|
||||
REAL(kind=DP) :: qr(msh), sin_qr(msh), cos_qr(msh)
|
||||
#endif
|
||||
!
|
||||
if (abs (q) < eps) then
|
||||
if (l == -1) then
|
||||
call errore ('sph_bes', 'j_{-1}(0) ?!?', 1)
|
||||
elseif (l == 0) then
|
||||
jl(:) = 1.d0
|
||||
else
|
||||
jl(:) = 0.d0
|
||||
endif
|
||||
else
|
||||
if (abs (q * r (1) ) > eps) then
|
||||
IF ( ABS( q ) < eps8 ) THEN
|
||||
!
|
||||
IF ( l == -1 ) THEN
|
||||
!
|
||||
CALL errore( 'sph_bes', 'j_{-1}(0) ?!?', 1 )
|
||||
!
|
||||
ELSE IF ( l == 0 ) THEN
|
||||
!
|
||||
jl(:) = 1.D0
|
||||
!
|
||||
ELSE
|
||||
!
|
||||
jl(:) = 0.D0
|
||||
!
|
||||
END IF
|
||||
!
|
||||
ELSE
|
||||
!
|
||||
IF ( ABS( q*r(1) ) > eps8 ) THEN
|
||||
!
|
||||
ir0 = 1
|
||||
else
|
||||
if (l == -1) then
|
||||
call errore ('sph_bes', 'j_{-1}(0) ?!?', 2)
|
||||
elseif (l == 0) then
|
||||
jl (1) = 1.d0
|
||||
else
|
||||
jl (1) = 0.d0
|
||||
endif
|
||||
!
|
||||
ELSE
|
||||
!
|
||||
IF (l == -1) THEN
|
||||
!
|
||||
CALL errore( 'sph_bes', 'j_{-1}(0) ?!?', 2 )
|
||||
!
|
||||
ELSE IF ( l == 0 ) THEN
|
||||
!
|
||||
jl(1) = 1.D0
|
||||
!
|
||||
ELSE
|
||||
!
|
||||
jl(1) = 0.D0
|
||||
!
|
||||
ENDIF
|
||||
!
|
||||
ir0 = 2
|
||||
endif
|
||||
|
||||
if (l == - 1) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
END IF
|
||||
!
|
||||
IF ( l == - 1 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = cos_qr(ir) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = cos (q * r (ir) ) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = COS (q * r (ir) ) / (q * r (ir) )
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 0) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 0 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vsin( sin_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = sin_qr(ir) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = sin (q * r (ir) ) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = SIN (q * r (ir) ) / (q * r (ir) )
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 1) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 1 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh)
|
||||
CALL vsin( sin_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = ( sin_qr(ir) / (q * r (ir) ) - cos_qr(ir) ) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = (sin (q * r (ir) ) / (q * r (ir) ) - cos (q * r ( &
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (SIN (q * r (ir) ) / (q * r (ir) ) - COS (q * r ( &
|
||||
ir) ) ) / (q * r (ir) )
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 2) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 2 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh)
|
||||
CALL vsin( sin_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = ( (3.d0 / (q * r (ir) ) - (q * r (ir) ) ) * sin_qr( &
|
||||
ir ) - 3.d0 * cos_qr(ir) ) / (q * r (ir) ) ** 2
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = ( (3.d0 / (q * r (ir) ) - (q * r (ir) ) ) * sin ( &
|
||||
q * r (ir) ) - 3.d0 * cos (q * r (ir) ) ) / (q * r (ir) ) ** &
|
||||
2
|
||||
enddo
|
||||
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = ( (3.d0 / (q * r (ir) ) - (q * r (ir) ) ) * SIN ( &
|
||||
q * r (ir) ) - 3.d0 * COS (q * r (ir) ) ) / (q * r (ir) ) ** 2
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 3) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 3 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh)
|
||||
CALL vsin( sin_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (sin_qr (ir) * (15.d0 / (q * r (ir) ) &
|
||||
- 6.d0 * (q * r (ir) ) ) + cos_qr (ir) * ( (q * r (ir) &
|
||||
) **2 - 15.d0) ) / (q * r (ir) ) **3
|
||||
enddo
|
||||
|
||||
ENDDO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = (sin (q * r (ir) ) * (15.d0 / (q * r (ir) ) &
|
||||
- 6.d0 * (q * r (ir) ) ) + cos (q * r (ir) ) * ( (q * r (ir) &
|
||||
) **2 - 15.d0) ) / (q * r (ir) ) **3
|
||||
enddo
|
||||
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (SIN (q * r (ir) ) * (15.d0 / (q * r (ir) ) &
|
||||
- 6.d0 * (q * r (ir) ) ) + COS (q * r (ir) ) * ( (q * r (ir) &
|
||||
) **2 - 15.d0) ) / (q * r (ir) ) ** 3
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 4) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 4 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh)
|
||||
CALL vsin( sin_qr, qr, msh)
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (sin_qr (ir) * (105.d0 - 45.d0 * (q * r (ir) &
|
||||
) **2 + (q * r (ir) ) **4) + cos_qr (ir) * (10.d0 * &
|
||||
(q * r (ir) ) **3 - 105.d0 * (q * r (ir) ) ) ) / (q * r (ir) &
|
||||
) **5
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = (sin (q * r (ir) ) * (105.d0 - 45.d0 * (q * r (ir) &
|
||||
) **2 + (q * r (ir) ) **4) + cos (q * r (ir) ) * (10.d0 * &
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (SIN (q * r (ir) ) * (105.d0 - 45.d0 * (q * r (ir) &
|
||||
) **2 + (q * r (ir) ) **4) + COS (q * r (ir) ) * (10.d0 * &
|
||||
(q * r (ir) ) **3 - 105.d0 * (q * r (ir) ) ) ) / (q * r (ir) &
|
||||
) **5
|
||||
enddo
|
||||
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 5) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 5 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh )
|
||||
CALL vsin( sin_qr, qr, msh )
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (-cos_qr(ir) - (945.d0*cos_qr(ir))/(q*r(ir)) ** 4 + &
|
||||
(105.d0*cos_qr(ir))/ (q*r(ir)) ** 2 + &
|
||||
(945.d0*sin_qr(ir))/ (q*r(ir)) ** 5 - &
|
||||
(420.d0*sin_qr(ir))/(q*r(ir)) ** 3 + &
|
||||
(15.d0*sin_qr(ir))/ (q*r(ir)))/(q*r(ir))
|
||||
end do
|
||||
|
||||
END DO
|
||||
!
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = (-cos(q*r(ir)) - (945.d0*cos(q*r(ir)))/(q*r(ir)) ** 4 + &
|
||||
(105.d0*cos(q*r(ir)))/ (q*r(ir)) ** 2 + &
|
||||
(945.d0*sin(q*r(ir)))/ (q*r(ir)) ** 5 - &
|
||||
(420.d0*sin(q*r(ir)))/(q*r(ir)) ** 3 + &
|
||||
(15.d0*sin(q*r(ir)))/ (q*r(ir)))/(q*r(ir))
|
||||
end do
|
||||
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = (-COS(q*r(ir)) - (945.d0*COS(q*r(ir)))/(q*r(ir)) ** 4 + &
|
||||
(105.d0*COS(q*r(ir)))/ (q*r(ir)) ** 2 + &
|
||||
(945.d0*SIN(q*r(ir)))/ (q*r(ir)) ** 5 - &
|
||||
(420.d0*SIN(q*r(ir)))/(q*r(ir)) ** 3 + &
|
||||
(15.d0*SIN(q*r(ir)))/ (q*r(ir)))/(q*r(ir))
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
|
||||
elseif (l == 6) then
|
||||
|
||||
#if defined __MASS
|
||||
|
||||
!
|
||||
ELSE IF ( l == 6 ) THEN
|
||||
!
|
||||
#if defined (__MASS)
|
||||
!
|
||||
qr = q * r
|
||||
call vcos( cos_qr, qr, msh)
|
||||
call vsin( sin_qr, qr, msh)
|
||||
do ir = ir0, msh
|
||||
!
|
||||
CALL vcos( cos_qr, qr, msh )
|
||||
CALL vsin( sin_qr, qr, msh )
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = ((-10395.d0*cos_qr(ir))/(q*r(ir)) ** 5 + &
|
||||
(1260.d0*cos_qr(ir))/(q*r(ir)) ** 3 - &
|
||||
(21.d0*cos_qr(ir))/ (q*r(ir)) - sin_qr(ir) + &
|
||||
(10395.d0*sin_qr(ir))/(q*r(ir)) ** 6 - &
|
||||
(4725.d0*sin_qr(ir))/ (q*r(ir)) ** 4 + &
|
||||
(210.d0*sin_qr(ir))/(q*r(ir)) ** 2)/(q*r(ir))
|
||||
end do
|
||||
|
||||
#else
|
||||
|
||||
do ir = ir0, msh
|
||||
jl (ir) = ((-10395.d0*cos(q*r(ir)))/(q*r(ir)) ** 5 + &
|
||||
(1260.d0*cos(q*r(ir)))/(q*r(ir)) ** 3 - &
|
||||
(21.d0*cos(q*r(ir)))/ (q*r(ir)) - sin(q*r(ir)) + &
|
||||
(10395.d0*sin(q*r(ir)))/(q*r(ir)) ** 6 - &
|
||||
(4725.d0*sin(q*r(ir)))/ (q*r(ir)) ** 4 + &
|
||||
(210.d0*sin(q*r(ir)))/(q*r(ir)) ** 2)/(q*r(ir))
|
||||
end do
|
||||
|
||||
#endif
|
||||
|
||||
else
|
||||
|
||||
call errore ('sph_bes', 'not implemented', l)
|
||||
|
||||
endif
|
||||
|
||||
endif
|
||||
END DO
|
||||
!
|
||||
return
|
||||
end subroutine sph_bes
|
||||
#else
|
||||
!
|
||||
DO ir = ir0, msh
|
||||
jl (ir) = ((-10395.d0*COS(q*r(ir)))/(q*r(ir)) ** 5 + &
|
||||
(1260.d0*COS(q*r(ir)))/(q*r(ir)) ** 3 - &
|
||||
(21.d0*COS(q*r(ir)))/ (q*r(ir)) - SIN(q*r(ir)) + &
|
||||
(10395.d0*SIN(q*r(ir)))/(q*r(ir)) ** 6 - &
|
||||
(4725.d0*SIN(q*r(ir)))/ (q*r(ir)) ** 4 + &
|
||||
(210.d0*SIN(q*r(ir)))/(q*r(ir)) ** 2)/(q*r(ir))
|
||||
END DO
|
||||
!
|
||||
#endif
|
||||
!
|
||||
ELSE
|
||||
!
|
||||
CALL errore ('sph_bes', 'not implemented', l)
|
||||
!
|
||||
END IF
|
||||
!
|
||||
END IF
|
||||
!
|
||||
RETURN
|
||||
!
|
||||
END SUBROUTINE sph_bes
|
||||
|
|
|
@ -0,0 +1,109 @@
|
|||
!
|
||||
! Copyright (C) 2001-2004 PWSCF-FPMD-CP90 group
|
||||
! This file is distributed under the terms of the
|
||||
! GNU General Public License. See the file `License'
|
||||
! in the root directory of the present distribution,
|
||||
! or http://www.gnu.org/copyleft/gpl.txt .
|
||||
!
|
||||
!----------------------------------------------------------------------------
|
||||
SUBROUTINE sph_dbes( MMAX, R, XG, L, DJL )
|
||||
!----------------------------------------------------------------------------
|
||||
!
|
||||
! ... calculates derivatives of spherical bessel functions j_l(Gr)
|
||||
! ... with respect to h_alpha,beta (without the factor GAGK(KK,IG)*HTM1)
|
||||
! ... i.e. -x * D(jl(x))/dx
|
||||
!
|
||||
USE kinds, ONLY : DP
|
||||
USE constants, ONLY : eps8
|
||||
!
|
||||
IMPLICIT NONE
|
||||
!
|
||||
INTEGER :: MMAX, L
|
||||
REAL(KIND=DP) :: XG
|
||||
REAL(KIND=DP) :: DJL(MMAX), R(MMAX)
|
||||
!
|
||||
INTEGER :: IR
|
||||
REAL(KIND=DP) :: XRG, XRG2
|
||||
!
|
||||
!
|
||||
IF ( L == 1 ) THEN ! S PART
|
||||
IF( XG < eps8 ) THEN
|
||||
DO IR=1,MMAX
|
||||
DJL(IR) = 0.D0
|
||||
END DO
|
||||
ELSE
|
||||
DJL(1) = 0.D0
|
||||
DO IR=2,MMAX
|
||||
XRG=R(IR)*XG
|
||||
DJL(IR) = SIN(XRG)/XRG-COS(XRG)
|
||||
END DO
|
||||
ENDIF
|
||||
ENDIF
|
||||
!
|
||||
IF ( L == 2 ) THEN ! P PART
|
||||
IF( XG < eps8 ) THEN
|
||||
DO IR=1,MMAX
|
||||
DJL(IR) = 0.D0
|
||||
END DO
|
||||
ELSE
|
||||
DJL(1) = 0.D0
|
||||
DO IR=2,MMAX
|
||||
XRG=R(IR)*XG
|
||||
DJL(IR) = 2.D0*(SIN(XRG)/XRG-COS(XRG))/XRG - SIN(XRG)
|
||||
END DO
|
||||
ENDIF
|
||||
ENDIF
|
||||
!
|
||||
IF ( L == 3 ) THEN ! D PART
|
||||
IF ( XG < eps8 ) THEN
|
||||
DO IR=1,MMAX
|
||||
DJL(IR) = 0.D0
|
||||
END DO
|
||||
ELSE
|
||||
DJL(1) = 0.D0
|
||||
DO IR=2,MMAX
|
||||
XRG=R(IR)*XG
|
||||
DJL(IR) = ( SIN(XRG)*(9.D0/(XRG*XRG)-4.D0) - &
|
||||
9.D0*COS(XRG)/XRG ) /XRG + COS(XRG)
|
||||
END DO
|
||||
END IF
|
||||
END IF
|
||||
!
|
||||
IF ( L == 4 ) THEN ! F PART
|
||||
IF ( XG < eps8 ) THEN
|
||||
DO IR=1,MMAX
|
||||
DJL(IR) = 0.D0
|
||||
END DO
|
||||
ELSE
|
||||
DJL(1) = 0.D0
|
||||
DO IR=2,MMAX
|
||||
XRG=R(IR)*XG
|
||||
XRG2=XRG*XRG
|
||||
DJL(IR) = SIN(XRG)*(60.D0/(XRG2*XRG2)-27.D0/XRG2+1.d0) - &
|
||||
COS(XRG)*(60.D0/XRG2-7.D0)/XRG
|
||||
END DO
|
||||
END IF
|
||||
END IF
|
||||
!
|
||||
IF ( L == 5 ) THEN ! G PART
|
||||
IF ( XG < eps8 ) THEN
|
||||
DO IR=1,MMAX
|
||||
DJL(IR) = 0.D0
|
||||
END DO
|
||||
ELSE
|
||||
DJL(1) = 0.D0
|
||||
DO IR=2,MMAX
|
||||
XRG=R(IR)*XG
|
||||
XRG2=XRG*XRG
|
||||
DJL(IR) = SIN(XRG)*(525.D0/(XRG2*XRG2)-240.D0/XRG2+11.D0)/XRG - &
|
||||
COS(XRG)*(525.D0/(XRG2*XRG2)-65.D0/XRG2+1.D0)
|
||||
END DO
|
||||
END IF
|
||||
END IF
|
||||
!
|
||||
IF ( L <= 0 .OR. L >= 6 ) &
|
||||
CALL errore( 'sph_dbes', ' L NOT PROGRAMMED, L= ',L )
|
||||
!
|
||||
RETURN
|
||||
!
|
||||
END SUBROUTINE sph_dbes
|
Loading…
Reference in New Issue