1%       Examples for the algorithmic calculation of formal
2%       Puiseux, Laurent and power series,
3%
4%       Wolfram Koepf, Freie Universitaet Berlin, Germany
5%       (taken from the original paper and adapted to REDUCE
6%        form by Winfried Neun, ZIB Berlin)
7
8% Formal Laurent series
9
10fps(E^x,x);
11
12
13              k
14             x
15infsum(--------------,k,0,infinity)
16        factorial(k)
17
18
19fps(E^x/(x^3),x);
20
21
22               k
23              x
24infsum(-----------------,k,0,infinity)
25                      3
26        factorial(k)*x
27
28
29fps(x * e^(x^4),x);
30
31
32            4*k
33           x   *x
34infsum(--------------,k,0,infinity)
35        factorial(k)
36
37
38fps(sin (x + y),x);
39
40
41         2*k       k
42        x   *( - 1) *cos(y)*x
43infsum(-----------------------,k,0,infinity)
44         factorial(2*k + 1)
45
46            2*k       k
47           x   *( - 1) *sin(y)
48 + infsum(---------------------,k,0,infinity)
49             factorial(2*k)
50
51
52simplede (sin x,x);
53
54
55df(y,x,2) + y
56 %find a DE for sin
57
58simplede (sin (x)^2,x,w);
59
60
61df(w,x,3) + 4*df(w,x)
62 % DE in w and x
63
64fps(asin x,x);
65
66
67            2*k
68           x   *factorial(2*k)*x
69infsum(------------------------------,k,0,infinity)
70         2*k             2
71        2   *factorial(k) *(2*k + 1)
72
73
74fps((asin x)^2,x);
75
76
77         2*k  2*k             2  2
78        x   *2   *factorial(k) *x
79infsum(----------------------------,k,0,infinity)
80        factorial(2*k + 1)*(k + 1)
81
82
83fps(e^(asin x),x);
84
85
86         2*k  k         2
87        x   *2 *prod(2*j  - 2*j + 1,j,1,k)*x
88infsum(--------------------------------------,k,0,infinity)
89                 factorial(2*k + 1)
90
91            2*k         2
92           x   *prod(4*j  - 8*j + 5,j,1,k)
93 + infsum(---------------------------------,k,0,infinity)
94                   factorial(2*k)
95
96
97fps(e^(asinh x),x);
98
99
100            2*k       k
101         - x   *( - 1) *factorial(2*k)
102infsum(--------------------------------,k,0,infinity) + x
103           k             2
104          4 *factorial(k) *(2*k - 1)
105
106
107fps((x + sqrt(1+x^2))^A,x);
108
109
110         2*k       k  2*k              - a                 a
111        x   *( - 1) *2   *pochhammer(------,k)*pochhammer(---,k)
112                                       2                   2
113infsum(----------------------------------------------------------,k,0,infinity)
114                             factorial(2*k)
115
116           2*k       k  2*k              - a + 1                 a + 1
117          x   *( - 1) *2   *pochhammer(----------,k)*pochhammer(-------,k)*a*x
118                                           2                       2
119+ infsum(----------------------------------------------------------------------,
120                                   factorial(2*k + 1)
121
122         k,0,infinity)
123
124
125fps(e^(x^2)*erf x,x);
126
127
128           2*k           2*k
129        2*x   *sqrt(pi)*2   *factorial(k)*x
130infsum(-------------------------------------,k,0,infinity)
131               factorial(2*k + 1)*pi
132
133
134fps(e^x - 2 e^(-x/2) * cos(sqrt(3) * x/2 -pi/3),x);
135
136
137           3*k  2
138        9*x   *x *(k + 1)
139infsum(--------------------,k,0,infinity)
140        factorial(3*k + 3)
141
142
143% fps(int(e^(-a^2*t^2) * cos(2*x*t),t,0,infinity),x)  % not yet
144
145% fps(4/x * int(e^(t^2)*erf(t),t,0,sqrt(x)/2),x);
146
147fps(sin x * e^x,x);
148
149
150         k  k/2      k*pi
151        x *2   *sin(------)
152                      4
153infsum(---------------------,k,0,infinity)
154           factorial(k)
155
156
157fps(cos x * e^(2*x),x);
158
159
160         k  k/2           1
161        x *5   *cos(atan(---)*k)
162                          2
163infsum(--------------------------,k,0,infinity)
164              factorial(k)
165
166
167fps(1/(x-x^3),x);
168
169
170          k       k    k
171         x *( - 1)  - x
172 infsum(-----------------,k,0,infinity)*x + 1
173                    k
174            2*( - 1)
175----------------------------------------------
176                      x
177
178
179fps(1/(x^2 + 3 x + 2),x);
180
181
182           k  k    k
183        2*x *2  - x
184infsum(--------------,k,0,infinity)
185                k  k
186        2*( - 1) *2
187
188
189fps(x/(1-x-x^2),x);
190
191
192          x
193fps(--------------,x,0)
194                2
195     (1 - x) - x
196
197
198% Logarithmic singularities and Puisieux series
199
200fps(sin sqrt x,x);
201
202
203         (2*k + 1)/2       k
204        x           *( - 1)
205infsum(----------------------,k,0,infinity)
206         factorial(2*k + 1)
207
208
209fps(((1 + sqrt x)/x)^(1/3),x);
210
211
212         (6*k + 1)/6             2
213        x           *pochhammer(---,2*k)
214                                 3
215infsum(----------------------------------,k,0,infinity)
216              3*factorial(2*k + 1)
217
218            k              - 1
219           x *pochhammer(------,2*k)
220                           3
221 + infsum(---------------------------,k,0,infinity)
222               1/3
223              x   *factorial(2*k)
224
225
226fps(asech x,x);
227
228
229
230% some more (Wolfram Koepf, priv. comm.)
231
232fps((1+x)^alpha,x);
233
234
235         k       k
236        x *( - 1) *pochhammer( - alpha,k)
237infsum(-----------------------------------,k,0,infinity)
238                  factorial(k)
239
240
241fps((1+sqrt(1+x))^beta,x);
242
243
244          k       k  beta
245         x *( - 1) *2    *pochhammer( - beta,2*k)
246infsum(---------------------------------------------,k,0,infinity)
247         2*k
248        2   *factorial(k)*pochhammer( - beta + 1,k)
249
250
251fps(sin(x)^2+cos(x)^2,x);
252
253
2541
255
256
257fps(sin(x)^2*cos(x)^2,x);
258
259
260            2*k       k  4*k  2
261           x   *( - 1) *2   *x
262infsum(----------------------------,k,0,infinity)
263        factorial(2*k + 1)*(k + 1)
264
265
266fps(sin(x)*cos(x^2),x);
267
268
269                2
270fps(sin(x)*cos(x ),x,0)
271
272
273fps((x-1)^(-1),x);
274
275
276           k
277infsum( - x ,k,0,infinity)
278
279
280fps(atan(x+y),x);
281
282
283fps(atan(x + y),x,0)
284
285
286fps((1-x^5)^6,x);
287
288
289 30      25       20       15       10      5
290x   - 6*x   + 15*x   - 20*x   + 15*x   - 6*x  + 1
291
292
293fps(asec x,x);
294
295
296
297fps(besseli(0,x),x);
298
299
300                2*k
301               x
302infsum(--------------------,k,0,infinity)
303         2*k             2
304        2   *factorial(k)
305
306fps(besseli(1,x),x);
307
308
309                        2*k
310                       x   *x
311infsum(--------------------------------------,k,0,infinity)
312           2*k
313        2*2   *factorial(k + 1)*factorial(k)
314
315
316fps(exp(x^(1/3)),x);
317
318
319            (3*k + 1)/3
320           x
321infsum(--------------------,k,0,infinity)
322        factorial(3*k + 1)
323
324                  k
325                 x
326 + infsum(----------------,k,0,infinity)
327           factorial(3*k)
328
329              (3*k + 2)/3
330           3*x           *(k + 1)
331 + infsum(------------------------,k,0,infinity)
332             factorial(3*k + 3)
333
334fps(log(1-x),x);
335
336
337            k
338         - x *x
339infsum(---------,k,0,infinity)
340         k + 1
341
342fps(exp x*sinh x,x);
343
344
345             k  k
346            x *2 *x
347infsum(------------------,k,0,infinity)
348        factorial(k + 1)
349
350fps(atan x,x);
351
352
353         2*k       k
354        x   *( - 1) *x
355infsum(----------------,k,0,infinity)
356           2*k + 1
357
358fps(sin x+sinh x,x);
359
360
361                4*k
362             2*x   *x
363infsum(--------------------,k,0,infinity)
364        factorial(4*k + 1)
365
366fps(sin x*sinh x,x);
367
368
369             4*k       k  2*k  2
370            x   *( - 1) *2   *x
371infsum(------------------------------,k,0,infinity)
372        factorial(4*k + 1)*(2*k + 1)
373
374fps(int(erf(x),x),x);
375
376
377            2*k                k
378         - x   *sqrt(pi)*( - 1)
379infsum(---------------------------,k,0,infinity)
380        factorial(k)*pi*(2*k - 1)
381
382fps(sqrt(2-x),x);
383
384
385            k
386         - x *sqrt(2)*factorial(2*k)
387infsum(------------------------------,k,0,infinity)
388          k             2
389         8 *factorial(k) *(2*k - 1)
390
391fps(sqrt(1+x)+sqrt(1-x),x);
392
393
394                 2*k
395            - 2*x   *factorial(4*k)
396infsum(--------------------------------,k,0,infinity)
397         2*k               2
398        4   *factorial(2*k) *(4*k - 1)
399
400fps(exp(a+b*x)*exp(c+d*x),x);
401
402
403         k  a + c        k
404        x *e     *(b + d)
405infsum(--------------------,k,0,infinity)
406           factorial(k)
407
408fps(1/cos(asin x),x);
409
410
411         2*k
412        x   *factorial(2*k)
413infsum(---------------------,k,0,infinity)
414         2*k             2
415        2   *factorial(k)
416
417fps(sqrt(1-x^2)+x*asin x,x);
418
419
420                2*k
421               x   *factorial(2*k)
422infsum(-----------------------------------,k,0,infinity)
423         k             2     2
424        4 *factorial(k) *(4*k  - 4*k + 1)
425
426fps(sqrt(1-sqrt(x)),x);
427
428
429                 (2*k + 1)/2
430              - x           *factorial(4*k)
431infsum(------------------------------------------,k,0,infinity)
432           4*k
433        2*2   *factorial(2*k + 1)*factorial(2*k)
434
435                    k
436                 - x *factorial(4*k)
437 + infsum(--------------------------------,k,0,infinity)
438            2*k               2
439           4   *factorial(2*k) *(4*k - 1)
440
441fps(cos(n*acos x),x);
442
443
444         2*k  2*k      n*pi                - n                 n
445        x   *2   *cos(------)*pochhammer(------,k)*pochhammer(---,k)
446                        2                  2                   2
447infsum(--------------------------------------------------------------,k,0,
448                               factorial(2*k)
449
450       infinity) + infsum(
451
452     2*k  2*k              - n + 1                 n + 1          n*pi
453    x   *2   *pochhammer(----------,k)*pochhammer(-------,k)*sin(------)*n*x
454                             2                       2             2
455   --------------------------------------------------------------------------,k,
456                               factorial(2*k + 1)
457
458   0,infinity)
459
460fps(cos x+I*sin x,x);
461
462
463            k  k
464           x *i
465infsum(--------------,k,0,infinity)
466        factorial(k)
467
468fps(cos(3*asinh x),x);
469
470
471         2*k       k         2
472        x   *( - 1) *prod(4*j  - 8*j + 13,j,1,k)
473infsum(------------------------------------------,k,0,infinity)
474                     factorial(2*k)
475
476fps(cos(n*asinh x),x);
477
478
479         2*k       k  2*k              - i*n                 i*n
480        x   *( - 1) *2   *pochhammer(--------,k)*pochhammer(-----,k)
481                                        2                     2
482infsum(--------------------------------------------------------------,k,0,
483                               factorial(2*k)
484
485       infinity)
486
487fps(sin(n*log(x+sqrt(1+x^2))),x);
488
489
490         2*k       k  2*k              - i*n + 1                 i*n + 1
491infsum((x   *( - 1) *2   *pochhammer(------------,k)*pochhammer(---------,k)*n*x
492                                          2                         2
493
494        )/factorial(2*k + 1),k,0,infinity)
495
496fps(sqrt(1+x^2)*asinh x-x,x);
497
498
499           2*k       k  2*k                                3
500        2*x   *( - 1) *2   *factorial(k + 1)*factorial(k)*x
501infsum(------------------------------------------------------,k,0,infinity)
502                         factorial(2*k + 3)
503
504
505fps(int(erf(x)/x,x),x);
506
507
508              2*k                k
509           2*x   *sqrt(pi)*( - 1) *x
510infsum(----------------------------------,k,0,infinity)
511                            2
512        factorial(k)*pi*(4*k  + 4*k + 1)
513
514                erf(x)
515 + sub(x=0,int(--------,x))
516                  x
517
518fps(asin(x)^2/x^4,x);
519
520
521            2*k  2*k             2
522           x   *2   *factorial(k)
523infsum(-------------------------------,k,0,infinity)
524                            2
525        factorial(2*k + 1)*x *(k + 1)
526
527
528
529% we had problems here:
530
531fps(cos(asin x),x);
532
533
534              2*k
535           - x   *factorial(2*k)
536infsum(----------------------------,k,0,infinity)
537         k             2
538        4 *factorial(k) *(2*k - 1)
539
540fps(sinh(log x),x);
541
542
543fps(sinh(log(x)),x,0)
544
545fps(atan(cot x),x);
546
547
548Could not find the limit of: atan(cot(x)),x,0
549
550
551% we can cure this one by defining the limit:
552
553let limit(atan(cot ~x),x,0) => pi/2;
554
555
556fps(atan(cot x),x);
557
558
559 pi - 2*x
560----------
561    2
562
563
564fps(exp(nnn*x)*cos(mmm*x),x);
565
566
567         k              2                                          2
568infsum((x *((impart(mmm)  + 2*impart(mmm)*repart(nnn) + impart(nnn)
569
570                                                       2              2    k
571              - 2*impart(nnn)*repart(mmm) + repart(mmm)  + repart(nnn) )**---
572                                                                           2
573
574                       impart(nnn) - repart(mmm)                    2
575            *cos(atan(---------------------------)*k) + (impart(mmm)
576                       impart(mmm) + repart(nnn)
577
578                                                         2
579                - 2*impart(mmm)*repart(nnn) + impart(nnn)
580
581                                                         2              2    k
582                + 2*impart(nnn)*repart(mmm) + repart(mmm)  + repart(nnn) )**---
583                                                                             2
584
585                       impart(nnn) + repart(mmm)                    2
586            *cos(atan(---------------------------)*k) - (impart(mmm)
587                       impart(mmm) - repart(nnn)
588
589                                                         2
590                - 2*impart(mmm)*repart(nnn) + impart(nnn)
591
592                                                         2              2    k
593                + 2*impart(nnn)*repart(mmm) + repart(mmm)  + repart(nnn) )**---
594                                                                             2
595
596                       impart(nnn) + repart(mmm)                      2
597            *sin(atan(---------------------------)*k)*i + (impart(mmm)
598                       impart(mmm) - repart(nnn)
599
600                                                         2
601                + 2*impart(mmm)*repart(nnn) + impart(nnn)
602
603                                                         2              2    k
604                - 2*impart(nnn)*repart(mmm) + repart(mmm)  + repart(nnn) )**---
605                                                                             2
606
607                       impart(nnn) - repart(mmm)
608            *sin(atan(---------------------------)*k)*i))/(2*factorial(k)),k,0,
609                       impart(mmm) + repart(nnn)
610
611       infinity)
612
613fps(sqrt(2-x^2),x);
614
615
616                   2*k
617              - 2*x   *factorial(2*k)
618infsum(------------------------------------,k,0,infinity)
619                 k             2
620        sqrt(2)*8 *factorial(k) *(2*k - 1)
621
622fps(ci x,x);
623
624
625                   2*k       k
626                  x   *( - 1) *infinity*x
627ci(0) + infsum(------------------------------,k,0,infinity)
628                factorial(2*k + 1)*(2*k + 1)
629
630fps(log(1-2*x*y+x^2),x);
631
632
633                     2
634fps(log(1 - 2*x*y + x ),x,0)
635
636
637FPS(sin x,x,pi);
638
639
640                   2*k       k
641        ( - pi + x)   *( - 1) *( - pi + x)
642infsum(------------------------------------,k,0,infinity)
643                factorial(2*k + 1)
644
645
646% detect logarithmic singularity and bail out
647
648fps(dilog x,x,0);
649
650
651fps(dilog(x),x,0)
652
653fps(dilog x,x,1);
654
655
656fps(dilog(x),x,1)
657
658fps(dilog x,x,2);
659
660
661fps(dilog(x),x,2)
662
663
664% This one takes ages :
665
666%fps(acos(cos(x)),x);
667
668fps_search_depth := 7;
669
670
671fps_search_depth := 7
672 % does not find aa DE with the default
673fps(sin(x^(1/3)),x);
674
675
676            2*k       k    k                         k
677infsum(( - x   *( - 1) *108 *factorial(k)*x)/(6*46656 *factorial(3*k + 1)
678
679                                          7                 5
680          *factorial(2*k + 1)*pochhammer(---,k)*pochhammer(---,k)),k,0,infinity)
681                                          6                 6
682
683            (6*k + 2)/3       k  2*k  3*k                              k
684 + infsum((x           *( - 1) *2   *3   *factorial(k + 1)*x)/(20*46656
685
686                                                                11
687             *factorial(3*k + 3)*factorial(2*k + 1)*pochhammer(----,k)
688                                                                6
689
690                          7
691             *pochhammer(---,k)),k,0,infinity) + infsum(
692                          6
693
694                   (6*k + 1)/3       k  2*k  3*k
695                  x           *( - 1) *2   *3   *factorial(k)
696   --------------------------------------------------------------------------,k,
697         k                                           7                 5
698    46656 *factorial(3*k)*factorial(2*k)*pochhammer(---,k)*pochhammer(---,k)
699                                                     6                 6
700
701   0,infinity)
702
703
704end;
705
706Tested on x86_64-pc-windows CSL
707Time (counter 1): 1793 ms  plus GC time: 157 ms
708
709End of Lisp run after 1.79+0.22 seconds
710real 2.24
711user 0.00
712sys 0.06
713