PROGRAM D7R2
C	Driver for routine RAN1
C	Calculates pi statistically using volume of unit n-sphere
	PARAMETER(PI=3.1415926)
	DIMENSION IY(3),YPROB(3)
	FNC(X1,X2,X3,X4)=SQRT(X1**2+X2**2+X3**2+X4**2)
	IDUM=-1
	DO 11 I=1,3
		IY(I)=0
11	CONTINUE
	WRITE(*,'(1X,/,T15,A)') 'Volume of unit n-sphere, n=2,3,4'
	WRITE(*,'(1X,/,T3,A,T17,A,T26,A,T37,A)')
     *		'# points','pi','(4/3)*pi','(1/2)*pi^2'
	DO 14 J=1,15
		DO 12 K=2**(J-1),2**J
			X1=RAN1(IDUM)
			X2=RAN1(IDUM)
			X3=RAN1(IDUM)
			X4=RAN1(IDUM)
			IF(FNC(X1,X2,0.0,0.0).LT.1.0) IY(1)=IY(1)+1
			IF(FNC(X1,X2,X3,0.0).LT.1.0) IY(2)=IY(2)+1
			IF(FNC(X1,X2,X3,X4).LT.1.0) IY(3)=IY(3)+1
12		CONTINUE
		DO 13 I=1,3
			YPROB(I)=1.0*(2**(I+1))*IY(I)/(2**J)
13		CONTINUE
		WRITE(*,'(1X,I8,3F12.6)') 2**J,(YPROB(I),I=1,3)
14	CONTINUE
	WRITE(*,'(1X,/,T4,A,3F12.6,/)') 'actual',PI,4.0*PI/3.0,0.5*(PI**2)
	END