InputMode:=Word

PrecisionDigits:=20

LOAD("CARLSON.MTH")

"Use numerical checks from Numerical Algorithms 10(1995)13-26."

"B.C. Carlson Numerical computation of real or complex elliptic integrals"

"page 22 3. Numerical checks"

mrf:=[[1,2,0,1.3110287771461],[#i,-#i,0,1.8540746773014],[0.5,1,0,1.854074677~
3014],[#i-1,#i,0,0.79612586584234-#i*(1.2138566698365)],[2,3,4,0.584082841677~
15],[#i,-#i,2,1.0441445654064],[#i-1,#i,1-#i,0.93912050218619-#i*(0.532962520~
18635)]]

APPROX(MAX(VECTOR(APPROX(ABS(R_F(mrf SUB i SUB 1,mrf SUB i SUB 2,mrf SUB i SU~
B 3)-mrf SUB i SUB 4)),i,1,DIMENSION(mrf))))

"Check R_C."

mrc:=[[0,1/4,pi],[9/4,2,LN(2)],[0,#i,(1-#i)*(1.1107207345396)],[-#i,#i,1.2260~
849569072-#i*(0.34471136988768)],[1/4,-2,LN(2)/3],[#i,-1,0.77778596920447+#i*~
(0.19832484993429)]]

APPROX(MAX(VECTOR(APPROX(ABS(R_C(mrc SUB i SUB 1,mrc SUB i SUB 2)-mrc SUB i S~
UB 3)),i,1,DIMENSION(mrc))))

"I can add some more exact solutions."

mrc:=[[0,1/4,pi],[9/4,2,LN(2)],[2,1,-LN(SQRT(2)-1)],[-#i,#i,-LN(SQRT(2)-1)/2+~
pi/4-#i*(LN(SQRT(2)-1)/2+pi/4)],[1/4,-2,LN(2)/3],[#i,-1,SQRT(SQRT(2)/4-1/4)*A~
TAN(SQRT(SQRT(2)-1))-SQRT(SQRT(2)/16+1/16)*LN(-SQRT(2*SQRT(2)+2)+SQRT(2)+1)+#~
i*(SQRT(SQRT(2)/4+1/4)*ATAN(SQRT(SQRT(2)-1))+SQRT(SQRT(2)/16-1/16)*LN(-SQRT(2~
*SQRT(2)+2)+SQRT(2)+1))],[0,1,pi/2],[#i,#i+1,pi/4+#i*LN(SQRT(2)-1)/2]]

APPROX(MAX(VECTOR(APPROX(ABS(R_C(mrc SUB i SUB 1,mrc SUB i SUB 2)-mrc SUB i S~
UB 3)),i,1,DIMENSION(mrc))))

"Check R_J."

mrj:=[[0,1,2,3,0.77688623778582],[2,3,4,5,0.14297579667157],[2,3,4,-1+#i,0.13~
613945827771-#i*(0.38207561624427)],[#i,-#i,0,2,1.6490011662711],[-1+#i,-1-#i~
,1,2,0.9414835884122],[#i,-#i,0,1-#i,1.8260115229009+#i*(1.2290661908643)],[-~
1+#i,-1-#i,1,-3+#i,-0.61127970812028-#i*(1.0684038390007)],[-1+#i,-2-#i,-#i,-~
1+#i,1.8249027393704-#i*(1.2218475784827)],[2,3,4,-0.5,0.24723819703052],[2,3~
,4,-5,-0.12711230042964]]

APPROX(MAX(VECTOR(APPROX(ABS(R_J(mrj SUB i SUB 1,mrj SUB i SUB 2,mrj SUB i SU~
B 3,mrj SUB i SUB 4)-mrj SUB i SUB 5)),i,1,DIMENSION(mrj))))

"page 23 Check R_D."

mrd:=[[0,2,1,1.7972103521034],[2,3,4,0.16510527294261],[#i,-#i,2,0.6593385415~
422],[0,#i,-#i,1.270819627191+#i*(2.7811120159521)],[0,#i-1,#i,-1.85772354392~
39-#i*(0.96193450888839)],[-2-#i,-#i,-1+#i,1.8249027393704-#i*(1.221847578482~
7)]]

APPROX(MAX(VECTOR(APPROX(ABS(R_D(mrd SUB i SUB 1,mrd SUB i SUB 2,mrd SUB i SU~
B 3)-mrd SUB i SUB 4)),i,1,DIMENSION(mrd))))

"Check R_G."

mrg:=[[0,16,16,pi],[2,3,4,1.7255030280692],[0,#i,-#i,0.42360654239699],[0,#i-~
1,#i,0.44660591677018+#i*(0.70768352357515)],[-#i,#i-1,#i,0.36023392184473+#i~
*(0.40348623401722)],[0,0.0796,4,1.0284758090288]]

APPROX(MAX(VECTOR(APPROX(ABS(R_G(mrg SUB i SUB 1,mrg SUB i SUB 2,mrg SUB i SU~
B 3)-mrg SUB i SUB 4)),i,1,DIMENSION(mrg))))

"Check the Legendre elliptic integrals using numerical integration."

APPROX(MAX(VECTOR(MAX(VECTOR(ABS(APPROX(LEGENDRE_F(phi,m)-ELLIPTIC_F(phi,m),6)),phi,~
pi/16,pi/2,pi/16)),m,1/8,7/8,1/8)))

APPROX(MAX(VECTOR(MAX(VECTOR(ABS(APPROX(LEGENDRE_E(phi,m)-ELLIPTIC_E(phi,m),6)),phi,~
pi/16,pi/2,pi/16)),m,1/8,7/8,1/8)))

"This will take a long time."

APPROX(MAX(VECTOR(MAX(VECTOR(MAX(VECTOR(ABS(APPROX(LEGENDRE_PI(phi,m,n)-ELLIPTIC_PI(~
phi,m,n),6)),phi,pi/16,pi/2,pi/16)),m,1/8,7/8,1/8)),n,1/8,7/8,1/8)))