Merge branch 'bug-4296-expop-assign'
[maxima.git] / tests / rtest_taylor.mac
blob72dd58045ec9f6adfffa4281416983bc8aaea3b8
1 (kill(all),0);
2 0$
4 /*--- taylorinfo ---*/
6 /* SF bug 867310 */
8 taylorinfo(x);
9 false$
11 taylorinfo(a = b);
12 false$
14 taylorinfo(taylor(x,x,0,5));
15 [[x,0,5]]$
17 taylorinfo(taylor(x,[x,y],0,5));
18 [[x,y],[0,0],[5,5]]$
20 taylorinfo(taylor(x,[x,y],0,5,z,a,7));
21 [[x,y],[0,0],[5,5],[z,a,7]]$
23 (tp : taylor(x, [x, a, 5, asympt]),0);
26 taylorinfo(tp);
27 [[x,a,5, asympt]]$
29 taylorinfo(taylor(taylor(x,x,a,7),y,b,9));
30 [[y,b,9],[x,a,7]]$
32 taylorinfo(taylor(x,x,0,inf));
33 [[x,0,'inf]]$
35 /*--- constants ----*/
37 ratdisrep(taylor(0,x,1,2));
40 ratdisrep(taylor(-42,x,1,2));
41 -42$
43 ratdisrep(taylor(6.023e23,x,1,2)), keepfloat ;
44 6.023e23$
46 ratdisrep(taylor(sqrt(5),x,1,2));
47 sqrt(5)$
49 ratdisrep(taylor(sqrt(5) - %i,x,1,2));
50 sqrt(5) - %i$
52 taylor(x[2],x[1],0,2);
53 x[2]$
55 /*---polynomials ---*/
57 ratdisrep(taylor(x,x,1,2));
60 ratdisrep(taylor(x + sqrt(y),x,1,2));
61 x + sqrt(y)$
63 ratdisrep(taylor(x*(x+1) + sqrt(y),x,1,2));
64 sqrt(y)+(x-1)^2+3*(x-1)+2$
66 ratdisrep(taylor(1 + sqrt(6) * x + %pi*x^2,x,-%phi,9));
67 %pi*(x+%phi)^2+(sqrt(6)-2*%phi*%pi)*(x+%phi)+(%phi+1)*%pi-sqrt(6)*%phi+1$
69 radcan(taylor((1 + x*y) * (1 - x),[x,y],sqrt(2),3));
70 (x-x^2)*y-x+1$
72 expand(taylor(x*y + x^2 * y,[x,y],5,3));
73 x*y + x^2 * y$
75 taylor(rat(x*(x-1)),x,0,5);
76 x^2 - x$
79 /* algebraic functions ---*/
81 taylor(abs(x),x,-sqrt(2),2);
82 -x$
84 (tp : taylor(sqrt(1 + x),x,0,5),0);
87 ratdisrep(tp^2);
88 1+x$
90 (tp : taylor((1 + x)^(3/4),x,0,5),0);
93 ratdisrep(tp^(4/3));
94 1+x$
96 (tp : taylor((1 + x)^mu,x,0,5),0);
99 ratdisrep(tp^(1/mu));
100 1+x$
102 (tp : taylor((1 - a * x)^mu,x,0,5),0);
105 ratdisrep(tp^(1/mu));
106 1-a * x$
108 (tp : taylor((1 - x / sqrt(3))^(1/5),x,0,5),0);
111 ratdisrep(tp^5 - (1-(sqrt(3)*x)/3));
114 (tp : taylor(sqrt(x*(x-1)),x,1,5),0);
117 ratdisrep(tp^2 - x * (x-1));
120 (tp : taylor(sqrt((1-x)/(1-y)),[x,y],0,5),0);
123 (1-y)*tp^2;
124 1-x$
126 /*---log and trig like ----*/
128 taylor(cos(x),x,0,4);
129 1 - x^2/2! + x^4 / 4!$
131 taylor(sin(x),x,0,4);
132 x - x^3 / 3!$
134 taylor(tan(x),x,0,4);
135 x + x^3/3$
137 taylor(cosh(x),x,0,5);
138 1 + ((x^2)/2) + ((x^4)/24)$
140 taylor(sinh(x),x,0,5);
141 x + ((x^3)/6) + ((x^5)/120)$
143 taylor(tanh(x),x,0,5);
144 x - ((x^3)/3) + ((2 * x^5)/15)$
146 taylor(sin(x) / cos(x),x,0,10) - taylor(tan(x),x,0,10);
149 taylor(sin(x) / cos(x),x, %pi/6,10) - taylor(tan(x),x, %pi/6,10);
152 trigsimp(taylor(sin(x + x^2) / cos(x + x^2),x, %pi/6,10) - taylor(tan(x + x^2),x, %pi/6,10));
155 taylor(sec(x),x,0,5) * taylor(cos(x),x,0,5);
158 taylor(sec(x),x,%pi/3,5) * taylor(cos(x),x,%pi/3,5);
161 taylor(csc(x),x,0,5) * taylor(sin(x),x,0,5);
164 taylor(csc(x),x,%pi/4,5) * taylor(sin(x),x,%pi/4,5);
167 taylor(cot(x),x,0,5) * taylor(tan(x),x,0,5);
170 taylor(cot(x),x,12,5) * taylor(tan(x),x,12,5);
173 trigsimp(taylor(sinh(x) / cosh(x),x,0,10) - taylor(tanh(x),x,0,10));
176 trigsimp(taylor(sinh(x) / cosh(x),x, 1/6,10) - taylor(tanh(x),x, 1/6,10));
179 trigsimp(taylor(sinh(x + x^2) / cosh(x + x^2),x, %pi/6,10) - taylor(tanh(x + x^2),x, %pi/6,10));
182 taylor(sech(x),x,0,5) * taylor(cosh(x),x,0,5);
185 taylor(sech(x),x,%pi/3,5) * taylor(cosh(x),x,%pi/3,5);
188 taylor(csch(x),x,0,5) * taylor(sinh(x),x,0,5);
191 taylor(csch(x),x,%pi/4,5) * taylor(sinh(x),x,%pi/4,5);
194 taylor(coth(x),x,0,5) * taylor(tanh(x),x,0,5);
197 taylor(coth(x),x,12,5) * taylor(tanh(x),x,12,5);
200 taylor(log(x),x,1,3) - (x-1-(x-1)^2/2+(x-1)^3/3);
203 taylor(exp(x),x,0,7);
204 ''(sum(x^k / k!,k,0,7))$
206 taylor(exp(-x),x,0,7);
207 ''(sum((-x)^k / k!,k,0,7))$
209 log(taylor(exp(x),x,0,10));
212 exp(taylor(log(x),x,1,10));
215 (tp : taylor(x,x,0,5),0);
218 is(equal(taylor(cos(x),x,0,5), cos(tp)));
219 true$
221 is(equal(taylor(sin(x),x,0,5), sin(tp)));
222 true$
224 is(equal(taylor(tan(x),x,0,5), tan(tp)));
225 true$
227 is(equal(taylor(exp(x),x,0,5), exp(tp)));
228 true$
230 /*--- asymptotic ---*/
232 taylor(x + 1/x + 2/x^2,[x,inf, 2, asympt]);
233 x+1/x + 2/x^2$
235 taylor(x + 1/x + 2/x^2,[x, inf, 1, asympt]);
236 x + 1/x$
238 /*--- jacobian elliptic ---*/
240 /*  G&R 8.125 1 through 3  */
242 taylor(jacobi_sn(x,k),x,0,5);
243 x - (1 + k) * x^3 / 6 + (1 + 14*k + k^2) * x^5 / 120;
245 taylor(jacobi_cn(x,k),x,0,4);
246 1 - x^2 / 2 + (1 + 4 * k) * x^4 / 24$
248 taylor(jacobi_dn(x,k),x,0,4);
249 1 - (k * x^2 / 2) + (k^2 + 4 * k) * x^4 / 24$
252 /*  G&R  8.154 1,2, 4, and 5 */
254 taylor(jacobi_sn(x,k)^2,x,0,19) - taylor((1 - jacobi_cn(2*x,k))/(1 + jacobi_dn(2*x,k)),x,0,19);
257 taylor(jacobi_cn(x,k)^2,x,0,19)-taylor((jacobi_cn(2*x,k) + jacobi_dn(2*x,k))/(1 + jacobi_dn(2*x,k)),x,0,19);
260 taylor(jacobi_sn(x,k)^2,x,0,15) + taylor(jacobi_cn(x,k)^2,x,0,15);
263 taylor(jacobi_dn(x,k)^2,x,0,15) + k * taylor(jacobi_sn(x,k)^2,x,0,15);
266 taylor(inverse_jacobi_sn(jacobi_sn(x,k),k),x,0,4);
269 taylor(jacobi_sn(inverse_jacobi_sn(x,k),k),x,0,4);
272 /* The Jacobian elliptic function cd(u,m) = cn(u,m)/dn(u,m). */
274 taylor(jacobi_cd(x,k),x,0,8) - taylor(jacobi_cn(x,k) / jacobi_dn(x,k),x,0,8);
278 /*--- Bessel ---*/
280 taylor(bessel_j(0,x),x,0,5);
281 1 - x^2/4 + x^4/64$
283 taylor(bessel_j(1,x),x,0,5);
284 x/2 - x^3/16 + x^5/384$
286 (assume(m > 0),0);
289 taylor(bessel_j(m,x),x,0,1);
290 x^m /(2^m * gamma(m + 1))$
292 taylor(bessel_i(0,x),x,0,5);
293 1 + x^2/4 + x^4/64$
295 taylor(bessel_i(1,x),x,0,5);
296 x/2 + x^3/16 + x^5/384$
298 taylor(bessel_i(m,x) / x^m,x,0,0);
299 1/(2^m * gamma(m + 1))$
301 /*--- gamma, polylog,  and zeta --- */
303 taylor(gamma(x),x,1,1);
304 1 - %gamma * (x - 1)$
306 taylor(log(gamma(x)),x,0,5);
307 (-log(x))+ -%gamma*x+(%pi^2*x^2)/12-(zeta(3)*x^3)/3+(%pi^4*x^4)/360-(zeta(5)*x^5)/5$
310  * The expansions for gamma_incomplete and gamma_incomplete_lower only works if
311  * gamma_expand is false.
312  */
313 taylor(gamma_incomplete(1/2,z), z, 0, 9);
314 sqrt(%pi)-2*sqrt(z)+(2*z^(3/2))/3-z^(5/2)/5+z^(7/2)/21-z^(9/2)/108
315                 +z^(11/2)/660-z^(13/2)/4680+z^(15/2)/37800-z^(17/2)/342720$
317 taylor(gamma_incomplete_lower(1/2,z), z, 0, 9);
318 2*sqrt(z)-(2*z^(3/2))/3+z^(5/2)/5-z^(7/2)/21+z^(9/2)/108-z^(11/2)/660
319                 +z^(13/2)/4680-z^(15/2)/37800+z^(17/2)/342720;
321 taylor(gamma_incomplete_lower(1/3,z), z, 0, 9);
322 3*z^(1/3)-(3*z^(4/3))/4+(3*z^(7/3))/14-z^(10/3)/20+z^(13/3)/104
323                 -z^(16/3)/640+z^(19/3)/4560-z^(22/3)/36960+z^(25/3)/336000;
325 taylor(zeta(x),x,0,1);
326  -1/2 - (log(%pi) + log(2)) * x / 2$
328 taylor(li[1](x),x,0,3);
329 x + ((x^2)/2) + ((x^3)/3)$
331 taylor(li[s](x),x,0,3);
332 x + ((x^2)/(2^s)) + ((x^3)/(3^s))$
334 errcatch(taylor(li[2](x), x, 1, 9));
337 /*--- expintegral ---*/
338 taylor(expintegral_si(z), z, 0, 9);
339 z-z^3/18+z^5/600-z^7/35280+z^9/3265920;
341 /*--- deftaylor ---*/
342 (deftaylor(unk(x), sum(x^k / k,k,1,inf)),0);
345 taylor(unk(s),s,0,3);
346 s + s^2/2 + s^3/3$
348 /*--- SF bug 1641487 ---*/
350 taylor(unk(x),x,1,2);
351 unk(1)+(at('diff(unk(x),x,1),x=1))*(x-1)$
353 (deftaylor(nam[n](x), sum(x^k/(k + n),k,0,inf)),0);
356 taylor(nam[5](x),x,0,4);
357 1/5+x/6+x^2/7+x^3/8+x^4/9$
359 taylor(nam[q](x),x,0,4);
360 1/q+x/(q+1)+x^2/(q+2)+x^3/(q+3)+x^4/(q+4)$
362 /*---other ---*/
364 (assume(t > 0),0);
367 taylor('integrate(exp(-x^4),x,0,t),t,0,9) - integrate(taylor(exp(-x^4),x,0,8),x,0,t);
370 (forget(t > 0),0);
373 is(equal(taylor(taylor(sqrt(1+x+y),[x,0,5]),[y,0,7]), taylor(sqrt(1+x+y),[y,0,7], [x,0,5])));
374 true$
376 taylor(sum(x^k/gamma(k + 1/2),k,0,inf),x,0,3);
377 sqrt(%pi)/%pi+(2*sqrt(%pi)*x)/%pi+(4*sqrt(%pi)*x^2)/(3*%pi)+(8*sqrt(%pi)*x^3)/(15*%pi)$
379 taylor(sum(x^k * a[k],k,0,inf),x,0,1);
380 a[0] + a[1] * x$
382 taylor(a^x+b^x,[a,b],0,1);
383 b^x + a^x$
385 /* SF bug 902757  */
387 taylor(x,[x],inf,2);
390 /* SF bug 836780 */
392 taylor(acosh(x),x,1,1) - taylor(logarc(acosh(x)),x,1,1);
395 /*  SF bug 774065 */
397 taylor(x/log(x),x,inf,1);
398 x/log(x)$
400 /* SF bug 788933 */
403 taylor(x^y,x,0,1,y,0,1);
404 1 + log(x) * y$ 
406 taylor(taylor(x^y,y,0,1),x,0,1);
407 1 + log(x) * y$ */
410 /* SF bug  701628*/
412 taylor(x^2,x,minf,2);
413 x^2$
415 /* SF bug 593449 -- fixed by gcd : spmod */
417 taylor(log(sqrt(y^e*y+1)+y),e,0,2);
418 log(sqrt(y+1)+y)+(sqrt(y+1)*log(y)*y*e)/(2*y^2+(2*sqrt(y+1)+2)*y+2*sqrt(y+1))+(((sqrt(y+1)+1)*log(y)^2*y^4+(2*sqrt(y+1)+5)*log(y)^2*y^3+(2*sqrt(y+1)+4)*log(y)^2*y^2+2*sqrt(y+1)*log(y)^2*y)*e^2)/(8*y^5+(24*sqrt(y+1)+40)*y^4+(56*sqrt(y+1)+80)*y^3+(48*sqrt(y+1)+72)*y^2+(24*sqrt(y+1)+24)*y+8*sqrt(y+1))$
420 /* SF #817516 */
422 taylor((log(n)-n)/(log(n)-1),n,inf,1);
423 ((-1/log(n)))*n+(1+1/log(n))$
425 /* SF #939022 */
427 (m: matrix([taylor(1+a*x,x,0,1),0], [0,taylor(1+d*x,x,0,1)]),0);
430 (mi : m^^-1,0);
433 is(mi[1,1] = taylor(1/(1+a*x),x,0,1));
434 true$
436 is(mi[2,2] = taylor(1/(1+d*x),x,0,1));
437 true$
439 taylorp(mi[1,1]);
440 true$
442 taylorp(mi[2,2]);
443 true$
445 /* SF #1203443 */
447 taylor(1/x^2,x,minf,2);
448 1/x^2$
450 taylor(1/x^2,x,-inf,2);
451 1/x^2$
453 /* SF #910008 */
455 taylor(exp(x),x,0,2) - taylor(exp(x),x,0,0);
458 /* SF #709138 */
460 errcatch(taylor(x^3,x,0,2)/taylor(x^3,x,0,2));
463 (remvalue(tp),0);
466 /* [ 696818 ] Taylor internal error (rat problem?) */
467 taylor(log(sqrt(e*%e^x+1)+e),x,0,2),gcd:subres;
468 log(sqrt(e+1)+e)+sqrt(e+1)*e*x/(2*e^2+(2*sqrt(e+1)+2)*e+2*sqrt(e+1))
469                       +((sqrt(e+1)+1)*e^5+(3*sqrt(e+1)+6)*e^4
470                                          +(4*sqrt(e+1)+9)*e^3
471                                          +(4*sqrt(e+1)+4)*e^2+2*sqrt(e+1)*e)
472                        *x^2
473                        /(8*e^6+(24*sqrt(e+1)+48)*e^5+(80*sqrt(e+1)+120)*e^4
474                               +(104*sqrt(e+1)+152)*e^3+(72*sqrt(e+1)+96)*e^2
475                               +(32*sqrt(e+1)+24)*e+8*sqrt(e+1))$
477 taylor(exp(x), x, minf, 1);
478 exp(x)$
480 /* [ 1682435 ] taylor of rat(factor(1+x)) gives internal error */
481 taylor(rat(factor(1+x)),x,0,1);
482 1+x$
484 /* [ 807275 ] Taylor Illegal log kernel: log(cos(th)) @ %pi/2 */
485 taylor(log(cos(th)),th,%pi/2,2);
486 log((2*th-%pi)/2)+log(-1)-(th-%pi/2)^2/6$
488 /* SF bug  752332 */
489 (assume(a < 1),0);
491 taylor(atan(1/sqrt(1-a))/sqrt(1-a),a,1,2);
492 ((%i * %pi)/(2 * sqrt(a - 1))) - 1 - ((a - 1)/3) - (((a - 1)^2)/5)$ 
494 /* [ 1358239 ] taylor(exp(sqrt(log(x)*log(log(x))))@inf > stack overflow */
495 taylor(exp(sqrt(log(x)*log(log(x)))), x, inf, 2);
496 exp(sqrt(log(x))*sqrt(log(log(x))));
498 /* [ 2010843 ] diff of Taylor poly */
499 diff(taylor(sqrt((f - z^2/4/f)^2 + (d-z)^2) + sqrt((f - z^2/4/f)^2 + (zp-z)^2),
500             [z, zp, d], 0, 4),
501      z);
502 (z-zp-d)/f+((3*zp+3*d)*z^2+(-5*zp^2-5*d^2)*z+2*zp^3+2*d^3)/(4*f^3);
504 /* [ 734851 ] pade interfered with by taylor */
505 ts:taylor(sin(x),x,0,3);
506 x-x^3/6;
508 taylor(x,x,0,3);
511 pade(ts,2,2);
512 [6*x/(x^2+6)];
514 /* [ 1729430 ] numberp of taylor poly */
515 numberp(taylor(x^2, x, 0, 1));
516 false;
518 /* Bug ID: 2907952 - diff of a taylor series
519  */
520 diff(taylor(1/(1-x^2),x,0,10),x);
521 2*x+4*x^3+6*x^5+8*x^7+10*x^9;
523 /* Bug ID: 660876 - taylor and keepfloat ERR
524  */
525 (old:keepfloat, keepfloat:true);
526 true;
527 (qq: taylor((x+1.0)^2,x,1,3),done);
528 done;
529 keepfloat:false;
530 false;
531 qq^2;
532 16+32*(x-1)+24*(x-1)^2+8*(x-1)^3;
533 (keepfloat:old,done);
534 done;
536 /* Bug ID: 2985866 - derivatives of functions of taylor polys
537  */
538 diff('(f(taylor(x,x,0,6))),x);
539 'diff(f(taylor(x,x,0,6)),x,1);
540 ev(%,nouns);
541 'diff(f(+x),x,1);
543 /* pade examples nicked from reference manual */
545 (taylor (1 + x + x^2 + x^3, x, 0, 3),
546  pade (%%, 1, 1));
547 [-1/(x - 1)];
549 (kill (t), 
550  t: taylor(-(83787*x^10 - 45552*x^9 - 187296*x^8
551                    + 387072*x^7 + 86016*x^6 - 1507328*x^5
552                    + 1966080*x^4 + 4194304*x^3 - 25165824*x^2
553                    + 67108864*x - 134217728)
554        /134217728, x, 0, 10),
555  pade (t, 4, 4));
558 (foo : [- (520256329*x^5  - 96719020632*x^4  - 489651410240*x^3
559  - 1619100813312*x^2  - 2176885157888*x - 2386516803584)
560 /(47041365435*x^5  + 381702613848*x^4  + 1360678489152*x^3
561  + 2856700692480*x^2  + 3370143559680*x + 2386516803584)],
562  pade (t, 5, 5), is (equal (%%, foo)));
563 true;
565 /* Bug #2116: lambda form for taylor_simplifier */
567 (kill (f), f : lambda ([s], s), 0);
570 /* This worked... */
571 block ([taylor_simplifier : 'f],
572   taylor (sin (x), x, 0, 4));
573 x - x^3 / 3!$
575 /* ...but this failed since the value of taylor_simplifier
576  * was previously required to be a symbol.
577  */
578 block ([taylor_simplifier : f],
579   taylor (sin (x), x, 0, 4));
580 x - x^3 / 3!$
582 (kill (f), 0);
585 /* mailing list 2020-04-30: "Bug in taylor?" */
587 kill (f, t, u);
588 done;
590 taylor( 'diff(f(t,u),t), t, 0, 1);
591 'at('diff(f(t,u),t,1),t = 0)+('at('diff(f(t,u),t,2),t = 0))*t;
593 taylor( 'diff(f(t,u),t), u, 0, 1);
594 'diff(f(t,0),t,1)+('at('diff(f(t,u),t,1,u,1),u = 0))*u;
596 taylor( 'diff(f(t,u),u), t, 0, 1);
597 'diff(f(0,u),u,1)+('at('diff(f(t,u),t,1,u,1),t = 0))*t;
599 taylor( 'diff(f(t,u),u), u, 0, 1);
600 'at('diff(f(t,u),u,1),u = 0)+('at('diff(f(t,u),u,2),u = 0))*u;
602 /*--- signum ---*/
603 taylor(signum(x), x, -2, 5)$
606 /* Bug #2167: taylor(rat(...)...) gives incorrect results
607  * Bug #3068: taylor of CRE fails
608  */
610 /* This used to return the bogus series 1/x */
611 taylor(rat(1/x),x,1,2);
612 1-(x-1)+(x-1)^2;
614 /* This used to signal an error about an unfamiliar singularity */
615 taylor(rat(1/z+1),z,0,1);
616 1/z+1;
618 /* This used to signal an error about an unfamiliar singularity */
619 taylor(rat(z+1),z,inf,0);
620 z+1;
622 /* This used to signal an error about an unfamiliar singularity */
623 taylor(rat(x/(1+x)),x,inf,2);
624 1-1/x+1/x^2;
626 /* Bug #1848: taytorat leaks internal gensyms from multivar expansions
627  * Bug #2446: horner of multivariate taylor gives junk
628  */
630 taytorat (taylor (x, [x], 0, 1));
631 ''(rat (x));
633 taytorat (taylor (taylor (x + y, [x], 0, 1), [y], 0, 1));
634 ''(rat (x + y));
636 horner (taylor (x, [x], 0, 1));
639 /* Bug #608: taylor(x^a,[x],0,1) unsimplified
640  * Bug #545: multivar taylor gives 1^2
641  */
643 /* This used to yield "+1^2*x^a" */
644 string(taylor(x^a,[x],0,1));
645 "+x^a";
647 /* This used to yield "+1^4*(x*y)^a" */
648 string(taylor((x*y)^a,[x,y],0,1));
649 "+(x*y)^a";