Mon Feb 22 22:09:22 MST 1999 aquarius% macsyma Starting Macsyma math engine with no window system... /opt/local/macsyma_422/macsyma.422 local This is Macsyma 422.0 for Sparc (Solaris 2.x) computers. Copyright (c) 1982 - 1998 Macsyma Inc. All rights reserved. Portions copyright (c) 1982 Massachusetts Institute of Technology. All rights reserved. Type "DESCRIBE(TRADE_SECRET);" to see important legal notices. Type "HELP();" for more information. /aquarius/data2/opt/local/macsyma_422/system/init.lsp being loaded. /aquarius/home/wester/macsyma-init.lsp being loaded. (c1) (c2) /* ----------[ M a c s y m a ]---------- */ /* ---------- Initialization ---------- */ symbol_display_case: lower_case$ Time= 0 msecs (c3) showtime: all$ Time= 0 msecs (c4) prederror: false$ Time= 0 msecs (c5) /* ---------- Number Theory ---------- */ /* Display the largest 6-digit prime and the smallest 7-digit prime => [999983, 1000003] */ [prime(78498), prime(78499)]; /aquarius/data2/opt/local/macsyma_422/library1/combin.so being loaded. Time= 760 msecs (d5) [999983, 1000003] (c6) /* Primitive root => 19 */ 191; Time= 0 msecs (d6) 191 (c7) /* (a + b)^p mod p => a^p + b^p for p prime and a, b in Z_p [Chris Hurlburt] See Thomas W. Hungerford, _Algebra_, Springer-Verlag, 1974, p. 121 for a more general simplification: (a +- b)^(p^n) => a^(p^n) +- b^(p^n) */ assume(equal(primep(p), true))$ Time= 40 msecs (c8) errcatch(mod((a + b)^p, p)); Improper value for MODULUS in call to MOD: p Time= 0 msecs (d8) [] (c9) forget(equal(primep(p), true))$ Time= 10 msecs (c10) /* Congruence equations. See Harold M. Stark, _An Introduction to Number Theory_, The MIT press, 1984. 9 x = 15 mod 21 => x = 4 mod 7 or {4, 11, 18} mod 21 [Stark, p. 68] */ modulus: 21$ Warning: MODULUS being set to 21, a non-prime. Time= 0 msecs (c11) errcatch(solve(9*x = 15, x)); Inverse of zero divisor? Time= 10 msecs (d11) [] (c12) /* 7 x = 22 mod 39 => x = 5 mod 13 or 31 mod 39 [Stark, p. 69] */ modulus: 39$ Warning: MODULUS being set to 39, a non-prime. Time= 0 msecs (c13) solve(7*x = 22, x); Time= 10 msecs (d13) [x = - 8] (c14) /* x^2 + x + 4 = 0 mod 8 => x = {3, 4} mod 8 [Stark, p. 97] */ modulus: 8$ Warning: MODULUS being set to 8, a non-prime. Time= 0 msecs (c15) errcatch(solve(x^2 + x + 4 = 0, x)); Inverse of zero divisor? Time= 10 msecs (d15) [] (c16) /* x^3 + 2 x^2 + 5 x + 6 = 0 mod 11 => x = 3 mod 11 [Stark, p. 97] */ modulus: 11$ Time= 0 msecs (c17) solve(x^3 + 2*x^2 + 5*x + 6 = 0, x); Time= 20 msecs (d17) [x = 3] (c18) multiplicities; Time= 0 msecs (d18) [3] (c19) /* {x = 7 mod 9, x = 13 mod 23, x = 1 mod 2} => x = 151 mod 414 [Stark, p. 76] */ /* {5 x + 4 y = 6 mod 7, 3 x - 2 y = 6 mod 7} => x = 1 mod 7, y = 2 mod 7 [Stark, p. 76] */ modulus: 7$ Time= 0 msecs (c20) solve([5*x + 4*y = 6, 3*x - 2*y = 6], [x, y]); Time= 20 msecs (d20) [[x = 1, y = 2]] (c21) /* 2 x + 3 y = 1 mod 5 => (x, y) = {(0, 2), (1, 3), (2, 4), (3, 0), (4, 1)} mod 5 */ modulus: 5$ Time= 0 msecs (c22) solve(2*x + 3*y = 1, [x, y]); Time= 10 msecs (d22) [[x = %r1 - 2, y = %r1]] (c23) /* 2 x + 3 y = 1 mod 6 => [Stark, p. 76] (x, y) = {(2, 1), (2, 3), (2, 5), (5, 1), (5, 3), (5, 5)} mod 6 */ modulus: 6$ Warning: MODULUS being set to 6, a non-prime. Time= 0 msecs (c24) errcatch(solve(2*x + 3*y = 1, [x, y])); Inverse of zero divisor? Time= 10 msecs (d24) [] (c25) modulus: false$ Time= 0 msecs (c26) /* Diophantine equations => x = 2, y = 5 (Wallis) [Stark, p. 147] */ declare([x, y], integer)$ Time= 10 msecs (c27) solve(x^4 + 9 = y^2, [x, y]); /aquarius/data2/opt/local/macsyma_422/library1/algsys.so being loaded. /aquarius/data2/opt/local/macsyma_422/library1/triangsy.so being loaded. Time= 230 msecs 2 1/4 2 1/4 (d27) [[y = %r3, x = (%r3 - 9) ], [y = %r3, x = - %i (%r3 - 9) ], 2 1/4 2 1/4 [y = %r3, x = - (%r3 - 9) ], [y = %r3, x = %i (%r3 - 9) ]] (c28) /* => x = 11, y = 5 (Fermat) [Stark, p. 147] */ solve(x^2 + 4 = y^3, [x, y]); Time= 10 msecs 3 2 (d28) [[x = %r4, y = root_of(y - %r4 - 4, y)]] (c29) /* => (x, y, t, z, w) = (3, 4, 5, 12, 13), (7, 24, 25, 312, 313), ... [Stark, p. 154] */ declare([t, z, w], integer)$ Time= 210 msecs (c30) system: [x^2 + y^2 = t^2, t^2 + z^2 = w^2]; Time= 0 msecs 2 2 2 2 2 2 (d30) [y + x = t , z + t = w ] (c31) solve(system, [x, y, t, z, w]); Time= 370 msecs 2 2 (d31) [[t = %r5, x = %r6, z = %r7, y = sqrt(%r5 - %r6 ), 2 2 2 2 w = sqrt(%r7 + %r5 )], [t = %r5, x = %r6, z = %r7, y = sqrt(%r5 - %r6 ), 2 2 2 2 w = - sqrt(%r7 + %r5 )], [t = %r5, x = %r6, z = %r7, y = - sqrt(%r5 - %r6 ), 2 2 2 2 w = sqrt(%r7 + %r5 )], [t = %r5, x = %r6, z = %r7, y = - sqrt(%r5 - %r6 ), 2 2 w = - sqrt(%r7 + %r5 )]] (c32) remvalue(system)$ Time= 0 msecs (c33) remove([t, w, x, t, z], integer)$ Time= 10 msecs (c34) /* Rational approximation of sqrt(3) with an error tolerance of 1/500 => 26/15 */ ev(rat(sqrt(3.)), ratepsilon = 1/500); RAT replaced 1.732051 by 26/15 = 1.733333 Time= 10 msecs 26 (d34)/R/ -- 15 (c35) /* Continued fractions => 3 + 1/(7 + 1/(15 + 1/(1 + 1/(292 + ... */ cf(3.1415926535d0); Time= 0 msecs (d35) [3, 7, 15, 1, 292] (c36) cfdisrep(%); Time= 0 msecs 1 (d36) 3 + ---------------- 1 7 + ------------ 1 15 + ------- 1 1 + --- 292 (c37) /* => 4 + 1/(1 + 1/(3 + 1/(1 + 1/(8 + 1/(1 + 1/(3 + 1/(1 + 1/(8 + ... [Stark, p. 340] */ cf(sqrt(23)); Time= 0 msecs (d37) [4, 1, 3, 1, 8] (c38) /* => 1 + 1/(1 + 1/(1 + 1/(1 + ... See Oskar Perron, _Die Lehre von den Kettenbr\"uchen_, Chelsea Publishing Company, 1950, p. 52. */ cf((1 + sqrt(5))/2); Time= 40 msecs (d38) [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2] (c39) /* => 1/(2 x + 1/(6 x + 1/(10 x + 1/(14 x + ... [Perron, p. 353] */ errcatch(cf((exp(1/x) - 1)/(exp(1/x) + 1))); Can't raise continued fraction to non-integral powers Time= 240 msecs (d39) [] (c40) /* => 1/(2 x + 1/(2 x + 1/(2 x + ... (Re x > 0) From Liyang Xu, ``Method Derived from Continued Fraction Approximations'', draft. */ errcatch(cf(sqrt(x^2 + 1) - x)); X - not a continued fraction Time= 10 msecs (d40) [] (c41) /* ---------- Quit ---------- */ quit(); Bye. real 5.92 user 2.76 sys 1.11