C
C
      SUBROUTINE  BESSY1 (X,BJ1,BY)
C
C     -------------------------------------------------------------
C
C       SUBROUTINE PURPOSE:
C       EVALUATION OF BESSEL FUNCTION OF 2ND KIND AND 1ST ORDER.
C
C     -------------------------------------------------------------
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
C
      PI = 3.14159265359D0
      GA = 0.57721566490D0
      IF (X.LT.12.0) THEN
         BY=-X/2.0
         BBY=DABS(BY)
         BM=0.0
         BK=0.0
         F=1.0
         DO 50 I=1,200
            BM=BM+1.0/I
            BBY=BBY*X*X/(4.0*I*(I+1.0))
            BY=BY+BBY*(2.0*BM+1.0/(I+1.0))*F
            IF (I.GT.2.AND.DABS(BK-BY).LT.1.D-9) THEN
               BBY=BY+2.0*((DLOG(X/2.0)+GA)*BJ1-1.0/X)
               BY=BBY/PI
               RETURN
            ELSE
               BK=BY
               F=-F
            ENDIF
50       CONTINUE
         BBY=BY+2.0*((DLOG(X/2.0)+GA)*BJ1-1.0/X)
         BY=BBY/PI
      ELSE
         BS1 = 1.0
         BS2 = 0.375/X
         BY1 = BS1
         BY2 = BS2
         BZ = DSQRT(2.0/(PI*X))
         F = -1.0
         DO 100 K = 1,100
            BS1 = BS1*(4*K+1)*(4*K-1)*(3-4*K)*(5-4*K)/
     .            (128.0*K*(2*K-1)*X*X)
            BY1 = BY1+F*BS1
            IF (ABS(BS1).LT.1.D-9) GO TO 150
            BS2 = BS2*(4*K+1)*(4*K+3)*(3-4*K)*(1-4*K)/
     .            (128.0*(2*K+1)*K*X*X)
            BY2 = BY2+F*BS2
            F = -F
100      CONTINUE
150      BY = BZ*(DCOS(X-0.75*PI)*BY2+DSIN(X-0.75*PI)*BY1)
      ENDIF
C
C     RETURN
      END
