PROGRAM D15R3
C	Driver for routine RKQC
	EXTERNAL DERIVS
	PARAMETER(N=4)
	DIMENSION Y(N),DYDX(N),YSCAL(N)
	X=1.0
	Y(1)=BESSJ0(X)
	Y(2)=BESSJ1(X)
	Y(3)=BESSJ(2,X)
	Y(4)=BESSJ(3,X)
	DYDX(1)=-Y(2)
	DYDX(2)=Y(1)-Y(2)
	DYDX(3)=Y(2)-2.0*Y(3)
	DYDX(4)=Y(3)-3.0*Y(4)
	DO 11 I=1,N
		YSCAL(I)=1.0
11	CONTINUE
	HTRY=0.1
	WRITE(*,'(/1X,T8,A,T19,A,T31,A,T43,A)')
     *		'eps','htry','hdid','hnext'
	DO 12 I=1,15
		EPS=EXP(-FLOAT(I))
		CALL RKQC(Y,DYDX,N,X,HTRY,EPS,YSCAL,HDID,HNEXT,DERIVS)
		WRITE(*,'(2X,E12.4,F8.2,2X,2F12.6)') EPS,HTRY,HDID,HNEXT
12	CONTINUE
	END

	SUBROUTINE DERIVS(X,Y,DYDX)
	DIMENSION Y(1),DYDX(1)
	DYDX(1)=-Y(2)
	DYDX(2)=Y(1)-(1.0/X)*Y(2)
	DYDX(3)=Y(2)-(2.0/X)*Y(3)
	DYDX(4)=Y(3)-(3.0/X)*Y(4)
	RETURN
	END