Use %%PRETTY-FNAME in more quadpack error messages
[maxima.git] / tests / wester_problems / test_number_theory.mac
blob6f7d18396ac826a0a6cddcd00a920dbf32fb67be
1 /* Original version of this file copyright 1999 by Michael Wester,
2  * and retrieved from http://www.math.unm.edu/~wester/demos/NumberTheory/problems.macsyma
3  * circa 2006-10-23.
4  *
5  * Released under the terms of the GNU General Public License, version 2,
6  * per message dated 2007-06-03 from Michael Wester to Robert Dodier
7  * (contained in the file wester-gpl-permission-message.txt).
8  *
9  * See: "A Critique of the Mathematical Abilities of CA Systems"
10  * by Michael Wester, pp 25--60 in
11  * "Computer Algebra Systems: A Practical Guide", edited by Michael J. Wester
12  * and published by John Wiley and Sons, Chichester, United Kingdom, 1999.
13  */
14 /* ----------[ M a c s y m a ]---------- */
15 /* ---------- Initialization ---------- */
16 showtime: all$
17 prederror: false$
18 /* ---------- Number Theory ---------- */
19 /* Display the largest 6-digit prime and the smallest 7-digit prime
20    => [999983, 1000003] */
21 [prime(78498), prime(78499)];
22 /* Primitive root => 19 */
23 191;
24 /* (a + b)^p mod p => a^p + b^p for p prime and a, b in Z_p   [Chris Hurlburt]
25    See Thomas W. Hungerford, _Algebra_, Springer-Verlag, 1974, p. 121 for a
26    more general simplification: (a +- b)^(p^n) => a^(p^n) +- b^(p^n) */
27 assume(equal(primep(p), true))$
28 errcatch(mod((a + b)^p, p));
29 forget(equal(primep(p), true))$
30 /* Congruence equations.  See Harold M. Stark, _An Introduction to Number
31    Theory_, The MIT press, 1984.
32    9 x = 15 mod 21 => x = 4 mod 7   or   {4, 11, 18} mod 21   [Stark, p. 68] */
33 modulus: 21$
34 errcatch(solve(9*x = 15, x));
35 /* 7 x = 22 mod 39 => x = 5 mod 13   or   31 mod 39   [Stark, p. 69] */
36 modulus: 39$
37 solve(7*x = 22, x);
38 /* x^2 + x + 4 = 0 mod 8 => x = {3, 4} mod 8   [Stark, p. 97] */
39 modulus: 8$
40 errcatch(solve(x^2 + x + 4 = 0, x));
41 /* x^3 + 2 x^2 + 5 x + 6 = 0 mod 11 => x = 3 mod 11   [Stark, p. 97] */
42 modulus: 11$
43 solve(x^3 + 2*x^2 + 5*x + 6 = 0, x);
44 multiplicities;
45 /* {x = 7 mod 9, x = 13 mod 23, x = 1 mod 2} => x = 151 mod 414   [Stark,
46    p. 76] */
47 /* {5 x + 4 y = 6 mod 7, 3 x - 2 y = 6 mod 7} => x = 1 mod 7, y = 2 mod 7
48    [Stark, p. 76] */
49 modulus: 7$
50 solve([5*x + 4*y = 6, 3*x - 2*y = 6], [x, y]);
51 /* 2 x + 3 y = 1 mod 5 =>
52    (x, y) = {(0, 2), (1, 3), (2, 4), (3, 0), (4, 1)} mod 5 */
53 modulus: 5$
54 solve(2*x + 3*y = 1, [x, y]);
55 /* 2 x + 3 y = 1 mod 6 =>   [Stark, p. 76]
56    (x, y) = {(2, 1), (2, 3), (2, 5), (5, 1), (5, 3), (5, 5)} mod 6 */
57 modulus: 6$
58 errcatch(solve(2*x + 3*y = 1, [x, y]));
59 modulus: false$
60 /* Diophantine equations => x = 2, y = 5 (Wallis)   [Stark, p. 147] */
61 declare([x, y], integer)$
62 solve(x^4 + 9 = y^2, [x, y]);
63 /* => x = 11, y = 5 (Fermat)   [Stark, p. 147] */
64 solve(x^2 + 4 = y^3, [x, y]);
65 /* => (x, y, t, z, w) = (3, 4, 5, 12, 13), (7, 24, 25, 312, 313), ...
66       [Stark, p. 154] */
67 declare([t, z, w], integer)$
68 system: [x^2 + y^2 = t^2, t^2 + z^2 = w^2];
69 solve(system, [x, y, t, z, w]);
70 remvalue(system)$
71 remove([t, w, x, t, z], integer)$
72 /* Rational approximation of sqrt(3) with an error tolerance of 1/500 => 26/15
73    */
74 ev(rat(sqrt(3.)), ratepsilon = 1/500);
75 /* Continued fractions => 3 + 1/(7 + 1/(15 + 1/(1 + 1/(292 + ... */
76 cf(3.1415926535d0);
77 cfdisrep(%);
78 /* => 4 + 1/(1 + 1/(3 + 1/(1 + 1/(8 + 1/(1 + 1/(3 + 1/(1 + 1/(8 + ...
79       [Stark, p. 340] */
80 cf(sqrt(23));
81 /* => 1 + 1/(1 + 1/(1 + 1/(1 + ...   See Oskar Perron, _Die Lehre von den
82       Kettenbr\"uchen_, Chelsea Publishing Company, 1950, p. 52. */
83 cf((1 + sqrt(5))/2);
84 /* => 1/(2 x + 1/(6 x + 1/(10 x + 1/(14 x + ...   [Perron, p. 353] */
85 errcatch(cf((exp(1/x) - 1)/(exp(1/x) + 1)));
86 /* => 1/(2 x + 1/(2 x + 1/(2 x + ...   (Re x > 0)   From Liyang Xu, ``Method
87       Derived from Continued Fraction Approximations'', draft. */
88 errcatch(cf(sqrt(x^2 + 1) - x));