Mon Feb 16 22:57:02 MST 1998 aquarius% mupad *----* MuPAD 1.4.0 -- Multi Processing Algebra Data Tool /| /| *----* | Copyright (c) 1997 - 98 by SciFace Software GmbH | *--|-* All rights reserved. |/ |/ *----* Licensed to: Michael Wester >> # ----------[ M u P A D ]---------- # >> # ---------- Initialization ---------- # >> TEXTWIDTH:= 80: >> read("../../Time.mupad"): >> # ---------- Series ---------- # >> # Taylor series---this first example comes from special relativity &> => 1 + 1/2 (v/c)^2 + 3/8 (v/c)^4 + 5/16 (v/c)^6 + O((v/c)^8) # >> 1/sqrt(1 - (v/c)^2); 1 ------------- / 2 \1/2 | v | | 1 - -- | | 2 | \ c / Time: 280 msec Type: "_power" >> series(%, v = 0, 7); 2 4 6 v 3 v 5 v 7 1 + ---- + ---- + ----- + O(v ) 2 4 6 2 c 8 c 16 c Time: 3580 msec Type: Puiseux >> 1/%^2; 2 v 7 1 - -- + O(v ) 2 c Time: 140 msec Type: Puiseux >> # Note: sin(x) = x - x^3/6 + x^5/120 - x^7/5040 + O(x^9) &> cos(x) = 1 - x^2/2 + x^4/24 - x^6/720 + O(x^8) &> tan(x) = x + x^3/3 + 2/15 x^5 + 17/315 x^7 + O(x^9) # >> tsin:= series(sin(x), x = 0, 8); 3 5 7 x x x 8 x - -- + --- - ---- + O(x ) 6 120 5040 Time: 320 msec Type: Puiseux >> tcos:= series(cos(x), x = 0, 7); 2 4 6 x x x 7 1 - -- + -- - --- + O(x ) 2 24 720 Time: 300 msec Type: Puiseux >> # Note that additional terms will be computed as needed # >> tsin/tcos; 3 5 7 x 2 x 17 x 8 x + -- + ---- + ----- + O(x ) 3 15 315 Time: 310 msec Type: Puiseux >> series(tan(x), x = 0, 8); 3 5 7 x 2 x 17 x 8 x + -- + ---- + ----- + O(x ) 3 15 315 Time: 350 msec Type: Puiseux >> tsin:= NIL: >> tcos:= NIL: >> # => -x^2/6 - x^4/180 - x^6/2835 - O(x^8) # >> series(ln(sin(x)/x), x = 0, 8); 2 4 6 x x x 7 - -- - --- - ---- + O(x ) 6 180 2835 Time: 1450 msec Type: Puiseux >> series(sin(x)/x, x = 0, 8); 2 4 6 x x x 7 1 - -- + --- - ---- + O(x ) 6 120 5040 Time: 170 msec Type: Puiseux >> series(ln(expr(%)), x = 0, 8); 2 4 6 x x x 8 - -- - --- - ---- + O(x ) 6 180 2835 Time: 360 msec Type: Puiseux >> # => [a f'(a d) + g(b d) + integrate(h(c y), y = 0..d)] &> + [a^2 f''(a d) + b g'(b d) + h(c d)] (x - d) # >> q:= diff(f(a*x), x) + g(b*x) + int(h(c*y), y = 0..x); g(b x) + a D(f)(a x) + int(h(c y), y = 0..x) Time: 19010 msec Type: "_plus" >> series(q, x = d, 1); (g(b d) + a D(f)(a d)) + O(- d + x ) Time: 650 msec Type: Puiseux >> series(q, x = d, 2); series(g(b x) + a D(f)(a x) + int(h(c y), y = 0..x), x = d, 2) Time: 4680 msec Type: "series" >> q:= NIL: >> # Taylor series of nonscalar objects (noncommutative multiplication) &> => (B A - A B) t^2/2 + O(t^3) [Stanly Steinberg] # >> E^((A + B)*t) - E^(A*t) * E^(B*t); t (A + B) A t B t exp(1) - exp(1) exp(1) Time: 440 msec Type: "_plus" >> simplify(%); exp(t (A + B)) - exp(A t) exp(B t) Time: 1170 msec Type: "_plus" >> series(E^((A + B)*t) - E^(A*t) * E^(B*t), t = 0, 4); 4 O(t ) Time: 560 msec Type: Puiseux >> # Laurent series: &> => sum( Bernoulli[k]/k! x^(k - 2), k = 1..infinity ) &> = 1/x^2 - 1/(2 x) + 1/12 - x^2/720 + x^4/30240 + O(x^6) &> [Levinson and Redheffer, p. 173] # >> series(1/(x*(exp(x) - 1)), x = 0, 9); 2 4 1 1 x x 6 -- - --- + 1/12 - --- + ----- + O(x ) 2 2 x 720 30240 x Time: 370 msec Type: Puiseux >> # Puiseux series (terms with fractional degree): &> => 1/sqrt(x - 3/2 pi) + (x - 3/2 pi)^(3/2) / 12 + O([x - 3/2 pi]^(7/2)) # >> series(sqrt(sec(x)), x = 3/2*PI, 4); / 3 PI \3/2 | x - ---- | 1 \ 2 / / / 3 PI \5/2 \ --------------- + --------------- + O| | x - ---- | | / 3 PI \1/2 12 \ \ 2 / / | x - ---- | \ 2 / Time: 580 msec Type: Puiseux >> # Generalized Taylor series => sum( [x log x]^k/k!, k = 0..infinity ) # >> series(x^x, x = 0, 4); 2 2 3 3 x ln(x) x ln(x) 4 1 + x ln(x) + --------- + --------- + O(x ) 2 6 Time: 460 msec Type: Puiseux >> # Compare the generalized Taylor series of two different formulations of a &> function => log(z) + log(cosh(w)) + tanh(w) z + O(z^2) # >> s1:= series(ln(sinh(z)) + ln(cosh(z + w)), z = 0, 1); Error: Illegal operand [_index]; in procedure 'funcattr(ln, "series")' >> s2:= series(ln(sinh(z) * cosh(z + w)), z = 0, 1); Error: Illegal operand [_index]; in procedure 'funcattr(ln, "series")' >> s1 - s2; s1 - s2 Time: 1060 msec Type: "_plus" >> s1:= series(ln(sinh(z)) + ln(cosh(z + w)), z = 0, 3); z sinh(w) 2 (ln(z) + ln(cosh(w))) + --------- + O(z ) cosh(w) Time: 550 msec Type: Puiseux >> s2:= series(ln(sinh(z) * cosh(z + w)), z = 0, 3); z sinh(w) 2 (ln(z) + ln(cosh(w))) + --------- + O(z ) cosh(w) Time: 470 msec Type: Puiseux >> s1 - s2; 0 Time: 230 msec Type: DOM_INT >> s1:= NIL: >> s2:= NIL: >> # Look at the generalized Taylor series around x = 1 &> => (x - 1)^a/e^b [1 - (a + 2 b) (x - 1) / 2 + O((x - 1)^2)] # >> ln(x)^a*exp(-b*x); a ln(x) exp(-b x) Time: 680 msec Type: "_mult" >> series(%, x = 1, 2); a series(ln(x) exp(-b x), x = 1, 2) Time: 5010 msec Type: "series" >> # Asymptotic expansions => sqrt(2) x + O(1/x) # >> series(sqrt(2*x^2 + 1), x = infinity, 1); Error: division by zero [Puiseux::_invert] >> series(sqrt(2*x^2 + 1), x = infinity, 3); 1/2 / 1 \ x 2 + O| - | \ x / Time: 700 msec Type: Puiseux >> # Wallis' product => 1/sqrt(pi n) + ... [Knopp, p. 385] # >> series(1/2^(2*n) * binomial(2*n, n), n = infinity, 3); Error: Division by zero [_power]; in procedure 'Series::unknown' >> # => 0!/x - 1!/x^2 + 2!/x^3 - 3!/x^4 + O(1/x^5) [Knopp, p. 544] # >> exp(x) * int(exp(-t)/t, t = x..infinity); / exp(-t) \ exp(x) int| -------, t = x..infinity | \ t / Time: 8910 msec Type: "_mult" >> series(%, x = infinity, 5); / / exp(-t) \ \ series| exp(x) int| -------, t = x..infinity |, x = infinity, 5 | \ \ t / / Time: 11100 msec Type: "series" >> # Multivariate Taylor series expansion => 1 - (x^2 + 2 x y + y^2)/2 + O(x^4) &> # >> series(series(cos(x + y), x = 0, 4), y = 0, 4); 2 4 3 2 2 (- 0.5 x + O(x ) + 1.0 ) + y (- 1.0 x + 0.1666666666 x ) + y (0.25 x - 0.5) 3 3 4 + y (0.1666666666 x - 0.02777777777 x ) + O(y ) Time: 720 msec Type: Puiseux >> # Power series (compute the general formula) # >> series(ln(sin(x)/x), x, infinity); / / sin(x) \ \ series| ln| ------ |, x, infinity | \ \ x / / Time: 490 msec Type: "series" >> series(exp(-x)*sin(x), x, infinity); / n n 1/2 / 3 n PI \ \ | x (2 ) sin| ------ | | | \ 4 / | sum| ------------------------, n = 0..infinity | \ fact(n) / Time: 4620 msec Type: "function" >> # Derive an explicit Taylor series solution of y as a function of x from the &> following implicit relation: &> y = x - 1 + (x - 1)^2/2 + 2/3 (x - 1)^3 + (x - 1)^4 + 17/10 (x - 1)^5 + ... # >> x = sin(y) + cos(y); x = cos(y) + sin(y) Time: 380 msec Type: "_equal" >> q:= revert(series(op(%, 2), y = 0)); 2 3 5 (y - 1) 2 (y - 1) 4 17 (y - 1) 6 (y - 1) + -------- + ---------- + (y - 1) + ----------- + O((y - 1) ) 2 3 10 Time: 620 msec Type: Puiseux >> y = subs(q, y = x); 2 3 5 (y - 1) 2 (y - 1) 4 17 (y - 1) 6 y = (y - 1) + -------- + ---------- + (y - 1) + ----------- + O((y - 1) ) 2 3 10 Time: 390 msec Type: "_equal" >> y = subsop(q, [5, 1] = x); 2 3 5 (x - 1) 2 (x - 1) 4 17 (x - 1) 6 y = (x - 1) + -------- + ---------- + (x - 1) + ----------- + O((x - 1) ) 2 3 10 Time: 380 msec Type: "_equal" >> q:= NIL: >> # Pade (rational function) approximation => (2 - x)/(2 + x) # >> Dom::Pade(exp(-x), x, 3); - x + 2 -------- x + 2 Time: 1180 msec Type: Dom::Pade >> # Fourier series of f(x) of period 2 p over the interval [-p, p] &> => - (2 p / pi) sum( (-1)^n sin(n pi x / p) / n, n = 1..infinity ) # >> x; x Time: 430 msec Type: DOM_IDENT >> # => p / 2 &> - (2 p / pi^2) sum( [1 - (-1)^n] cos(n pi x / p) / n^2, n = 1..infinity ) # >> abs(x); abs(x) Time: 360 msec Type: "abs" >> # ---------- Quit ---------- # >> quit real 77.25 user 75.69 sys 1.40