C	TPOLY
	DIMENSION PL(26,26),PROD(26)
9000	FORMAT(' GENERATE LEGENDRE POLYNOMIALS')
9001	FORMAT(//' WHAT IS MAXIMUM DEGREE? ',$)
9002	FORMAT(I)
9003	FORMAT(/'  P(',I2,')')
9004	FORMAT(/' PRODUCT OF P(',I2,') BY P(',
	1I2,') INTEGRATES TO ',F10.6)
9005	FORMAT('  THEORY SAYS INTEGRAL IS',F12.8)
C
	TYPE 9000
1	TYPE 9001
	ACCEPT 9002,L
	IF(L) 9,9,3
3	M=L+1
	CALL LEGEND(PL,L)
	DO 2 I=1,M
	K=I-1
	TYPE 9003,K
	CALL PRNT(PL(1,I),K)
	NP=K+L
	IF(NP-25)4,4,2
4	CALL MTALG(PL(1,I),K,PL(1,M),L,PROD)
	AREA=PINT(PROD,NP,-1.,+1.)
	TYPE 9004,K,L,AREA
	IF(I-M)2,5,5
5	EL=L
	THEORY=2./(1.+2.*EL)
	TYPE 9005,THEORY
2	CONTINUE
	GO TO 1
9	CALL EXIT
	END
C
	FUNCTION PINT(C,NC,XLOW,XHIGH)
	DIMENSION C(26)
	NPLUS=NC+1
	J=NPLUS
	DO 101 I=1,NPLUS
	D=J
	C(J+1)=C(J)/D
101	J=J-1
	C(1)=0.
	PINT=CLPLY(XHIGH,C,NPLUS)-CLPLY(XLOW,C,NPLUS)
	RETURN
	END
C
	SUBROUTINE LEGEND(CL,LL)
C	COMPUTE THE COEFFICIENTS OF THE LEGENDRE
C	POLYNOMIALS BY RECURSION
	DIMENSION CL(26,26),TA(26),TB(26),TC(26),TD(26)
	CL(1,1)=1.
	CL(1,2)=0.
	CL(2,2)=1.
	IF(LL-1) 207,207,201
201	L=0
	ND=1
	TD(1)=0.
	NC=0
	N=LL+1
	DO 202 I=1,N
	TA(I)=0.
202	TB(I)=0.
203	K=L
	L=L+1
	IF(L-LL) 204,207,207
204	M=L+1
	N=L+2
	EL=L
	TD(2)=2.*EL+1.
	TC(1)=EL
	CALL MTALG(TC,NC,CL(1,L),K,TB)
	CALL MTALG(TD,ND,CL(1,M),L,TA)
	DO 205 I=1,N
205	CL(I,N)=(TA(I)-TB(I))/(EL+1.)
	GO TO 203
207	RETURN
	END
C
	SUBROUTINE PRNT(A,NA)
	DIMENSION A(26)
1300	FORMAT('  POLYNOMIAL IS')
1301	FORMAT(1X,F12.6,5H  X**,I3)
	TYPE 1300
	NPLUS=NA+1
	J=NPLUS
	K=NA
	DO 303 I=1,NPLUS
	IF(A(J)) 301,302,301
301	TYPE 1301,A(J),K
302	J=J-1
303	K=K-1
	RETURN
	END

