1float2bf : true$
2true;
3
4/*  Rodrigue's formulae */
5
6(jacobi_p_rod(n,a,b,x) := block([ an, rho, g ],
7    an : (-1)^n * 2^n * n!,
8    rho : (1-t)^a * (1 + t)^b,
9    g : 1-t^2,
10    rat(subst(x,t, diff(rho * g^n,t,n) / (an * rho)))),
11
12 gen_laguerre_rod(n, a, x) := block([an, rho, g],
13   an : n!,
14   rho : exp(-t) * t^a,
15   g : t,
16  rat(subst(x,t, diff(rho * g^n,t,n) / (an * rho)))),
17
18 hermite_rod(n,x) := block([an, rho, g],
19   an : (-1)^n,
20   rho : exp(-t^2),
21   g : 1,
22  rat(subst(x,t, diff(rho * g^n,t,n) / (an * rho)))),
23
24/*See A&S 10.1.25 page 439. */
25
26 spherical_bessel_j_rod(n,x) := block([sofar,k],
27  sofar : sin(x)/x,
28  for k : 1 thru n do (
29     sofar : -diff(sofar,x) / x
30  ),
31  x^n * sofar),
32
33 spherical_bessel_y_rod(n,x) := block([sofar,k],
34  sofar : cos(x) / x,
35  for k : 1 thru n do (
36     sofar : -diff(sofar,x) / x
37  ),
38  -x^n * sofar),
39
40 all_functions :   [jacobi_p,
41  ultraspherical,
42  assoc_legendre_p,
43  legendre_q,
44  assoc_legendre_q,
45  chebyshev_t,
46  chebyshev_u,
47  laguerre,
48  gen_laguerre,
49  legendre_p,
50  hermite,
51  spherical_hankel2,
52  spherical_hankel1,
53  spherical_bessel_j,
54  spherical_bessel_y,
55  assoc_leg_cos,
56  spherical_harmonic],
57
58
59 errors_found : [ ],
60 tests_pass : [ ],
61
62 zerop(e) := is(0=e),
63
64 check_zero_list(e) := block([ k, okay,n],
65   kill(labels),
66   okay : true,
67   k : 0,
68   n : length(e),
69   while okay and k < n do (
70      k : k + 1,
71      if not zerop(e[ k ]) then (
72          okay : false
73      )
74  ),
75  if okay then (
76      tests_pass : endcons(test_name, tests_pass),
77      print("okay:  ", test_name)
78  ) else (
79         print("error: ", test_name),
80         print("should vanish = ", e[ k ]),
81         errors_found : endcons(test_name, errors_found)
82  )
83 ),
84
85 check_true_list(e) := block([ ],
86   if member(false, e) then  (
87       print("error:  ",  test_name),
88       errors_found : endcons(test_name, errors_found)
89   )  else (
90      print("okay:  ", test_name),
91      tests_pass : endcons(test_name, tests_pass)
92   )
93 ),
94
950);
960;
97
98/* Jacobi Rodrigues test*/
99
100(print (test_name :  "Jacobi Rodrigues test"),
101 makelist(jacobi_p(k,a,b,x) - jacobi_p_rod(k,a,b,x), k,0,5),
102 rat(%%));
103''(makelist (0, k, 0, 5));
104
105/*  gen_laguerre Rodrigues*/
106
107(print (test_name : "gen_laguerre Rodrigues"),
108 makelist(gen_laguerre(k,a,x) - gen_laguerre_rod(k,a,x), k,0,5),
109 rat(%%));
110''(makelist (0, k, 0, 5));
111
112/*  hermite Rodrigues*/
113
114(print (test_name : "hermite Rodrigues"),
115 makelist(hermite(k,x) - hermite_rod(k,x), k,0,5),
116 rat(%%));
117''(makelist (0, k, 0, 5));
118
119/*  spherical_bessel_j Rodrigues*/
120
121(print (test_name : "spherical_bessel_j Rodrigues"),
122 makelist(spherical_bessel_j(k,x) - spherical_bessel_j_rod(k,x), k,0,5),
123 ratsimp(exponentialize(%%)));
124''(makelist (0, k, 0, 5));
125
126/*  spherical_bessel_y Rodrigues*/
127
128(print (test_name : "spherical_bessel_y Rodrigues"),
129 makelist(spherical_bessel_y(k,x) - spherical_bessel_y_rod(k,x), k,0,5),
130 ratsimp(exponentialize(%%)));
131''(makelist (0, k, 0, 5));
132
133/*---------------*/
134
135(print (test_name : "spherical harmonic orthogonality"),
136
137 f(n1,m1,n2,m2) := defint(defint(trigreduce(spherical_harmonic(n1,m1,th,ph) *
138   conjugate(spherical_harmonic(n2,m2,th,ph)) * sin(th)),th,0,%pi),ph,0,2*%pi),
139 0);
1400;
141
142block ([sofar : [ ]],
143       for n1 : 0 thru 2 do (
144          for m1 : -n1 thru n1 do (
145               for n2 : 0 thru 2 do (
146	           for  m2 : - n2 thru n2 do (
147	                sofar : cons(f(n1,m1,n2,m2) - kron_delta(n1,n2) * kron_delta(m1,m2), sofar))))),
148       if every (zerop, sofar) then true else sofar);
149true;
150
151/*----------------*/
152
153/* See A&S 22.3.14 page 776 and 22.5.4 page 777.
154*/
155
156(print (test_name : "A&S 22.3.14"),
157 makelist(ultraspherical(n,a,cos(x))/a - 2*cos(n*x)/n,n,1, 3),
158 limit(%%,a,0),
159 ratsimp(trigreduce(%%))),
160''(makelist (0, 3));
161
162/* See A&S 22.3.14 page 776
163*/
164
165(print (test_name : "A&S 22.3.14"),
166 makelist(chebyshev_t(n,cos(x)) - cos(n*x),n,0, 3),
167 ratsimp(trigreduce(%%)))$
168''(makelist (0, 4));
169
170/* See A&S 22.3.15 page 776
171*/
172
173(print (test_name : "A&S 22.3.15"),
174 makelist(sin(x) * chebyshev_u(n,cos(x)) - sin((n+1)*x),n,0, 4),
175 ratsimp(trigreduce(%%)))$
176''(makelist (0, 5));
177
178/* See A&S 22.7.15  page 782.
179*/
180
181(print (test_name : "A&S 22.7.15"),
182 jacobi_rec(n) := (n + a/2+b/2+1)*(1-x)*jacobi_p(n,a+1,b,x)
183      - (n+a+1)*jacobi_p(n,a,b,x) + (n+1)*jacobi_p(n+1,a,b,x),
184 makelist(ratsimp(jacobi_rec(n)),n,0,  7));
185''(makelist (0, 8));
186
187/* See  A&S 22.7.16 page 782
188*/
189
190(print (test_name : "A&S 22.7.16"),
191 jacobi_rec(n) := (n + a/2+b/2+1)*(1+x)*jacobi_p(n,a,b+1,x)
192      - (n+b+1)*jacobi_p(n,a,b,x) - (n+1)*jacobi_p(n+1,a,b,x),
193 makelist(ratsimp(jacobi_rec(n)),n,0, 7));
194''(makelist (0, 8));
195
196/* See A&S 22.7.17 page 782
197*/
198
199(print (test_name : "A&S 22.7.17"),
200 jacobi_rec(n) := (1-x)*jacobi_p(n,a+1,b,x)
201      + (1+x)*jacobi_p(n,a,b+1,x) - 2*jacobi_p(n,a,b,x),
202 makelist(ratsimp(jacobi_rec(n)),n,0, 7));
203''(makelist (0, 8));
204
205/* See A&S 22.7.18 page 782
206*/
207
208test_name : "A&S 22.7.18"$
209
210jacobi_rec(n) := (2*n+a+b)*jacobi_p(n,a-1,b,x)
211      - (n+a+b)*jacobi_p(n,a,b,x) +  (n+b)*jacobi_p(n-1,a,b,x)$
212
213check_zero_list(makelist(ratsimp(jacobi_rec(n)),n,1, 7))$
214
215/* See A&S 22.7.19 page 782
216*/
217
218test_name : "A&S 22.7.19"$
219
220jacobi_rec(n) := (2*n+a+b)*jacobi_p(n,a,b-1,x)
221      - (n+a+b)*jacobi_p(n,a,b,x) -  (n+a)*jacobi_p(n-1,a,b,x)$
222
223check_zero_list(makelist(ratsimp(jacobi_rec(n)),n, 1, 7))$
224
225/* See A&S 22.7.20 page 782
226*/
227
228test_name : "A&S 22.7.20"$
229
230jacobi_rec(n) := jacobi_p(n,a,b-1,x)
231      - jacobi_p(n,a-1,b,x) -  jacobi_p(n-1,a,b,x)$
232
233check_zero_list(makelist(ratsimp(jacobi_rec(n)),n,1, 2))$
234
235/* See A&S 22.7.21 page 782
236*/
237
238test_name : "A&S 22.7.21"$
239
240ultraspherical_rec(n) := 2*a*(1-x^2)*ultraspherical(n-1,a+1,x)
241- (2*a+n-1) * ultraspherical(n-1,a,x) + n*x*ultraspherical(n,a,x)$
242
243foo : makelist(ratsimp(ultraspherical_rec(n)),n,1, 7)$
244foo : append(ev(foo,a=4/5), ev(foo, a= -2/3), ev(foo, a=8))$
245check_zero_list(foo)$
246
247
248/* See A&S 22.7.23 page 782;  funny 22.7.22 is missing a lhs.
249    Maxima lacks simplification rules to simplify linear
250    combinations of gamma functions.
251*/
252
253test_name : "A&S 22.7.23"$
254
255ultraspherical_rec(n) := (n+a)*ultraspherical(n+1,a-1,x)
256- (a-1) * (ultraspherical(n+1,a,x) - ultraspherical(n-1,a,x))$
257
258foo : makelist(ratsimp(ultraspherical_rec(n)),n,1, 7)$
259foo : append(ev(foo,a=4/5,rat), ev(foo, a=8,rat))$
260check_zero_list(foo)$
261
262
263/* See A&S 22.7.24 page 782
264*/
265
266test_name : "A&S 22.7.24"$
267
268cheb_rec(n,m) := 2 * chebyshev_t(m,x) * chebyshev_t(n,x)
269- chebyshev_t(n+m,x) - chebyshev_t(n-m,x)$
270
271foo : makelist(makelist(cheb_rec(n, m), m,1,n),n,1, 9)$
272foo : rat(foo)$
273foo : apply(append,foo)$
274check_zero_list(foo)$
275
276
277/* See A&S 22.7.25 page 782
278*/
279
280test_name : "A&S 22.7.25"$
281
282cheb_rec(n,m) := 2*(x^2-1)*chebyshev_u(m-1,x) * chebyshev_u(n-1,x)
283 - chebyshev_t(n+m,x) + chebyshev_t(n-m,x)$
284
285foo : makelist(makelist(cheb_rec(n, m), m,1,n),n,1, 9)$
286foo : rat(foo)$
287foo : apply(append,foo)$
288check_zero_list(foo)$
289
290
291/* See A&S 22.7.26 page 782
292*/
293
294test_name : "A&S 22.7.26"$
295
296cheb_rec(n,m) := 2*chebyshev_t(m,x) * chebyshev_u(n-1,x)
297   -chebyshev_u(n+m-1,x) - chebyshev_u(n-m-1,x)$
298
299foo : makelist(makelist(cheb_rec(n, m), m, 0, n-1),n,1, 9)$
300foo : rat(foo)$
301foo : apply(append,foo)$
302check_zero_list(foo)$
303
304/* See A&S 22.7.27 page 782
305*/
306
307test_name : "A&S 22.7.27"$
308
309cheb_rec(n,m) := 2*chebyshev_t(n,x) * chebyshev_u(m-1,x)
310   -chebyshev_u(n+m-1,x) + chebyshev_u(n-m-1,x)$
311
312foo : makelist(makelist(cheb_rec(n, m), m,1,n-1),n,  2,  10)$
313foo : rat(foo)$
314foo : apply(append,foo)$
315check_zero_list(foo)$
316
317/* See A&S 22.7.28 page 782
318*/
319
320test_name : "A&S 22.7.28"$
321
322cheb_rec(n) := 2*chebyshev_t(n,x) * chebyshev_u(n-1,x)
323   -chebyshev_u(2*n-1,x)$
324
325foo : makelist(cheb_rec(n), n,2,10)$
326foo : rat(foo)$
327check_zero_list(foo)$
328
329/* See A&S 22.7.29 page 783
330*/
331
332test_name : "A&S 22.7.29"$
333
334lag_rec(n) := gen_laguerre(n,1/3+1,x) -((x-n) * gen_laguerre(n,1/3,x) + (1/3+n) * gen_laguerre(n-1,1/3,x))/x$
335check_zero_list(ratsimp(makelist(lag_rec(n),n,1, 10)))$
336
337lag_rec(n) := gen_laguerre(n,5+1,x) -((x-n) * gen_laguerre(n,5,x) + (5+n) * gen_laguerre(n-1,5,x))/x$
338check_zero_list(ratsimp(makelist(lag_rec(n),n,1,10)))$
339
340lag_rec(n) := gen_laguerre(n,1,x) -((x-n) * gen_laguerre(n,0,x) + (n) * gen_laguerre(n-1,0,x))/x$
341check_zero_list(ratsimp(makelist(lag_rec(n),n,1,10)))$
342
343
344/* See A&S 22.7.30 page 783
345*/
346
347test_name : "A&S 22.7.30"$
348
349declare(a,integer)$
350
351lag_rec(n) := gen_laguerre(n,a-1,x) -  gen_laguerre(n,a,x) +  gen_laguerre(n-1,a,x)$
352
353check_zero_list(ratsimp(makelist(lag_rec(n),n,1, 10)))$
354remove(a,integer)$
355
356 /* See A&S 22.7.31 page 783
357 */
358
359test_name : "A&S 22.7.31"$
360
361lag_rec(n) := gen_laguerre(n,a+1,x) -  (n+a+1)*gen_laguerre(n,a,x) / x +  (n+1)*gen_laguerre(n+1,a,x)/x$
362check_zero_list(ratsimp(makelist(lag_rec(n),n,0, 7)))$
363
364 /* See A&S 22.7.32 page 783
365 */
366
367test_name : "A&S 22.7.32"$
368
369lag_rec(n) := gen_laguerre(n,a-1,x) -  (n+1)*gen_laguerre(n+1,a,x) / (n+a) +  (n+1-x)*gen_laguerre(n,a,x)/(n+a)$
370check_zero_list(ratsimp(makelist(lag_rec(n),n,0,7)))$
371
372/* See A&S 22.2.1 page 774
373*/
374
375test_name : "A&S 22.2.1"$
376
377f(a,b,n,m) := integrate((1-x)^a*(1+x)^b * jacobi_p(n,a,b,x)
378                                        * jacobi_p(m,a,b,x),x,-1,1)$
379foo : makelist(makelist(f(1/2,-1/2,n,m), n, 0, m - 1), m, 1, 5)$
380foo : apply(append,foo)$
381check_zero_list(foo)$
382
383foo : makelist(makelist(f(1/3,2/3,n,m), n, 0, m - 1), m, 1, 5)$
384foo : apply(append,foo)$
385check_zero_list(foo)$
386
387
388test_name : "jacobi normalization"$
389
390stand(n) := expand(jacobi_p(n,a,b,1) - binomial(n+a,n))$
391check_zero_list(makelist(stand(n),n,0, 7))$
392
393test_name : "A&S 22.2.3"$
394
395assume(a > 1/2)$
396baz(n,m) := 'integrate((1-x^2)^(a-1/2) * ultraspherical(n,a,x)* ultraspherical(m,a,x),x,-1,1)$
397foo : makelist(makelist(ev(baz(n,m),integrate),n,0,m-1),m,1, 2)$
398foo : ev(foo, rat)$
399forget(a > 1/2)$
400foo : apply(append,foo)$
401check_zero_list(foo)$
402
403/* See A&S 22.2.3 page 774; Maxima doesn't know enough about
404the gamma functions to simplify these expressions to zero.
405*/
406
407test_name : "A&S 22.2.3"$
408
409stand(n) := ultraspherical(n,a,1) - binomial(n+2*a-1,n)$
410foo : makelist(stand(n),n,0, 7)$
411foo : append(ev(foo,a=2/3), ev(foo, a=7), ev(foo,a=1/3))$
412foo : rat(foo)$
413check_zero_list(foo)$
414
415/* See A&S 22.2.10 page 774
416*/
417
418test_name : "A&S 22.2.10"$
419
420f(n,m) := integrate(legendre_p(n,(rat(x))) * legendre_p(m,rat(x)),x,-1,1) - kron_delta(n,m) * 2 /(2*n+1)$
421foo : makelist(makelist(f(n,m),n,0,3),m,0,3)$
422foo : apply(append, foo)$
423check_zero_list(foo)$
424
425test_name : "legendre poly normalization"$
426stand(n) := legendre_p(n,1) -1$
427foo : makelist(stand(n),n,0,35)$
428check_zero_list(foo)$
429
430/* See A&S 22.2.12 page 774
431*/
432
433test_name : "A&S 22.2.12"$
434baz(n,m) := 'integrate(gen_laguerre(n,a,x) * gen_laguerre(m,a,x) * exp(-x)*x^a,x,0,inf) - kron_delta(n,m) * gamma(a+n+1)/n!$
435foo : makelist(makelist(baz(n,m),n,0,5),m,0,5)$
436foo : apply(append, ev(foo,a=2/3, integrate))$
437check_zero_list(foo)$
438
439
440/* See A&S 22.2.13 page 775
441*/
442
443test_name : "A&S 22.2.13"$
444baz(n,m) := 'integrate(laguerre(n,x) * laguerre(m,x) * exp(-x),x,0,inf) - kron_delta(n,m)$
445foo : makelist(makelist(ev(baz(n,m), integrate), n,0,4),m,0,4)$
446foo : rat(foo)$
447foo : apply(append, foo)$
448check_zero_list(foo)$
449
450/* See A&S 22.2.14 page 775
451*/
452
453test_name : "A&S 22.2.14"$
454baz(n,m) := 'integrate(hermite(n,x) * hermite(m,x) * exp(-x^2),x,-inf,inf) - kron_delta(n,m) * sqrt(%pi) * 2^n * n!$
455foo : makelist(makelist(ev(baz(n,m), integrate), n,0,4),m,0,4)$
456foo : apply(append, foo)$
457check_zero_list(foo)$
458
459
460/* Some Christoffel-Darboux sum formulae
461*/
462
463/* See A&S 22.12.2 page 785
464*/
465
466test_name : "A&S 22.12.2"$
467baz(n) := sum(chebyshev_t(2*k,x) ,k,0,n)  - (1 + chebyshev_u(2*n,x))/2$
468foo : makelist(rat(baz(n)),n,0,7)$
469check_zero_list(foo)$
470
471/* See A&S 22.12.3 page 785
472*/
473
474test_name : "A&S 22.12.3"$
475baz(n) := sum(chebyshev_t(2*k+1,x) ,k,0,n-1)  - chebyshev_u(2*n-1,x)/2$
476foo : makelist(rat(baz(n)),n,1, 7)$
477check_zero_list(foo)$
478
479/* See A&S 22.12.4 page 785
480*/
481
482test_name : "A&S 22.12.4"$
483baz(n) := sum(chebyshev_u(2*k,x) ,k,0,n)  - (1-chebyshev_t(2*n+2,x))/(2 * (1-x^2))$
484foo : makelist(rat(baz(n)),n,1, 11)$
485check_zero_list(foo)$
486
487/* See A&S 22.2.12.5 page 785
488*/
489
490test_name : "A&S 22.12.5"$
491baz(n) := sum(chebyshev_u(2*k+1,x) ,k,0,n-1)  - (x-chebyshev_t(2*n+1,x))/(2 * (1-x^2))$
492foo : makelist(rat(baz(n)),n,1,7)$
493check_zero_list(foo)$
494
495/* See A&S 22.12.6 page 785
496*/
497
498test_name : "A&S 22.12.6"$
499baz(n) := gen_laguerre(n,a+b+1,x+y) - sum(gen_laguerre(k,a,x) * gen_laguerre(n-k,b,y),k,0,n)$
500foo : makelist(rat(baz(n)),n,1, 7)$
501check_zero_list(foo)$
502
503/* See A&S 22.12.7 page 785
504*/
505
506test_name : "A&S 22.12.7"$
507baz(n) := gen_laguerre(n,a,x*y) -
508sum(binomial(n+a,m) * y^(n-m) * (1-y)^m * gen_laguerre(n-m,a,x),m,0,n)$
509foo : makelist(rat(baz(n)),n,1,3)$
510check_zero_list(foo)$
511
512/* See A&S 22.12.8 page 785
513*/
514
515test_name : "A&S 22.12.8"$
516baz(n) := hermite(n,x+y) - sum(binomial(n, k) * hermite(k,sqrt(2)*x) * hermite(n-k,sqrt(2)*y),k,0,n) / 2^(n/2)$
517foo : makelist(rat(baz(n)),n,0, 7)$
518check_zero_list(foo)$
519
520/* See A&S 22.5.17 page 778
521*/
522
523test_name : "A&S 22.5.17"$
524baz(n,m) := gen_laguerre(n,m,x) - (-1)^m * diff(laguerre(n+m,x),x,m)$
525foo : makelist(makelist(baz(i,j),j,0,i),i,0,7)$
526foo : rat(foo)$
527foo : apply(append,foo)$
528check_zero_list(foo)$
529
530/* See A&S 22.7.29 page 783
531*/
532
533test_name : "A&S 22.7.29"$
534baz(n) := gen_laguerre(n,a+1,x) - ((x-n)*gen_laguerre(n,a,x) + (a+n)*gen_laguerre(n-1,a,x))/x$
535foo : makelist(baz(i),i,1,8)$
536foo : rat(foo)$
537check_zero_list(foo)$
538
539/* float tests
540*/
541
542test_name : "A&S 22.12.2 - float"$
543orthopoly_returns_intervals : true$
544listarith : true$
545
546xargs(e) := if intervalp(e) then args(e) else [e,0]$
547
548f(n,x) := block([p,q],
549   p : sum(xargs(chebyshev_t(2*i,x)),i,0,n),
550   q : xargs(chebyshev_u(2*n,x)) / 2,
551   [part(p,1) - part(q,1), part(p,2) + part(q,2)])$
552
553gomer : []$
554for i : -50 thru 50 do (
555  foo : f(10, 1.25 * i),
556  e : part(foo,2),
557  foo : part(foo,1),
558  if (abs(foo - 0.5) <= e) then foo : 0 else foo : 1,
559  gomer : cons(foo, gomer))$
560check_zero_list(gomer)$
561
562gomer : []$
563for i : -100 thru 100 do (
564  foo : f(35, 0.01 * i),
565  e : part(foo,2),
566  foo : part(foo,1),
567  if (abs(foo - 0.5) <= e) then foo : 0 else foo : 1,
568  gomer : cons(foo, gomer))$
569check_zero_list(gomer)$
570
571gomer : []$
572for i : 1 thru 100 do (
573  foo : f(35, -0.01 * i),
574  e : part(foo,2),
575  foo : part(foo,1),
576  if (abs(foo - 0.5) <= e) then foo : 0 else foo : 1,
577  gomer : cons(foo, gomer))$
578check_zero_list(gomer)$
579
580gomer : []$
581for i : -100 thru 100 do (
582  foo : f(35, 0.01 * i + %i * 0.2),
583  e : part(foo,2),
584  foo : part(foo,1),
585  if (abs(foo - 0.5) <= e) then foo : 0 else foo : 1,
586  gomer : cons(foo, gomer))$
587check_zero_list(gomer)$
588
589gomer : []$
590for i : -100 thru 100 do (
591  foo : f(36, -0.01 * i + %i * 0.2),
592  e : part(foo,2),
593  foo : part(foo,1),
594  if (abs(foo - 0.5) <= e) then foo : 0 else foo : 1,
595  gomer : cons(foo, gomer))$
596check_zero_list(gomer)$
597
598test_name : "A&S 22.12.3 - float"$
599
600f(n,x) := block([p,q],
601   p : sum(xargs(chebyshev_t(2*i+1,x)),i,0,n-1),
602   q : xargs(chebyshev_u(2*n-1,x)) / 2,
603   [part(p,1) - part(q,1), part(p,2) + part(q,2)])$
604
605
606gomer : []$
607for i : -50 thru 50 do (
608  foo : f(10, 1.25 * i),
609  e : part(foo,2),
610  foo : part(foo,1),
611  if (abs(foo) <= e) then foo : 0 else foo : 1,
612  gomer : cons(foo, gomer))$
613check_zero_list(gomer)$
614
615gomer : []$
616for i : -100 thru 100 do (
617  foo : f(35, 0.01 * i),
618  e : part(foo,2),
619  foo : part(foo,1),
620  if (abs(foo) <= e) then foo : 0 else foo : 1,
621  gomer : cons(foo, gomer))$
622check_zero_list(gomer)$
623
624gomer : []$
625for i : 1 thru 100 do (
626  foo : f(35, -0.01 * i),
627  e : part(foo,2),
628  foo : part(foo,1),
629  if (abs(foo) <= e) then foo : 0 else foo : 1,
630  gomer : cons(foo, gomer))$
631check_zero_list(gomer)$
632
633gomer : []$
634for i : -100 thru 100 do (
635  foo : f(35, 0.01 * i + %i * 0.2),
636  e : part(foo,2),
637  foo : part(foo,1),
638  if (abs(foo) <= e) then foo : 0 else foo : 1,
639  gomer : cons(foo, gomer))$
640check_zero_list(gomer)$
641
642gomer : []$
643for i : -100 thru 100 do (
644  foo : f(36, -0.01 * i + %i * 0.2),
645  e : part(foo,2),
646  foo : part(foo,1),
647  if (abs(foo) <= e) then foo : 0 else foo : 1,
648  gomer : cons(foo, gomer))$
649check_zero_list(gomer)$
650
651
652test_name : "G&R 8.974-3-float"$
653
654f(n,a,x) := block([p,q],
655   p : sum(xargs(gen_laguerre(i,a,x)),i,0,n),
656   q : xargs(gen_laguerre(n,a+1,x)),
657   [part(p,1) - part(q,1), part(p,2) + part(q,2)])$
658
659gomer : []$
660for i : -100 thru 100 do (
661  foo : f(11,0.4, float(0.01 * i)),
662  e : part(foo,2),
663  foo : part(foo,1),
664  if (abs(foo) <= e) then foo : 0 else foo : 1,
665  gomer : cons(foo,gomer))$
666check_zero_list(gomer)$
667
668gomer : []$
669for i : -100 thru 100 do (
670  foo : f(12,0.4, float(0.01 * i)),
671  e : part(foo,2),
672  foo : part(foo,1),
673  if (abs(foo) <= e) then foo : 0 else foo : 1,
674  gomer : cons(foo,gomer))$
675check_zero_list(gomer)$
676
677gomer : []$
678for i : -100 thru 100 do (
679  foo : f(12,-0.4, float(0.01 * i)),
680  e : part(foo,2),
681  foo : part(foo,1),
682  if (abs(foo) <= e) then foo : 0 else foo : 1,
683  gomer : cons(foo,gomer))$
684check_zero_list(gomer)$
685
686gomer : []$
687for i : -100 thru 100 do (
688  foo : f(12,-41.0, float(0.01 * i)),
689  e : part(foo,2),
690  foo : part(foo,1),
691  if (abs(foo) <= e) then foo : 0 else foo : 1,
692  gomer : cons(foo,gomer))$
693check_zero_list(gomer)$
694
695gomer : []$
696for i : -100 thru 100 do (
697  foo : f(12,-4.1, float(0.01 * i * %i)),
698  e : part(foo,2),
699  foo : part(foo,1),
700  if (abs(foo) <= e) then foo : 0 else foo : 1,
701  gomer : cons(foo,gomer))$
702check_zero_list(gomer)$
703
704test_name : "random jacobi float"$
705
706gomer : [ ]$
707for i : 1 thru 500 do (
708   n : random(10),
709   a : random(11) / (1 + random(9)),
710   b : random(11) / (1 + random(9)),
711   x : (random(11) - 5) / (1 + random(10)),
712   exact : float(jacobi_p(n,a,b,x)),
713   approx : jacobi_p(n,float(a),float(b),float(x)),
714   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
715   if (foo = 1) then print ("n = ",n, " a = ",a," b = ",b, "x = ",x),
716   gomer : cons(foo,gomer))$
717check_zero_list(gomer)$
718
719test_name : "random ultraspherical float"$
720
721gomer : [ ]$
722for i : 1 thru 500 do (
723   n : random(10),
724   a : random(11) / (1 + random(9)),
725   x : (random(11) - 5) / (1 + random(10)),
726   exact : float(ultraspherical(n,a,x)),
727   approx : ultraspherical(n,float(a),float(x)),
728   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
729   if (foo = 1) then print ("n = ",n, " a = ",a,"x = ",x),
730   gomer : cons(foo,gomer))$
731check_zero_list(gomer)$
732
733test_name : "random chebyshev_t float"$
734gomer : [ ]$
735for i : 1 thru 500 do (
736   n : random(10),
737   x : (random(11) - 5) / (1 + random(10)),
738   exact : float(chebyshev_t(n,x)),
739   approx : chebyshev_t(n,float(x)),
740   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
741   if (foo = 1) then print("chebyshev_t float error; n = ",n, "x = ",x),
742   gomer : cons(foo,gomer))$
743check_zero_list(gomer)$
744
745test_name : "random chebyshev_u float"$
746gomer : [ ]$
747for i : 1 thru 500 do (
748   n : random(10),
749   x : (random(11) - 5) / (1 + random(10)),
750   exact : float(chebyshev_u(n,x)),
751   approx : chebyshev_u(n,float(x)),
752   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
753   if (foo = 1) then print("chebyshev_t float error; n = ",n, "x = ",x),
754   gomer : cons(foo,gomer))$
755check_zero_list(gomer)$
756
757
758test_name : "random legendre float"$
759gomer : [ ]$
760for i : 1 thru 500 do (
761   n : random(10),
762   x : (random(11) - 5) / (1 + random(10)),
763   exact : float(legendre_p(n,x)),
764   approx : legendre_p(n,float(x)),
765   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
766   if (foo = 1) then print("legendre float error; n = ",n, "x = ",x),
767   gomer : cons(foo,gomer))$
768check_zero_list(gomer)$
769
770test_name : "random hermite float"$
771gomer : [ ]$
772for i : 1 thru 500 do (
773   n : random(10),
774   x : (random(11) - 5) / (1 + random(10)),
775   exact : float(hermite(n,x)),
776   approx : hermite(n,float(x)),
777   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
778   if (foo = 1) then print("hermite float error; n = ",n, "x = ",x),
779   gomer : cons(foo,gomer))$
780check_zero_list(gomer)$
781
782
783test_name : "random gen_laguerre float"$
784
785gomer : [ ]$
786for i : 1 thru 500 do (
787   n : random(10),
788   a : random(11) / (1 + random(9)),
789   x : (random(11) - 5) / (1 + random(10)),
790   exact : float(gen_laguerre(n,a,x)),
791   approx : gen_laguerre(n,float(a),float(x)),
792   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
793   if (foo = 1) then print ("n = ",n, " a = ",a," b = ",b, "x = ",x),
794   gomer : cons(foo,gomer))$
795check_zero_list(gomer)$
796
797test_name : "random assoc_legendre_p float"$
798gomer : [ ]$
799for i : 1 thru 500 do (
800   n : random(10),
801   m : random(10),
802   x : (random(11) - 5)/ (1 + random(10)),
803   exact : float(assoc_legendre_p(n,m,x)),
804   approx : assoc_legendre_p(n,m,float(x)),
805   if not(intervalp(approx)) then approx : interval(approx,0),
806   if (abs(exact-part(approx,1)) <= part(approx,2)) then foo : 0 else foo : 1,
807   if (foo = 1) then print ("n =", n, " m = ", m, " x = ",x),
808   gomer : cons(foo,gomer))$
809check_zero_list(gomer)$
810
811remvalue(n)$
812remvalue(x)$
813remvalue(k)$
814q : []$
815
816test_name : "simple spherical_bessel_j test"$
817
818fpprec : 75$
819foo : ratsimp(spherical_bessel_j(5,x))$
820q : makelist(ev(foo / bfloat(spherical_bessel_j(5,0.1 * k)) - 1.0b0,
821   x = 0.1b0 * k),k,1,10)$
822q : map(lambda([s], if abs(s) < 5.0b-11 then 0 else 1), q)$
823check_zero_list(q)$
824
825foo : ratsimp(spherical_bessel_j(5,x))$
826q : makelist(ev(foo / bfloat(spherical_bessel_j(5,0.1 * k)) - 1.0b0,
827   x = 0.1b0 * k),k,-10,-1)$
828q : map(lambda([s], if abs(s) < 5.0b-11 then 0 else 1), q)$
829check_zero_list(q)$
830
831test_name : "simple spherical_bessel_y test"$
832
833q : []$
834foo : expand(spherical_bessel_y(5,x))$
835q : makelist(ev(foo / bfloat(spherical_bessel_y(5,0.1 * k)) - 1.0b0,
836   x = 0.1b0*k),k,1,10)$
837q : map(lambda([s], if abs(s) < 5.0b-11 then 0 else 1), q)$
838check_zero_list(q)$
839
840q : []$
841foo : spherical_bessel_y(5,x)$
842q : makelist(ev(foo / bfloat(spherical_bessel_y(5,0.1 * k)) - 1.0b0,
843   x = 0.1b0*k),k,-10,-1)$
844q : map(lambda([s], if abs(s) < 5.0b-11 then 0 else 1), q)$
845check_zero_list(q)$
846
847
848remvalue(n)$
849remvalue(i)$
850remvalue(x)$
851remvalue(a)$
852remvalue(b)$
853remvalue(gomer)$
854remvalue(exact)$
855remvalue(approx)$
856remfunction(f)$
857
858/* See A & S 8.5.4
859*/
860
861test_name : "A&S 8.5.3"$
862
863foo(n) := (n + 1)*legendre_q(n+1,x) - (2*n+1)*x*legendre_q(n,x)
864    + n*legendre_q(n-1,x)$
865gomer  : makelist(rat(foo(i)),i,1,15)$
866check_zero_list(gomer)$
867
868test_name : "A&S 8.5.3"$
869foo(n) :=  (n + 1)*legendre_p(n+1,x) - (2*n+1)*x*legendre_p(n,x)
870    + n*legendre_p(n-1,x)$
871gomer : makelist(rat(foo(i)),i,1,15)$
872check_zero_list(gomer)$
873
874test_name : "A&S 8.5.3"$
875
876foo(n,m) := (n - m + 1) * assoc_legendre_q(n+1,m,x) -
877   (2*n + 1) * x * assoc_legendre_q(n,m,x) + (n + m) * assoc_legendre_q(n-1,m,x)$
878
879gomer : makelist(makelist(foo(n,m),m,-n+1,n-1),n,1,5)$
880gomer : map(ratsimp,apply(append, gomer))$
881check_zero_list(gomer)$
882
883
884
885/* See A&S 8.6.7 page 334
886*/
887
888test_name : "A&S 8.6.7"$
889foo : makelist(makelist(assoc_legendre_q(n,m,x) - (-1)^m * (1-x^2)^(m/2) * diff(legendre_q(n,x),x,m),m,0,n),n,0,4)$
890foo : apply(append,foo)$
891foo : ratsimp(foo)$
892check_zero_list(foo)$
893
894
895/* See G&R 8.810 page 1014
896*/
897
898test_name : "G&R 8.810"$
899foo : makelist(makelist(assoc_legendre_p(n,m,x) - (-1)^m * (1-x^2)^(m/2) * diff(legendre_p(n,x),x,m),n,0,5),m,0,5)$
900foo : apply(append,foo)$
901foo : rat(foo)$
902check_zero_list(foo)$
903
904
905/* See G&R 8.813 page 1015
906*/
907
908test_name : "G&R 8.813 (1-6)"$
909
910foo : [assoc_legendre_p(1,1,x) + sqrt(1-x^2),
911         assoc_legendre_p(2,1,x) + 3 * x * sqrt(1-x^2),
912         assoc_legendre_p(2,2,x) - 3 *(1-x^2),
913         assoc_legendre_p(3,1,x) + 3 * sqrt(1-x^2) *(5*x^2-1) / 2,
914         assoc_legendre_p(3,2,x) - 15*x*(1-x^2),
915         assoc_legendre_p(3,3,x) + 15 * (1-x^2)^(3/2)]$
916
917foo :  rat(foo)$
918
919check_zero_list(foo)$
920
921/* See G&R 8.950 (1) page 1033
922*/
923
924test_name : "G&R 8.950 (1)"$
925foo : makelist(hermite(n,x) - (-1)^n * exp(x^2) * diff(exp(-x^2),x,n),n,0,9)$
926foo : rat(foo)$
927check_zero_list(foo)$
928
929
930/* See G&R 8.952 (1) page 1033
931*/
932
933test_name : "G&R 8.952 (1)"$
934foo : makelist(diff(hermite(n,x),x) - 2 * n * hermite(n-1,x),n,1,7)$
935foo : rat(foo)$
936check_zero_list(foo)$
937
938/* See G&R 8.952 (2) page 1033
939*/
940test_name : "G&R 8.952 (2)"$
941foo : makelist(hermite(n+1,x) - 2 * x * hermite(n,x) + 2 * n * hermite(n-1,x),n,1,7)$
942foo : rat(foo)$
943check_zero_list(foo)$
944
945/* See G&R 8.956 (1-3) page 1034
946*/
947
948test_name : "G&R 8.956 (1-5)"$
949
950foo : [hermite(0,x) - 1,
951         hermite(1,x) - 2*x,
952	 hermite(2,x) - (4*x^2 -2),
953	 hermite(3,x) - (8*x^3-12*x),
954	 hermite(4,x) - (16*x^4 - 48 * x^2 + 12)]$
955foo : rat(foo)$
956check_zero_list(foo)$
957
958/* See A&S 10.1.19 page 439
959*/
960
961test_name : "A&S 10.1.19 spherical_hankel2"$
962baz(n) := spherical_hankel2(n-1,x) + spherical_hankel2(n+1,x) - (2*n+1) * spherical_hankel2(n,x) /x$
963foo : makelist(baz(k),k,-7,7)$
964foo : rat(expand(foo))$
965check_zero_list(foo)$
966
967test_name : "A&S 10.1.19 spherical_hankel1"$
968baz(n) := spherical_hankel1(n-1,x) + spherical_hankel1(n+1,x) - (2*n+1) * spherical_hankel1(n,x) /x$
969foo : makelist(baz(k),k,-7,7)$
970foo : rat(expand(foo))$
971check_zero_list(foo)$
972
973test_name : "A&S 10.1.19 spherical_bessel_j"$
974baz(n) := spherical_bessel_j(n-1,x) + spherical_bessel_j(n+1,x) - (2*n+1) * spherical_bessel_j(n,x) /x$
975foo : makelist(baz(k),k,-7,7)$
976foo : rat(foo)$
977check_zero_list(foo)$
978
979test_name : "A&S 10.1.19 spherical_bessel_y"$
980baz(n) := spherical_bessel_y(n-1,x) + spherical_bessel_y(n+1,x) - (2*n+1) * spherical_bessel_y(n,x) /x$
981foo : makelist(baz(k),k,-7,7)$
982foo : rat(foo)$
983check_zero_list(foo)$
984
985/*--------------*/
986
987/* See A&S 10.1.20 page 439
988*/
989
990test_name : "A&S 10.1.20 spherical_hankel1"$
991remvalue(q)$
992baz(n) := n * spherical_hankel1(n-1,q) - (n+1)*spherical_hankel1(n+1,q)
993 -(2*n+1) * diff(spherical_hankel1(n,q),q)$
994foo : makelist(baz(k),k,-7,7)$
995foo : ratsimp(foo)$
996check_zero_list(foo)$
997
998test_name : "A&S 10.1.20 spherical_hankel2"$
999baz(n) := n * spherical_hankel2(n-1,q) - (n+1)*spherical_hankel2(n+1,q)
1000 -(2*n+1) * diff(spherical_hankel2(n,q),q)$
1001foo : makelist(baz(k),k,-7,7)$
1002foo : ratsimp(foo)$
1003check_zero_list(foo)$
1004
1005test_name : "A&S 10.1.20 spherical_bessel_j"$
1006baz(n) := n * spherical_bessel_j(n-1,q) - (n+1)*spherical_bessel_j(n+1,q)
1007 -(2*n+1) * diff(spherical_bessel_j(n,q),q)$
1008foo : makelist(baz(k),k,-7,7)$
1009foo : ratsimp(foo)$
1010check_zero_list(foo)$
1011
1012test_name : "A&S 10.1.20 spherical_bessel_y"$
1013baz(n) := n * spherical_bessel_y(n-1,q) - (n+1)*spherical_bessel_y(n+1,q)
1014 -(2*n+1) * diff(spherical_bessel_y(n,q),q)$
1015foo : makelist(baz(k),k,-7,7)$
1016foo : ratsimp(foo)$
1017check_zero_list(foo)$
1018
1019/*--------------*/
1020
1021kill(labels)$
1022
1023test_name : "A&S 10.1.31"$
1024f(n) := spherical_bessel_j(n,t) *spherical_bessel_y(n-1,t) -
1025    spherical_bessel_j(n-1,t) *spherical_bessel_y(n,t) - 1/t^2$
1026
1027foo : makelist(f(k),k,-3, 3)$
1028foo : trigreduce(ratsimp(foo))$
1029check_zero_list(foo)$
1030
1031/*--------------*/
1032
1033remove(n,integer)$
1034remove(k,integer)$
1035remove(i,integer)$
1036remvalue(a,b,x,n,m,i,j)$
1037test_name : "jacobi_p gradef test"$
1038q : [ ]$
1039for i : 0 thru 15 do (
1040   foo : diff(jacobi_p(n,a,b,x^2) - jacobi_p(i,a,b,x^2),x),
1041   foo : ev(foo, n=i),
1042   foo : rat(foo),
1043   q : cons(foo,q))$
1044check_zero_list(q)$
1045
1046/*--------------*/
1047
1048test_name : "ultraspherical gradef test"$
1049q : [ ]$
1050for i : 0 thru 15 do (
1051   foo : diff(ultraspherical(n,a,x^2) - ultraspherical(i,a,x^2),x),
1052   foo : ev(foo, n = i),
1053   foo : rat(foo),
1054   q : cons(foo, q))$
1055check_zero_list(q)$
1056
1057/*--------------*/
1058
1059test_name : "assoc_legendre_p gradef test"$
1060remvalue(q,i,j,n,m,x,foo)$
1061q : []$
1062for i : 0 thru 15 do (
1063   for j : -15 thru 15 do (
1064      foo : diff(assoc_legendre_p(n,m,x) - assoc_legendre_p(i,j,x),x),
1065      foo : ev(foo, n=i,m=j),
1066      foo : radcan(foo),
1067      q : cons(foo,q)))$
1068check_zero_list(q)$
1069
1070/*--------------*/
1071
1072test_name : "assoc_legendre_q gradef test"$
1073remvalue(q,i,j,n,m,x,foo,w)$
1074q : []$
1075for i : 0 thru 5 do (
1076   for j : -i thru i do (
1077      w : assoc_legendre_q(i,j,x),
1078      w : diff(w,x),
1079      foo : diff(assoc_legendre_q(n, m, x),x) - w,
1080      foo : ev(foo, n=i,m=j),
1081      foo : ratsimp(expand(foo)),
1082      q : cons(foo,q)))$
1083check_zero_list(q)$
1084
1085/*--------------*/
1086
1087test_name : "chebyshev_t gradef test"$
1088q : []$
1089for i : 0 thru 15 do (
1090   foo : diff(chebyshev_t(n,x^2) - chebyshev_t(i,x^2),x),
1091   foo : ev(foo, n = i),
1092   foo : rat(foo),
1093   q : cons(foo,q))$
1094check_zero_list(q)$
1095
1096/*--------------*/
1097
1098test_name : "chebyshev_u gradef test"$
1099q : []$
1100for i : 0 thru 15 do (
1101    foo : diff(chebyshev_u(n,x^2) - chebyshev_u(i,x^2),x),
1102    foo : ev(foo, n=i),
1103    foo : rat(foo),
1104    q : cons(foo,q))$
1105check_zero_list(q)$
1106
1107/*--------------*/
1108
1109test_name : "laguerre gradef test"$
1110q : [ ]$
1111for i : 0 thru 15 do (
1112    foo : diff(laguerre(n,x^2) - laguerre(i,x^2),x),
1113    foo : ev(foo, n=i),
1114    foo : rat(foo),
1115    q : cons(foo,q))$
1116check_zero_list(q)$
1117
1118/*--------------*/
1119
1120test_name : "gen_laguerre gradef test"$
1121q : [ ]$
1122for i : 0 thru 15 do (
1123    foo : diff(gen_laguerre(n,a,x^2) - gen_laguerre(i,a,x^2),x),
1124    foo : ev(foo, n=i),
1125    foo : rat(foo),
1126    q : cons(foo,q))$
1127check_zero_list(q)$
1128
1129/*--------------*/
1130
1131test_name : "legendre_p gradef test"$
1132q : [ ]$
1133for i : 0 thru 15 do (
1134     foo : diff(legendre_p(n,x^2) - legendre_p(i,x^2),x),
1135     foo : ev(foo, n=i),
1136     foo : rat(foo),
1137     q : cons(foo,q))$
1138check_zero_list(q)$
1139
1140/*--------------*/
1141
1142test_name : "legendre_q gradef test"$
1143q : [ ]$
1144for i : 0 thru 15 do (
1145     foo : diff(legendre_q(n,x) - legendre_q(i,x),x),
1146     foo : ev(foo, n=i),
1147     foo : rat(foo),
1148     q : cons(foo,q))$
1149check_zero_list(q)$
1150
1151
1152test_name : "hermite gradef test"$
1153q : []$
1154for i : 0 thru 15 do (
1155   foo : diff(hermite(n,x^2) - hermite(i,x^2),x),
1156   foo : ev(foo, n=i),
1157   foo : rat(foo),
1158   q : cons(foo,q))$
1159check_zero_list(q)$
1160
1161/*--------------*/
1162
1163test_name : "spherical_hankel2 gradef test"$
1164q : []$
1165for i : 0 thru 15 do (
1166   foo : diff(spherical_hankel2(n,x^2) - spherical_hankel2(i,x^2),x),
1167   foo : ev(foo, n=i),
1168   foo : rat(foo),
1169   q : cons(foo,q))$
1170check_zero_list(q)$
1171
1172/*--------------*/
1173
1174test_name : "spherical_hankel1 gradef test"$
1175q : [ ]$
1176for i : 0 thru 15 do (
1177    foo : diff(spherical_hankel1(n,x^2) - spherical_hankel1(i,x^2),x),
1178    foo : ev(foo, n=i),
1179    foo : ratsimp(foo),
1180    q : cons(foo,q))$
1181check_zero_list(q)$
1182
1183/*--------------*/
1184
1185test_name : "spherical_bessel_j gradef test"$
1186q : []$
1187for i : 0 thru 15 do (
1188   foo : diff(spherical_bessel_j(n,x^2) - spherical_bessel_j(i,x^2),x),
1189   foo : ev(foo, n=i),
1190   foo : ratsimp(foo),
1191   q : cons(foo,q))$
1192check_zero_list(q)$
1193
1194/*--------------*/
1195
1196test_name : "spherical_bessel_y gradef test"$
1197q : []$
1198for i : 0 thru 15 do (
1199   foo : diff(spherical_bessel_y(n,x^2) - spherical_bessel_y(i,x^2),x),
1200   foo : ev(foo, n=i),
1201   foo : rat(foo),
1202   q : cons(foo,q))$
1203check_zero_list(q)$
1204
1205
1206/*--------------*/
1207
1208test_name : "spherical_harmonic gradef test"$
1209remvalue(q,i,j,n,m,x,foo)$
1210q : []$
1211for i : 0 thru 5 do (
1212   for j : -i thru i do (
1213 foo : [diff(spherical_harmonic(n, m, x, y) - spherical_harmonic(i, j, x, y),x),
1214          diff(spherical_harmonic(n, m, x, y) - spherical_harmonic(i,j, x, y),y)],
1215     foo : ev(foo, n=i, m=j),
1216     foo : radcan(foo),
1217     foo : trigreduce(rat(foo)),
1218     q : append(foo,q)))$
1219check_zero_list(q)$
1220remvalue(q)$
1221
1222/*--------------*/
1223
1224test_name : "jacobi sum representation"$
1225declare(n,integer)$
1226foo : jacobi_p(n,p,q,2/3)$
1227foo : ev(foo, sum, n=7)$
1228foo : rat(foo - jacobi_p(7,p,q,2/3))$
1229check_zero_list([foo])$
1230
1231/*--------------*/
1232
1233test_name : "ultraspherical sum representation"$
1234foo : ultraspherical(n,p,-2/3)$
1235foo : ev(foo, sum, n=2)$
1236foo : rat(foo - ultraspherical(2,p,-2/3))$
1237check_zero_list([foo])$
1238
1239/*--------------*/
1240test_name : "legendre_p sum representation"$
1241foo : legendre_p(n,1/8)$
1242foo : ev(foo, sum, n=8)$
1243foo : rat(foo - legendre_p(8,1/8))$
1244check_zero_list([foo])$
1245
1246/*----------------*/
1247
1248test_name : "chebyshev_t sum representation"$
1249foo : chebyshev_t(n,2)$
1250foo : ev(foo, sum, n=5)$
1251foo : rat(foo - chebyshev_t(5,2))$
1252check_zero_list([foo])$
1253
1254/*----------------*/
1255
1256test_name : "chebyshev_u sum representation"$
1257foo : chebyshev_u(n,-1/4)$
1258foo : ev(foo, sum, n=15)$
1259foo : rat(foo - chebyshev_u(15,-1/4))$
1260check_zero_list([foo])$
1261
1262/*---------------*/
1263
1264test_name : "laguerre sum representation"$
1265foo : laguerre(n,2/3)$
1266foo : ev(foo,sum,n=4)$
1267foo : rat(foo - laguerre(4,2/3))$
1268check_zero_list([foo])$
1269
1270/*---------------*/
1271
1272test_name : "generalized laguerre sum representation"$
1273foo : gen_laguerre(n,a,-2/3)$
1274foo : ev(foo,sum,n=4)$
1275foo : rat(foo - gen_laguerre(4,a,-2/3))$
1276check_zero_list([foo])$
1277
1278
1279remvalue([x,n,k])$
1280
1281test_name : "pochhammer-1"$
1282check_zero_list(ratsimp([pochhammer(x,0) - 1, pochhammer(x,1) - x,
1283pochhammer(x,2) - x*(x+1), pochhammer(x,5)/pochhammer(x,4) - (x+4)]))$
1284
1285test_name : "pochhammer-2"$
1286check_zero_list(makelist(ratsimp(pochhammer(x,-k) * pochhammer(1-x,k) - (-1)^k),k,-5,5))$
1287
1288test_name : "pochhammer-grad"$
1289foo : pochhammer(x,n)$
1290dfoo : diff(foo,x)$
1291goober : makelist(diff(ev(foo,n = k),x) - ev(dfoo,n = k),k,-5,5)$
1292check_zero_list(expand(subst([x = 7/2], goober)));
1293check_zero_list(expand(subst([x = -7/2], goober)));
1294
1295/*-----------------*/
1296
1297test_name : "unit_step"$
1298check_zero_list([unit_step(-2), unit_step(-1/9), unit_step(-1.2),
1299unit_step(-1.5b-2), unit_step(0), unit_step(0.0), unit_step(0.0b0),
1300unit_step(2)-1, unit_step(1/9)-1, unit_step(4.5)-1, unit_step(8.23b3)-1,
1301unit_step(x^2+1)-1,unit_step(exp(x))-1])$
1302
1303/*----------------*/
1304
1305test_name : "kron_delta"$
1306
1307check_zero_list([kron_delta(1,2), kron_delta(1,1)-1, kron_delta(x,rat(x)) - 1,
1308kron_delta(x,y)-kron_delta(y,x), kron_delta(1.0,1.0)-1])$
1309
1310
1311
1312/*----------------*/
1313
1314test_name : "jacobi_p recursion"$
1315
1316foo : orthopoly_recur(jacobi_p,[m,a,b,x])$
1317foo : rhs(foo)-lhs(foo)$
1318foo : makelist(ev(foo,m=k),k,1,6)$
1319foo : rat(foo)$
1320check_zero_list(foo)$
1321
1322/*-----------------*/
1323test_name : "ultraspherical recursion"$
1324
1325foo : orthopoly_recur(ultraspherical,[m,a,x])$
1326foo : rhs(foo)-lhs(foo)$
1327foo : makelist(ev(foo,m=k),k,1,6)$
1328foo : rat(foo)$
1329check_zero_list(foo)$
1330
1331/*-----------------*/
1332test_name : "chebyshev_t_recursion"$
1333
1334foo : orthopoly_recur(chebyshev_t,[m,x])$
1335foo : rhs(foo)-lhs(foo)$
1336foo : makelist(ev(foo,m=k),k,1,6)$
1337foo : rat(foo)$
1338check_zero_list(foo)$
1339
1340/*-----------------*/
1341test_name : "chebyshev_u recursion"$
1342
1343foo : orthopoly_recur(chebyshev_u,[m,x])$
1344foo : rhs(foo)-lhs(foo)$
1345foo : makelist(ev(foo,m=k),k,1,6)$
1346foo : rat(foo)$
1347check_zero_list(foo)$
1348
1349/*-----------------*/
1350test_name : "legendre_p recursion"$
1351
1352foo : orthopoly_recur(legendre_p,[m,x])$
1353foo : rhs(foo)-lhs(foo)$
1354foo : makelist(ev(foo,m=k),k,1,6)$
1355foo : rat(foo)$
1356check_zero_list(foo)$
1357
1358/*-----------------*/
1359test_name : "legendre_q recursion"$
1360foo : orthopoly_recur(legendre_q,[m,x])$
1361foo : rhs(foo)-lhs(foo)$
1362foo : makelist(ev(foo,m=k),k,1,6)$
1363foo : rat(foo)$
1364check_zero_list(foo)$
1365
1366/*-----------------*/
1367test_name : "assoc_legendre_p recursion"$
1368foo : orthopoly_recur(assoc_legendre_p,[m,1,x])$
1369foo : rhs(foo)-lhs(foo)$
1370foo : makelist(ev(foo,m=k),k,2,6)$
1371foo : rat(foo)$
1372check_zero_list(foo)$
1373
1374/*-----------------*/
1375test_name : "assoc_legendre_p recursion"$
1376foo : orthopoly_recur(assoc_legendre_p,[m,2,x])$
1377foo : rhs(foo)-lhs(foo)$
1378foo : makelist(ev(foo,m=k),k,3,6)$
1379foo : rat(foo)$
1380check_zero_list(foo)$
1381
1382/*-----------------*/
1383test_name : "assoc_legendre_q recursion"$
1384foo : orthopoly_recur(assoc_legendre_q,[m,1,x])$
1385foo : rhs(foo)-lhs(foo)$
1386foo : makelist(ev(foo,m=k),k,2,6)$
1387foo : rat(foo)$
1388check_zero_list(foo)$
1389
1390/*----------------*/
1391
1392test_name : "laguerre recursion"$
1393
1394foo : orthopoly_recur(laguerre,[m,x])$
1395foo : rhs(foo)-lhs(foo)$
1396foo : makelist(ev(foo,m=k),k,1,6)$
1397foo : rat(foo)$
1398check_zero_list(foo)$
1399
1400/*----------------*/
1401
1402test_name : "gen_laguerre recursion"$
1403
1404foo : orthopoly_recur(gen_laguerre,[m,a,x])$
1405foo : rhs(foo)-lhs(foo)$
1406foo : makelist(ev(foo,m=k),k,1,6)$
1407foo : rat(foo)$
1408check_zero_list(foo)$
1409
1410/*-----------------*/
1411test_name : "hermite recursion"$
1412
1413foo : orthopoly_recur(hermite,[m,x])$
1414foo : rhs(foo)-lhs(foo)$
1415foo : makelist(ev(foo,m=k),k,1,6)$
1416foo : rat(foo)$
1417check_zero_list(foo)$
1418
1419/*----------------*/
1420test_name : "spherical_bessel_j recursion"$
1421foo : orthopoly_recur(spherical_bessel_j,[m,x])$
1422foo : rhs(foo)-lhs(foo)$
1423foo : makelist(ev(foo,m=k),k,1,6)$
1424foo : rat(foo)$
1425check_zero_list(foo)$
1426
1427
1428/*---------------------------------------*/
1429test_name : "spherical_bessel_j recursion"$
1430foo : orthopoly_recur(spherical_bessel_y,[m,x])$
1431foo : rhs(foo)-lhs(foo)$
1432foo : makelist(ev(foo,m=k),k,1,6)$
1433foo : rat(foo)$
1434check_zero_list(foo)$
1435
1436/*---------------------------------------*/
1437test_name : "spherical_hankel1 recursion"$
1438foo : orthopoly_recur(spherical_hankel1,[m,x])$
1439foo : rhs(foo)-lhs(foo)$
1440foo : makelist(ev(foo,m=k),k,1,6)$
1441foo : rat(foo)$
1442check_zero_list(foo)$
1443
1444/*---------------------------------------*/
1445test_name : "spherical_hankel2 recursion"$
1446foo : orthopoly_recur(spherical_hankel2,[m,x])$
1447foo : rhs(foo)-lhs(foo)$
1448foo : makelist(ev(foo,m=k),k,1,6)$
1449foo : rat(foo)$
1450check_zero_list(foo)$
1451
1452/*---------------------------------------*/
1453
1454test_name : "jacobi_weight"$
1455foo : orthopoly_weight(jacobi_p,[n,2,3,x])$
1456foo : integrate(jacobi_p(5,2,3,x) * jacobi_p(4,2,3,x) * foo[1],x,foo[2],foo[3])$
1457check_zero_list([foo])$
1458
1459/*---------------------------------------*/
1460
1461test_name : "ultraspherical_weight"$
1462foo : orthopoly_weight(ultraspherical,[n,2,x])$
1463foo : integrate(ultraspherical(5,2,x) * ultraspherical(4,2,x)
1464	* foo[1],x,foo[2],foo[3])$
1465check_zero_list([foo])$
1466
1467/*---------------------------------------*/
1468
1469test_name : "chebyshev_t_weight"$
1470foo : orthopoly_weight(chebyshev_t,[n,x])$
1471foo : integrate(chebyshev_t(5,x) * chebyshev_t(4,x)
1472	* foo[1],x,foo[2],foo[3])$
1473check_zero_list([foo])$
1474
1475/*---------------------------------------*/
1476
1477test_name : "chebyshev_u_weight"$
1478foo : orthopoly_weight(chebyshev_u,[n,x])$
1479foo : integrate(chebyshev_u(5,x) * chebyshev_u(4,x)
1480	* foo[1],x,foo[2],foo[3])$
1481check_zero_list([foo])$
1482
1483/*---------------------------------------*/
1484
1485test_name : "legendre_p_weight"$
1486foo : orthopoly_weight(legendre_p,[n,x])$
1487foo : integrate(legendre_p(5,x) * legendre_p(4,x)
1488	* foo[1],x,foo[2],foo[3])$
1489check_zero_list([foo])$
1490
1491/*---------------------------------------*/
1492
1493test_name : "laguerre_weight"$
1494foo : orthopoly_weight(laguerre,[n,x])$
1495foo : integrate(laguerre(5,x) * laguerre(4,x)
1496	* foo[1],x,foo[2],foo[3])$
1497check_zero_list([foo])$
1498
1499/*---------------------------------------*/
1500
1501test_name : "gen_laguerre_weight"$
1502foo : orthopoly_weight(gen_laguerre,[n,1/2,x])$
1503foo : integrate(gen_laguerre(5,1/2,x) * gen_laguerre(4,1/2,x)
1504	* foo[1],x,foo[2],foo[3])$
1505check_zero_list([foo])$
1506
1507/*---------------------------------------*/
1508
1509test_name : "hermite_weight"$
1510foo : orthopoly_weight(hermite,[n,x])$
1511foo : integrate(hermite(5,x) * hermite(4,x)
1512	* foo[1],x,foo[2],foo[3])$
1513check_zero_list([foo])$
1514
1515
1516/*---------------------------------------*/
1517
1518orthopoly_returns_intervals : false$
1519
1520test_name : "legendre_p negative degree--symbolic argument"$
1521foo1 : makelist (legendre_p (-k, u), k, 1, 8);
1522foo2 : makelist (legendre_p (k - 1, u), k, 1, 8);
1523check_zero_list(foo1 - foo2)$
1524
1525test_name : "legendre_p negative degree--rational argument"$
1526foo1 : makelist (legendre_p (-k, 11/7), k, 1, 8);
1527foo2 : makelist (legendre_p (k - 1, 11/7), k, 1, 8);
1528check_zero_list(foo1 - foo2)$
1529
1530test_name : "legendre_p negative degree--float argument"$
1531foo1 : makelist (legendre_p (-k, float (17/16)), k, 1, 8);
1532foo2 : makelist (legendre_p (k - 1, float (17/16)), k, 1, 8);
1533foo : map (lambda ([a, b], is (a = b)), foo1, foo2);
1534check_true_list (foo);
1535
1536/*---------------------------------------*/
1537
1538test_name : "assoc_legendre_p negative degree--symbolic argument"$
1539foo1 : makelist (assoc_legendre_p (-k, 1, u), k, 1, 8);
1540foo2 : makelist (assoc_legendre_p (k - 1, 1, u), k, 1, 8);
1541check_zero_list(foo1 - foo2)$
1542
1543test_name : "assoc_legendre_p negative degree--rational argument"$
1544foo1 : makelist (assoc_legendre_p (-k, 1, 11/7), k, 1, 8);
1545foo2 : makelist (assoc_legendre_p (k - 1, 1, 11/7), k, 1, 8);
1546check_zero_list(foo1 - foo2)$
1547
1548test_name : "assoc_legendre_p negative degree--float argument"$
1549foo1 : makelist (assoc_legendre_p (-k, 1, float (17/16)), k, 1, 8);
1550foo2 : makelist (assoc_legendre_p (k - 1, 1, float (17/16)), k, 1, 8);
1551foo : map (lambda ([a, b], is (a = b)), foo1, foo2);
1552check_true_list (foo);
1553
1554reset (orthopoly_returns_intervals);
1555
1556
1557print("orthopoly version = ", get('orthopoly,'version))$
1558print("errors found = ", errors_found)$
1559print("number of tests passed =", length(tests_pass))$
1560
1561
1562/* Generate A&S Figure 22.4 page 776
1563*/
1564
1565
1566orthopoly_returns_intervals : false$
1567
1568foo : ['(ultraspherical(2,0.5,x)),
1569'(ultraspherical(3,0.5,x)),
1570'(ultraspherical(4,0.5,x)),
1571'(ultraspherical(5,0.5,x))]$
1572
1573plot2d(foo,[x,-1,1])$
1574
1575/* Generate A&S Figure 22.5 page 777
1576*/
1577
1578
1579foo : ['(ultraspherical(5, 0.2, x)),
1580'(ultraspherical(5, 0.4, x)),
1581'(ultraspherical(5, 0.6, x)),
1582'(ultraspherical(5, 0.8, x)),
1583'(ultraspherical(5, 1.0, x))]$
1584
1585plot2d(foo, [x,-0.8,0.8])$
1586
1587/* Generate A&S Figure 22.6 page 778
1588*/
1589
1590foo : ['(chebyshev_t(1,x)),
1591'(chebyshev_t(2,x)),
1592'(chebyshev_t(3,x)),
1593'(chebyshev_t(4,x)),
1594'(chebyshev_t(5,x))]$
1595
1596
1597plot2d(foo,[x,-1.0,1.0])$
1598
1599/* Generate A&S Figure 22.7 page 779
1600*/
1601
1602
1603foo : ['(chebyshev_u(1,x)),
1604'(chebyshev_u(2,x)),
1605'(chebyshev_u(3,x)),
1606'(chebyshev_u(4,x)),
1607'(chebyshev_u(5,x))]$
1608
1609
1610plot2d(foo,[x,-1.0,1.0])$
1611
1612reset (float2bf);
1613[float2bf];
1614