Line data    Source code

       1             :       MODULE m_dsphbs
       2             :       IMPLICIT NONE
       3             : !********************************************************************
       4             : !     calculates the derivative of the spherical bessel functions
       5             : !     dfj(l) = d jl(x)/dx
       6             : !     for l=0,lmax and argument x
       7             : !     note that the spherical bessel functions fj(l), l=0,lmax are
       8             : !     needed (call sphbes to generate these)
       9             : !                                            m. weinert
      10             : !********************************************************************
      11             :       CONTAINS
      12     4625450 :       SUBROUTINE dsphbs(
      13     4625450 :      >                  lmax,x,fj,
      14     4625450 :      <                  dfj)
      15             : !$acc routine
      16             : !     ..
      17             : !     .. Arguments ..
      18             :       INTEGER, INTENT  (IN) :: lmax
      19             :       REAL,    INTENT  (IN) :: x
      20             :       REAL,    INTENT  (IN) ::  fj(0:lmax)
      21             :       REAL,    INTENT (OUT) :: dfj(0:lmax)
      22             : !
      23             : !     .. Parameters ..
      24             :       REAL,    PARAMETER :: xlim = 1.0e-04
      25             : !     ..
      26             : !     .. Local Scalars ..
      27             :       REAL fac,x2
      28             :       INTEGER l
      29             : 
      30     4625450 :       dfj(0) = -fj(1)
      31             : !--->    small x limit
      32     4625450 :       IF (x.LT.xlim) THEN
      33        2022 :          x2 = 0.5*x*x
      34        2022 :          fac = 1./3.
      35       18514 :          DO l = 1,lmax
      36       16492 :             dfj(l) = fac* (l-x2* (l+2)/ (2*l+3))
      37       18514 :             fac = x*fac/ (2*l+3)
      38             :          ENDDO
      39             :       ELSE
      40             : !--->    obtain dfj using recurrence relationship
      41    42237776 :          DO l = 1,lmax
      42    42237776 :             dfj(l) = fj(l-1) - (l+1)*fj(l)/x
      43             :          ENDDO
      44             :       END IF
      45             : 
      46     4625450 :       RETURN
      47             :       END SUBROUTINE dsphbs
      48             :       END MODULE m_dsphbs