| blob | d8c7f0be259ec7b8624092657c07dee52714c8fb |
1 @menu
2 * Functions and Variables for Number Theory::
3 @end menu
5 @c -----------------------------------------------------------------------------
6 @node Functions and Variables for Number Theory, , Number Theory, Number Theory
7 @section Functions and Variables for Number Theory
8 @c -----------------------------------------------------------------------------
10 @c -----------------------------------------------------------------------------
11 @anchor{bern}
12 @deffn {Function} bern (@var{n})
14 Returns the @var{n}'th Bernoulli number for integer @var{n}.
15 @c WELL, ACTUALLY bern SIMPLIFIES, LIKE FACTORIAL -- DO WE WANT TO GET INTO THAT ???
16 @c OR JUST PRETEND IT'S "RETURNED" ???
17 Bernoulli numbers equal to zero are suppressed if @code{zerobern} is
18 @code{false}.
20 See also @mrefdot{burn}
22 @example
23 (%i1) zerobern: true$
24 (%i2) map (bern, [0, 1, 2, 3, 4, 5, 6, 7, 8]);
25 1 1 1 1 1
26 (%o2) [1, - -, -, 0, - --, 0, --, 0, - --]
27 2 6 30 42 30
28 (%i3) zerobern: false$
29 (%i4) map (bern, [0, 1, 2, 3, 4, 5, 6, 7, 8]);
30 1 1 1 1 1 5 691 7
31 (%o4) [1, - -, -, - --, --, - --, --, - ----, -]
32 2 6 30 42 30 66 2730 6
33 @end example
35 @opencatbox{Categories:}
36 @category{Number theory}
37 @closecatbox
38 @end deffn
40 @c -----------------------------------------------------------------------------
41 @anchor{bernpoly}
42 @deffn {Function} bernpoly (@var{x}, @var{n})
44 Returns the @var{n}'th Bernoulli polynomial in the
45 variable @var{x}.
47 @opencatbox{Categories:}
48 @category{Number theory}
49 @closecatbox
50 @end deffn
52 @c -----------------------------------------------------------------------------
53 @anchor{bfzeta}
54 @deffn {Function} bfzeta (@var{s}, @var{n})
56 Returns the Riemann zeta function for the argument @var{s}.
57 The return value is a big float (bfloat);
58 @var{n} is the number of digits in the return value.
60 @opencatbox{Categories:}
61 @category{Number theory}
62 @category{Numerical evaluation}
63 @closecatbox
64 @end deffn
66 @c -----------------------------------------------------------------------------
67 @anchor{bfhzeta}
68 @deffn {Function} bfhzeta (@var{s}, @var{h}, @var{n})
70 Returns the Hurwitz zeta function for the arguments @var{s} and @var{h}.
71 The return value is a big float (bfloat);
72 @var{n} is the number of digits in the return value.
74 The Hurwitz zeta function is defined as
76 @tex
77 $$\zeta \left(s,h\right) = \sum_{k=0}^\infty {1 \over \left(k+h\right)^{s}}$$
78 @end tex
79 @ifnottex
80 @example
81 inf
82 ====
83 \ 1
84 zeta (s,h) = > --------
85 / s
86 ==== (k + h)
87 k = 0
88 @end example
89 @end ifnottex
91 @code{load ("bffac")} loads this function.
93 @opencatbox{Categories:}
94 @category{Number theory}
95 @category{Numerical evaluation}
96 @closecatbox
97 @end deffn
99 @c -----------------------------------------------------------------------------
100 @anchor{burn}
101 @deffn {Function} burn (@var{n})
103 Returns a rational number, which is an approximation of the @var{n}'th Bernoulli
104 number for integer @var{n}. @code{burn} exploits the observation that
105 (rational) Bernoulli numbers can be approximated by (transcendental) zetas with
106 tolerable efficiency:
108 @example
109 n - 1 1 - 2 n
110 (- 1) 2 zeta(2 n) (2 n)!
111 B(2 n) = ------------------------------------
112 2 n
113 %pi
114 @end example
116 @code{burn} may be more efficient than @mref{bern} for large, isolated @var{n}
117 as @mref{bern} computes all the Bernoulli numbers up to index @var{n} before
118 returning. @code{burn} invokes the approximation for even integers @code{n >
119 255}. For odd integers and @code{n <= 255} the function @mref{bern} is called.
121 @code{load ("bffac")} loads this function. See also @mrefdot{bern}
123 @opencatbox{Categories:}
124 @category{Number theory}
125 @closecatbox
126 @end deffn
128 @c -----------------------------------------------------------------------------
129 @anchor{chinese}
130 @deffn {Function} chinese ([@var{r_1}, @dots{}, @var{r_n}], [@var{m_1}, @dots{}, @var{m_n}])
132 Solves the system of congruences @code{x = r_1 mod m_1}, @dots{}, @code{x = r_n mod m_n}.
133 The remainders @var{r_n} may be arbitrary integers while the moduli @var{m_n} have to be
134 positive and pairwise coprime integers.
136 @example
137 (%i1) mods : [1000, 1001, 1003, 1007];
138 (%o1) [1000, 1001, 1003, 1007]
139 (%i2) lreduce('gcd, mods);
140 (%o2) 1
141 (%i3) x : random(apply("*", mods));
142 (%o3) 685124877004
143 (%i4) rems : map(lambda([z], mod(x, z)), mods);
144 (%o4) [4, 568, 54, 624]
145 (%i5) chinese(rems, mods);
146 (%o5) 685124877004
147 (%i6) chinese([1, 2], [3, n]);
148 (%o6) chinese([1, 2], [3, n])
149 (%i7) %, n = 4;
150 (%o7) 10
151 @end example
153 @opencatbox{Categories:}
154 @category{Number theory}
155 @closecatbox
156 @end deffn
158 @c -----------------------------------------------------------------------------
159 @anchor{cf}
160 @deffn {Function} cf (@var{expr})
162 Computes a continued fraction approximation.
163 @var{expr} is an expression comprising continued fractions,
164 square roots of integers, and literal real numbers
165 (integers, rational numbers, ordinary floats, and bigfloats).
166 @code{cf} computes exact expansions for rational numbers,
167 but expansions are truncated at @code{ratepsilon} for ordinary floats
168 and @code{10^(-fpprec)} for bigfloats.
170 Operands in the expression may be combined with arithmetic operators.
171 Maxima does not know about operations on continued fractions
172 outside of @code{cf}.
174 @code{cf} evaluates its arguments after binding @code{listarith} to
175 @code{false}. @code{cf} returns a continued fraction, represented as a list.
177 A continued fraction @code{a + 1/(b + 1/(c + ...))} is represented by the list
178 @code{[a, b, c, ...]}. The list elements @code{a}, @code{b}, @code{c}, @dots{}
179 must evaluate to integers. @var{expr} may also contain @code{sqrt (n)} where
180 @code{n} is an integer. In this case @code{cf} will give as many terms of the
181 continued fraction as the value of the variable @mref{cflength} times the
182 period.
184 A continued fraction can be evaluated to a number by evaluating the arithmetic
185 representation returned by @mrefdot{cfdisrep} See also @mref{cfexpand} for
186 another way to evaluate a continued fraction.
188 See also @mrefcomma{cfdisrep} @mrefcomma{cfexpand} and @mrefdot{cflength}
190 Examples:
192 @itemize @bullet
193 @item
194 @var{expr} is an expression comprising continued fractions and square roots of
195 integers.
197 @example
198 (%i1) cf ([5, 3, 1]*[11, 9, 7] + [3, 7]/[4, 3, 2]);
199 (%o1) [59, 17, 2, 1, 1, 1, 27]
200 (%i2) cf ((3/17)*[1, -2, 5]/sqrt(11) + (8/13));
201 (%o2) [0, 1, 1, 1, 3, 2, 1, 4, 1, 9, 1, 9, 2]
202 @end example
204 @item
205 @code{cflength} controls how many periods of the continued fraction
206 are computed for algebraic, irrational numbers.
208 @example
209 (%i1) cflength: 1$
210 (%i2) cf ((1 + sqrt(5))/2);
211 (%o2) [1, 1, 1, 1, 2]
212 (%i3) cflength: 2$
213 (%i4) cf ((1 + sqrt(5))/2);
214 (%o4) [1, 1, 1, 1, 1, 1, 1, 2]
215 (%i5) cflength: 3$
216 (%i6) cf ((1 + sqrt(5))/2);
217 (%o6) [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2]
218 @end example
220 @item
221 A continued fraction can be evaluated by evaluating the arithmetic
222 representation returned by @mrefdot{cfdisrep}
224 @example
225 (%i1) cflength: 3$
226 (%i2) cfdisrep (cf (sqrt (3)))$
227 (%i3) ev (%, numer);
228 (%o3) 1.731707317073171
229 @end example
231 @item
232 Maxima does not know about operations on continued fractions outside of
233 @code{cf}.
235 @example
236 (%i1) cf ([1,1,1,1,1,2] * 3);
237 (%o1) [4, 1, 5, 2]
238 (%i2) cf ([1,1,1,1,1,2]) * 3;
239 (%o2) [3, 3, 3, 3, 3, 6]
240 @end example
242 @end itemize
244 @opencatbox{Categories:}
245 @category{Continued fractions}
246 @closecatbox
247 @end deffn
249 @c NEEDS CLARIFICATION -- MAKE EXPLICIT HOW list IS RELATED TO a, b, c, ...
250 @c ALSO, CAN list CONTAIN ANYTHING OTHER THAN LITERAL INTEGERS ??
252 @c -----------------------------------------------------------------------------
253 @anchor{cfdisrep}
254 @deffn {Function} cfdisrep (@var{list})
256 Constructs and returns an ordinary arithmetic expression
257 of the form @code{a + 1/(b + 1/(c + ...))}
258 from the list representation of a continued fraction @code{[a, b, c, ...]}.
260 @example
261 (%i1) cf ([1, 2, -3] + [1, -2, 1]);
262 (%o1) [1, 1, 1, 2]
263 (%i2) cfdisrep (%);
264 1
265 (%o2) 1 + ---------
266 1
267 1 + -----
268 1
269 1 + -
270 2
271 @end example
273 @opencatbox{Categories:}
274 @category{Continued fractions}
275 @closecatbox
276 @end deffn
278 @c -----------------------------------------------------------------------------
279 @anchor{cfexpand}
280 @deffn {Function} cfexpand (@var{x})
282 Returns a matrix of the numerators and denominators of the last (column 1) and
283 next-to-last (column 2) convergents of the continued fraction @var{x}.
285 @example
286 (%i1) cf (rat (ev (%pi, numer)));
288 `rat' replaced 3.141592653589793 by 103993/33102 =3.141592653011902
289 (%o1) [3, 7, 15, 1, 292]
290 (%i2) cfexpand (%);
291 [ 103993 355 ]
292 (%o2) [ ]
293 [ 33102 113 ]
294 (%i3) %[1,1]/%[2,1], numer;
295 (%o3) 3.141592653011902
296 @end example
298 @opencatbox{Categories:}
299 @category{Continued fractions}
300 @closecatbox
301 @end deffn
303 @c -----------------------------------------------------------------------------
304 @anchor{cflength}
305 @defvr {Option variable} cflength
306 Default value: 1
308 @code{cflength} controls the number of terms of the continued fraction the
309 function @code{cf} will give, as the value @code{cflength} times the period.
310 Thus the default is to give one period.
312 @example
313 (%i1) cflength: 1$
314 (%i2) cf ((1 + sqrt(5))/2);
315 (%o2) [1, 1, 1, 1, 2]
316 (%i3) cflength: 2$
317 (%i4) cf ((1 + sqrt(5))/2);
318 (%o4) [1, 1, 1, 1, 1, 1, 1, 2]
319 (%i5) cflength: 3$
320 (%i6) cf ((1 + sqrt(5))/2);
321 (%o6) [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2]
322 @end example
324 @opencatbox{Categories:}
325 @category{Continued fractions}
326 @closecatbox
327 @end defvr
329 @c -----------------------------------------------------------------------------
330 @anchor{divsum}
331 @deffn {Function} divsum @
332 @fname{divsum} (@var{n}, @var{k}) @
333 @fname{divsum} (@var{n})
335 @code{divsum (@var{n}, @var{k})} returns the sum of the divisors of @var{n}
336 raised to the @var{k}'th power.
338 @code{divsum (@var{n})} returns the sum of the divisors of @var{n}.
340 @example
341 (%i1) divsum (12);
342 (%o1) 28
343 (%i2) 1 + 2 + 3 + 4 + 6 + 12;
344 (%o2) 28
345 (%i3) divsum (12, 2);
346 (%o3) 210
347 (%i4) 1^2 + 2^2 + 3^2 + 4^2 + 6^2 + 12^2;
348 (%o4) 210
349 @end example
351 @opencatbox{Categories:}
352 @category{Number theory}
353 @closecatbox
354 @end deffn
356 @c -----------------------------------------------------------------------------
357 @anchor{euler}
358 @deffn {Function} euler (@var{n})
360 Returns the @var{n}'th Euler number for nonnegative integer @var{n}.
361 Euler numbers equal to zero are suppressed if @code{zerobern} is
362 @code{false}.
364 For the Euler-Mascheroni constant, see @code{%gamma}.
366 @example
367 (%i1) zerobern: true$
368 (%i2) map (euler, [0, 1, 2, 3, 4, 5, 6]);
369 (%o2) [1, 0, - 1, 0, 5, 0, - 61]
370 (%i3) zerobern: false$
371 (%i4) map (euler, [0, 1, 2, 3, 4, 5, 6]);
372 (%o4) [1, - 1, 5, - 61, 1385, - 50521, 2702765]
373 @end example
375 @opencatbox{Categories:}
376 @category{Number theory}
377 @closecatbox
378 @end deffn
380 @c -----------------------------------------------------------------------------
381 @anchor{factors_only}
382 @defvr {Option variable} factors_only
383 Default value: @code{false}
385 Controls the value returned by @mrefdot{ifactors} The default @code{false}
386 causes @code{ifactors} to provide information about multiplicities of the
387 computed prime factors. If @code{factors_only} is set to @code{true},
388 @code{ifactors} returns nothing more than a list of prime factors.
390 Example: See @mrefdot{ifactors}
392 @opencatbox{Categories:}
393 @category{Number theory}
394 @closecatbox
395 @end defvr
397 @c -----------------------------------------------------------------------------
398 @anchor{fib}
399 @deffn {Function} fib (@var{n})
401 Returns the @var{n}'th Fibonacci number.
402 @code{fib(0)} is equal to 0 and @code{fib(1)} equal to 1, and
403 @code{fib (-@var{n})} equal to @code{(-1)^(@var{n} + 1) * fib(@var{n})}.
405 After calling @code{fib},
406 @code{prevfib} is equal to @code{fib(@var{n} - 1)},
407 the Fibonacci number preceding the last one computed.
409 @example
410 (%i1) map (fib, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8]);
411 (%o1) [- 3, 2, - 1, 1, 0, 1, 1, 2, 3, 5, 8, 13, 21]
412 @end example
414 @opencatbox{Categories:}
415 @category{Number theory}
416 @closecatbox
417 @end deffn
419 @c -----------------------------------------------------------------------------
420 @anchor{fibtophi}
421 @deffn {Function} fibtophi (@var{expr})
423 Expresses Fibonacci numbers in @var{expr} in terms of the constant @code{%phi},
424 which is @code{(1 + sqrt(5))/2}, approximately 1.61803399.
426 Examples:
428 @c ===beg===
429 @c fibtophi (fib (n));
430 @c fib (n-1) + fib (n) - fib (n+1);
431 @c fibtophi (%);
432 @c ratsimp (%);
433 @c ===end===
434 @example
435 (%i1) fibtophi (fib (n));
436 n n
437 %phi - (1 - %phi)
438 (%o1) -------------------
439 2 %phi - 1
440 (%i2) fib (n-1) + fib (n) - fib (n+1);
441 (%o2) - fib(n + 1) + fib(n) + fib(n - 1)
442 (%i3) fibtophi (%);
443 n + 1 n + 1 n n
444 %phi - (1 - %phi) %phi - (1 - %phi)
445 (%o3) - --------------------------- + -------------------
446 2 %phi - 1 2 %phi - 1
447 n - 1 n - 1
448 %phi - (1 - %phi)
449 + ---------------------------
450 2 %phi - 1
451 (%i4) ratsimp (%);
452 (%o4) 0
453 @end example
455 @opencatbox{Categories:}
456 @category{Number theory}
457 @closecatbox
458 @end deffn
460 @c -----------------------------------------------------------------------------
461 @anchor{ifactors}
462 @deffn {Function} ifactors (@var{n})
464 For a positive integer @var{n} returns the factorization of @var{n}. If
465 @code{n=p1^e1..pk^nk} is the decomposition of @var{n} into prime
466 factors, ifactors returns @code{[[p1, e1], ... , [pk, ek]]}.
468 Factorization methods used are trial divisions by primes up to 9973,
469 Pollard's rho and p-1 method and elliptic curves.
471 If the variable @code{ifactor_verbose} is set to @code{true}
472 ifactor produces detailed output about what it is doing including
473 immediate feedback as soon as a factor has been found.
475 The value returned by @code{ifactors} is controlled by the option variable @mrefdot{factors_only}
476 The default @code{false} causes @code{ifactors} to provide information about
477 the multiplicities of the computed prime factors. If @code{factors_only}
478 is set to @code{true}, @code{ifactors} simply returns the list of
479 prime factors.
481 @example
482 (%i1) ifactors(51575319651600);
483 (%o1) [[2, 4], [3, 2], [5, 2], [1583, 1], [9050207, 1]]
484 (%i2) apply("*", map(lambda([u], u[1]^u[2]), %));
485 (%o2) 51575319651600
486 (%i3) ifactors(51575319651600), factors_only : true;
487 (%o3) [2, 3, 5, 1583, 9050207]
488 @end example
490 @opencatbox{Categories:}
491 @category{Number theory}
492 @closecatbox
493 @end deffn
495 @c -----------------------------------------------------------------------------
496 @anchor{igcdex}
497 @deffn {Function} igcdex (@var{n}, @var{k})
499 Returns a list @code{[@var{a}, @var{b}, @var{u}]} where @var{u} is the greatest
500 common divisor of @var{n} and @var{k}, and @var{u} is equal to
501 @code{@var{a} @var{n} + @var{b} @var{k}}. The arguments @var{n} and @var{k}
502 must be integers.
504 @code{igcdex} implements the Euclidean algorithm. See also @mrefdot{gcdex}
506 The command @code{load("gcdex")} loads the function.
508 Examples:
510 @example
511 (%i1) load("gcdex")$
513 (%i2) igcdex(30,18);
514 (%o2) [- 1, 2, 6]
515 (%i3) igcdex(1526757668, 7835626735736);
516 (%o3) [845922341123, - 164826435, 4]
517 (%i4) igcdex(fib(20), fib(21));
518 (%o4) [4181, - 2584, 1]
519 @end example
521 @opencatbox{Categories:}
522 @category{Number theory}
523 @closecatbox
524 @end deffn
526 @c -----------------------------------------------------------------------------
527 @anchor{inrt}
528 @deffn {Function} inrt (@var{x}, @var{n})
530 Returns the integer @var{n}'th root of the absolute value of @var{x}.
532 @example
533 (%i1) l: [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]$
534 (%i2) map (lambda ([a], inrt (10^a, 3)), l);
535 (%o2) [2, 4, 10, 21, 46, 100, 215, 464, 1000, 2154, 4641, 10000]
536 @end example
538 @opencatbox{Categories:}
539 @category{Number theory}
540 @closecatbox
541 @end deffn
543 @c -----------------------------------------------------------------------------
544 @anchor{inv_mod}
545 @deffn {Function} inv_mod (@var{n}, @var{m})
547 Computes the inverse of @var{n} modulo @var{m}.
548 @code{inv_mod (n,m)} returns @code{false},
549 if @var{n} is a zero divisor modulo @var{m}.
551 @example
552 (%i1) inv_mod(3, 41);
553 (%o1) 14
554 (%i2) ratsimp(3^-1), modulus = 41;
555 (%o2) 14
556 (%i3) inv_mod(3, 42);
557 (%o3) false
558 @end example
560 @opencatbox{Categories:}
561 @category{Number theory}
562 @closecatbox
563 @end deffn
565 @c -----------------------------------------------------------------------------
566 @anchor{isqrt}
567 @deffn {Function} isqrt (@var{x})
569 Returns the "integer square root" of the absolute value of @var{x}, which is an
570 integer.
572 @opencatbox{Categories:}
573 @category{Mathematical functions}
574 @closecatbox
575 @end deffn
577 @c -----------------------------------------------------------------------------
578 @anchor{jacobi}
579 @deffn {Function} jacobi (@var{p}, @var{q})
581 Returns the Jacobi symbol of @var{p} and @var{q}.
583 @example
584 (%i1) l: [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]$
585 (%i2) map (lambda ([a], jacobi (a, 9)), l);
586 (%o2) [1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0]
587 @end example
589 @opencatbox{Categories:}
590 @category{Number theory}
591 @closecatbox
592 @end deffn
594 @c -----------------------------------------------------------------------------
595 @anchor{lcm}
596 @deffn {Function} lcm (@var{expr_1}, @dots{}, @var{expr_n})
598 Returns the least common multiple of its arguments.
599 The arguments may be general expressions as well as integers.
601 @code{load ("functs")} loads this function.
603 @opencatbox{Categories:}
604 @category{Number theory}
605 @closecatbox
606 @end deffn
608 @c -----------------------------------------------------------------------------
609 @anchor{lucas}
610 @deffn {Function} lucas (@var{n})
612 Returns the @var{n}'th Lucas number.
613 @code{lucas(0)} is equal to 2 and @code{lucas(1)} equal to 1, and
614 @code{lucas(-@var{n})} equal to @code{(-1)^(-@var{n}) * lucas(@var{n})}.
616 @example
617 (%i1) map (lucas, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8]);
618 (%o1) [7, - 4, 3, - 1, 2, 1, 3, 4, 7, 11, 18, 29, 47]
619 @end example
621 After calling @code{lucas}, the global variable
622 @code{next_lucas} is equal to @code{lucas (@var{n} + 1)},
623 the Lucas number following the last returned. The example shows
624 how Fibonacci numbers can be computed via @code{lucas} and @code{next_lucas}.
626 @example
627 (%i1) fib_via_lucas(n) :=
628 block([lucas : lucas(n)],
629 signum(n) * (2*next_lucas - lucas)/5 )$
630 (%i2) map (fib_via_lucas, [-4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7,
631 8]);
632 (%o2) [- 3, 2, - 1, 1, 0, 1, 1, 2, 3, 5, 8, 13, 21]
633 @end example
635 @opencatbox{Categories:}
636 @category{Number theory}
637 @closecatbox
638 @end deffn
640 @c -----------------------------------------------------------------------------
641 @anchor{mod}
642 @deffn {Function} mod (@var{x}, @var{y})
644 If @var{x} and @var{y} are real numbers and @var{y} is nonzero, return
645 @code{@var{x} - @var{y} * floor(@var{x} / @var{y})}. Further for all real
646 @var{x}, we have @code{mod (@var{x}, 0) = @var{x}}. For a discussion of the
647 definition @code{mod (@var{x}, 0) = @var{x}}, see Section 3.4, of
648 "Concrete Mathematics," by Graham, Knuth, and Patashnik. The function
649 @code{mod (@var{x}, 1)} is a sawtooth function with period 1 with
650 @code{mod (1, 1) = 0} and @code{mod (0, 1) = 0}.
652 To find the principal argument (a number in the interval @code{(-%pi, %pi]}) of
653 a complex number, use the function
654 @code{@var{x} |-> %pi - mod (%pi - @var{x}, 2*%pi)}, where @var{x} is an
655 argument.
657 When @var{x} and @var{y} are constant expressions (@code{10 * %pi}, for
658 example), @code{mod} uses the same big float evaluation scheme that @code{floor}
659 and @code{ceiling} uses. Again, it's possible, although unlikely, that
660 @code{mod} could return an erroneous value in such cases.
662 For nonnumerical arguments @var{x} or @var{y}, @code{mod} knows several
663 simplification rules:
665 @c ===beg===
666 @c mod (x, 0);
667 @c mod (a*x, a*y);
668 @c mod (0, x);
669 @c ===end===
670 @example
671 (%i1) mod (x, 0);
672 (%o1) x
673 (%i2) mod (a*x, a*y);
674 (%o2) a mod(x, y)
675 (%i3) mod (0, x);
676 (%o3) 0
677 @end example
679 @opencatbox{Categories:}
680 @category{Mathematical functions}
681 @closecatbox
682 @end deffn
684 @c -----------------------------------------------------------------------------
685 @anchor{next_prime}
686 @deffn {Function} next_prime (@var{n})
688 Returns the smallest prime bigger than @var{n}.
690 @example
691 (%i1) next_prime(27);
692 (%o1) 29
693 @end example
695 @opencatbox{Categories:}
696 @category{Number theory}
697 @closecatbox
698 @end deffn
700 @c -----------------------------------------------------------------------------
701 @anchor{partfrac}
702 @deffn {Function} partfrac (@var{expr}, @var{var})
704 Expands the expression @var{expr} in partial fractions
705 with respect to the main variable @var{var}. @code{partfrac} does a complete
706 partial fraction decomposition. The algorithm employed is based on
707 the fact that the denominators of the partial fraction expansion (the
708 factors of the original denominator) are relatively prime. The
709 numerators can be written as linear combinations of denominators, and
710 the expansion falls out.
712 @code{partfrac} ignores the value @code{true} of the option variable
713 @code{keepfloat}.
715 @example
716 (%i1) 1/(1+x)^2 - 2/(1+x) + 2/(2+x);
717 2 2 1
718 (%o1) ----- - ----- + --------
719 x + 2 x + 1 2
720 (x + 1)
721 (%i2) ratsimp (%);
722 x
723 (%o2) - -------------------
724 3 2
725 x + 4 x + 5 x + 2
726 (%i3) partfrac (%, x);
727 2 2 1
728 (%o3) ----- - ----- + --------
729 x + 2 x + 1 2
730 (x + 1)
731 @end example
732 @end deffn
734 @c -----------------------------------------------------------------------------
735 @anchor{power_mod}
736 @deffn {Function} power_mod (@var{a}, @var{n}, @var{m})
738 Uses a modular algorithm to compute @code{a^n mod m}
739 where @var{a} and @var{n} are integers and @var{m} is a positive integer.
740 If @var{n} is negative, @code{inv_mod} is used to find the modular inverse.
742 @example
743 (%i1) power_mod(3, 15, 5);
744 (%o1) 2
745 (%i2) mod(3^15,5);
746 (%o2) 2
747 (%i3) power_mod(2, -1, 5);
748 (%o3) 3
749 (%i4) inv_mod(2,5);
750 (%o4) 3
751 @end example
753 @opencatbox{Categories:}
754 @category{Number theory}
755 @closecatbox
756 @end deffn
758 @c -----------------------------------------------------------------------------
759 @anchor{primep}
760 @deffn {Function} primep (@var{n})
762 Primality test. If @code{primep (@var{n})} returns @code{false}, @var{n} is a
763 composite number and if it returns @code{true}, @var{n} is a prime number
764 with very high probability.
766 For @var{n} less than 341550071728321 a deterministic version of
767 Miller-Rabin's test is used. If @code{primep (@var{n})} returns
768 @code{true}, then @var{n} is a prime number.
770 For @var{n} bigger than 341550071728321 @code{primep} uses
771 @code{primep_number_of_tests} Miller-Rabin's pseudo-primality tests and one
772 Lucas pseudo-primality test. The probability that a non-prime @var{n} will
773 pass one Miller-Rabin test is less than 1/4. Using the default value 25 for
774 @code{primep_number_of_tests}, the probability of @var{n} being
775 composite is much smaller that 10^-15.
777 @opencatbox{Categories:}
778 @category{Predicate functions}
779 @category{Number theory}
780 @closecatbox
781 @end deffn
783 @c -----------------------------------------------------------------------------
784 @anchor{primep_number_of_tests}
785 @defvr {Option variable} primep_number_of_tests
786 Default value: 25
788 Number of Miller-Rabin's tests used in @code{primep}.
790 @opencatbox{Categories:}
791 @category{Number theory}
792 @closecatbox
793 @end defvr
795 @c -----------------------------------------------------------------------------
796 @anchor{primes}
797 @deffn {Function} primes (@var{start}, @var{end})
799 Returns the list of all primes from @var{start} to @var{end}.
801 @example
802 (%i1) primes(3, 7);
803 (%o1) [3, 5, 7]
804 @end example
806 @opencatbox{Categories:}
807 @category{Number theory}
808 @closecatbox
809 @end deffn
811 @c -----------------------------------------------------------------------------
812 @anchor{prev_time}
813 @deffn {Function} prev_prime (@var{n})
815 Returns the greatest prime smaller than @var{n}.
817 @example
818 (%i1) prev_prime(27);
819 (%o1) 23
820 @end example
822 @opencatbox{Categories:}
823 @category{Number theory}
824 @closecatbox
825 @end deffn
827 @c -----------------------------------------------------------------------------
828 @anchor{qunit}
829 @deffn {Function} qunit (@var{n})
831 Returns the principal unit of the real quadratic number field
832 @code{sqrt (@var{n})} where @var{n} is an integer,
833 i.e., the element whose norm is unity.
834 This amounts to solving Pell's equation @code{a^2 - @var{n} b^2 = 1}.
836 @example
837 (%i1) qunit (17);
838 (%o1) sqrt(17) + 4
839 (%i2) expand (% * (sqrt(17) - 4));
840 (%o2) 1
841 @end example
843 @opencatbox{Categories:}
844 @category{Number theory}
845 @closecatbox
846 @end deffn
848 @c -----------------------------------------------------------------------------
849 @anchor{totient}
850 @deffn {Function} totient (@var{n})
852 Returns the number of integers less than or equal to @var{n} which
853 are relatively prime to @var{n}.
855 @opencatbox{Categories:}
856 @category{Number theory}
857 @closecatbox
858 @end deffn
860 @c -----------------------------------------------------------------------------
861 @defvr {Option variable} zerobern
862 Default value: @code{true}
864 When @code{zerobern} is @code{false}, @code{bern} excludes the Bernoulli numbers
865 and @code{euler} excludes the Euler numbers which are equal to zero.
866 See @code{bern} and @code{euler}.
868 @opencatbox{Categories:}
869 @category{Number theory}
870 @closecatbox
871 @end defvr
873 @c -----------------------------------------------------------------------------
874 @anchor{zeta}
875 @deffn {Function} zeta (@var{n})
877 Returns the Riemann zeta function. If @var{n} is a negative integer, 0, or a
878 positive even integer, the Riemann zeta function simplifies to an exact value.
879 For a positive even integer the option variable @code{zeta%pi} has to be
880 @code{true} in addition (See @code{zeta%pi}). For a floating point or bigfloat
881 number the Riemann zeta function is evaluated numerically. Maxima returns a
882 noun form @code{zeta (@var{n})} for all other arguments, including rational
883 noninteger, and complex arguments, or for even integers, if @code{zeta%pi} has
884 the value @code{false}.
886 @code{zeta(1)} is undefined, but Maxima knows the limit
887 @code{limit(zeta(x), x, 1)} from above and below.
889 The Riemann zeta function distributes over lists, matrices, and equations.
891 See also @mref{bfzeta} and @mrefdot{zeta%pi}
893 Examples:
895 @c ===beg===
896 @c zeta([-2, -1, 0, 0.5, 2, 3,1+%i]);
897 @c limit(zeta(x),x,1,plus);
898 @c limit(zeta(x),x,1,minus);
899 @c ===end===
900 @example
901 (%i1) zeta([-2, -1, 0, 0.5, 2, 3, 1+%i]);
902 2
903 1 1 %pi
904 (%o1) [0, - --, - -, - 1.460354508809586, ----, zeta(3),
905 12 2 6
906 zeta(%i + 1)]
907 (%i2) limit(zeta(x),x,1,plus);
908 (%o2) inf
909 (%i3) limit(zeta(x),x,1,minus);
910 (%o3) minf
911 @end example
913 @opencatbox{Categories:}
914 @category{Number theory}
915 @closecatbox
916 @end deffn
918 @c -----------------------------------------------------------------------------
919 @anchor{zeta%pi}
920 @defvr {Option variable} zeta%pi
921 Default value: @mref{true}
923 When @code{zeta%pi} is @code{true}, @code{zeta} returns an expression
924 proportional to @code{%pi^n} for even integer @code{n}. Otherwise, @code{zeta}
925 returns a noun form @code{zeta (n)} for even integer @code{n}.
927 Examples:
929 @c ===beg===
930 @c zeta%pi: true$
931 @c zeta (4);
932 @c zeta%pi: false$
933 @c zeta (4);
934 @c ===end===
935 @example
936 (%i1) zeta%pi: true$
937 (%i2) zeta (4);
938 4
939 %pi
940 (%o2) ----
941 90
942 (%i3) zeta%pi: false$
943 (%i4) zeta (4);
944 (%o4) zeta(4)
945 @end example
947 @opencatbox{Categories:}
948 @category{Number theory}
949 @closecatbox
950 @end defvr
952 @c -----------------------------------------------------------------------------
953 @anchor{zn_add_table}
954 @deffn {Function} zn_add_table (@var{n})
956 Shows an addition table of all elements in (Z/@var{n}Z).
958 See also @mrefcomma{zn_mult_table} @mrefdot{zn_power_table}
960 @opencatbox{Categories:}
961 @category{Number theory}
962 @closecatbox
963 @end deffn
965 @c -----------------------------------------------------------------------------
966 @anchor{zn_characteristic_factors}
967 @deffn {Function} zn_characteristic_factors (@var{n})
969 Returns a list containing the characteristic factors of the totient of @var{n}.
971 Using the characteristic factors a multiplication group modulo @var{n}
972 can be expressed as a group direct product of cyclic subgroups.
974 In case the group itself is cyclic the list only contains the totient
975 and using @code{zn_primroot} a generator can be computed.
976 If the totient splits into more than one characteristic factors
977 @code{zn_factor_generators} finds generators of the corresponding subgroups.
979 Each of the @code{r} factors in the list divides the right following factors.
980 For the last factor @code{f_r} therefore holds @code{a^f_r = 1 (mod n)}
981 for all @code{a} coprime to @var{n}.
982 This factor is also known as Carmichael function or Carmichael lambda.
984 If @code{n > 2}, then @code{totient(n)/2^r} is the number of quadratic residues,
985 and each of these has @code{2^r} square roots.
987 See also @mrefcomma{totient} @mrefcomma{zn_primroot} @mrefdot{zn_factor_generators}
989 Examples:
991 The multiplication group modulo @code{14} is cyclic and its @code{6} elements
992 can be generated by a primitive root.
994 @example
995 (%i1) [zn_characteristic_factors(14), phi: totient(14)];
996 (%o1) [[6], 6]
997 (%i2) [zn_factor_generators(14), g: zn_primroot(14)];
998 (%o2) [[3], 3]
999 (%i3) M14: makelist(power_mod(g,i,14), i,0,phi-1);
1000 (%o3) [1, 3, 9, 13, 11, 5]
1001 @end example
1003 The multiplication group modulo @code{15} is not cyclic and its @code{8} elements
1004 can be generated by two factor generators.
1006 @example
1007 (%i1) [[f1,f2]: zn_characteristic_factors(15), totient(15)];
1008 (%o1) [[2, 4], 8]
1009 (%i2) [[g1,g2]: zn_factor_generators(15), zn_primroot(15)];
1010 (%o2) [[11, 7], false]
1011 (%i3) UG1: makelist(power_mod(g1,i,15), i,0,f1-1);
1012 (%o3) [1, 11]
1013 (%i4) UG2: makelist(power_mod(g2,i,15), i,0,f2-1);
1014 (%o4) [1, 7, 4, 13]
1015 (%i5) M15: create_list(mod(i*j,15), i,UG1, j,UG2);
1016 (%o5) [1, 7, 4, 13, 11, 2, 14, 8]
1017 @end example
1019 For the last characteristic factor @code{4} it holds that @code{a^4 = 1 (mod 15)}
1020 for all @code{a} in @code{M15}.
1022 @code{M15} has two characteristic factors and therefore @code{8/2^2} quadratic residues,
1023 and each of these has @code{2^2} square roots.
1025 @example
1026 (%i6) zn_power_table(15);
1027 [ 1 1 1 1 ]
1028 [ ]
1029 [ 2 4 8 1 ]
1030 [ ]
1031 [ 4 1 4 1 ]
1032 [ ]
1033 [ 7 4 13 1 ]
1034 (%o6) [ ]
1035 [ 8 4 2 1 ]
1036 [ ]
1037 [ 11 1 11 1 ]
1038 [ ]
1039 [ 13 4 7 1 ]
1040 [ ]
1041 [ 14 1 14 1 ]
1042 (%i7) map(lambda([i], zn_nth_root(i,2,15)), [1,4]);
1043 (%o7) [[1, 4, 11, 14], [2, 7, 8, 13]]
1044 @end example
1046 @opencatbox{Categories:}
1047 @category{Number theory}
1048 @closecatbox
1049 @end deffn
1051 @c -----------------------------------------------------------------------------
1052 @anchor{zn_carmichael_lambda}
1053 @deffn {Function} zn_carmichael_lambda (@var{n})
1055 Returns @code{1} if @var{n} is @code{1} and otherwise
1056 the greatest characteristic factor of the totient of @var{n}.
1058 For remarks and examples see @mrefdot{zn_characteristic_factors}
1060 @opencatbox{Categories:}
1061 @category{Number theory}
1062 @closecatbox
1063 @end deffn
1065 @c -----------------------------------------------------------------------------
1066 @anchor{zn_determinant}
1067 @deffn {Function} zn_determinant (@var{matrix}, @var{p})
1069 Uses the technique of LU-decomposition to compute the determinant of @var{matrix}
1070 over (Z/@var{p}Z). @var{p} must be a prime.
1072 However if the determinant is equal to zero the LU-decomposition might fail.
1073 In that case @code{zn_determinant} computes the determinant non-modular
1074 and reduces thereafter.
1076 See also @mrefdot{zn_invert_by_lu}
1078 Examples:
1080 @example
1081 (%i1) m : matrix([1,3],[2,4]);
1082 [ 1 3 ]
1083 (%o1) [ ]
1084 [ 2 4 ]
1085 (%i2) zn_determinant(m, 5);
1086 (%o2) 3
1087 (%i3) m : matrix([2,4,1],[3,1,4],[4,3,2]);
1088 [ 2 4 1 ]
1089 [ ]
1090 (%o3) [ 3 1 4 ]
1091 [ ]
1092 [ 4 3 2 ]
1093 (%i4) zn_determinant(m, 5);
1094 (%o4) 0
1095 @end example
1097 @opencatbox{Categories:}
1098 @category{Number theory}
1099 @closecatbox
1100 @end deffn
1102 @c -----------------------------------------------------------------------------
1103 @anchor{zn_factor_generators}
1104 @deffn {Function} zn_factor_generators (@var{n})
1106 Returns a list containing factor generators corresponding to the
1107 characteristic factors of the totient of @var{n}.
1109 For remarks and examples see @mrefdot{zn_characteristic_factors}
1111 @opencatbox{Categories:}
1112 @category{Number theory}
1113 @closecatbox
1114 @end deffn
1116 @c -----------------------------------------------------------------------------
1117 @anchor{zn_invert_by_lu}
1118 @deffn {Function} zn_invert_by_lu (@var{matrix}, @var{p})
1120 Uses the technique of LU-decomposition to compute the modular inverse of
1121 @var{matrix} over (Z/@var{p}Z). @var{p} must be a prime and @var{matrix}
1122 invertible. @code{zn_invert_by_lu} returns @code{false} if @var{matrix}
1123 is not invertible.
1125 See also @mrefdot{zn_determinant}
1127 Example:
1129 @example
1130 (%i1) m : matrix([1,3],[2,4]);
1131 [ 1 3 ]
1132 (%o1) [ ]
1133 [ 2 4 ]
1134 (%i2) zn_determinant(m, 5);
1135 (%o2) 3
1136 (%i3) mi : zn_invert_by_lu(m, 5);
1137 [ 3 4 ]
1138 (%o3) [ ]
1139 [ 1 2 ]
1140 (%i4) matrixmap(lambda([a], mod(a, 5)), m . mi);
1141 [ 1 0 ]
1142 (%o4) [ ]
1143 [ 0 1 ]
1144 @end example
1146 @opencatbox{Categories:}
1147 @category{Number theory}
1148 @closecatbox
1149 @end deffn
1151 @c -----------------------------------------------------------------------------
1152 @anchor{zn_log}
1153 @deffn {Function} zn_log @
1154 @fname{zn_log} (@var{a}, @var{g}, @var{n}) @
1155 @fname{zn_log} (@var{a}, @var{g}, @var{n}, [[@var{p1}, @var{e1}], @dots{}, [@var{pk}, @var{ek}]])
1157 Computes the discrete logarithm. Let (Z/@var{n}Z)* be a cyclic group, @var{g} a
1158 primitive root modulo @var{n} or a generator of a subgroup of (Z/@var{n}Z)*
1159 and let @var{a} be a member of this group.
1160 @code{zn_log (a, g, n)} then solves the congruence @code{g^x = a mod n}.
1161 Please note that if @var{a} is not a power of @var{g} modulo @var{n},
1162 @code{zn_log} will not terminate.
1164 The applied algorithm needs a prime factorization of @code{zn_order(g)} resp. @code{totient(n)}
1165 in case @var{g} is a primitive root modulo @var{n}.
1166 A precomputed list of factors of @code{zn_order(g)} might be used as the optional fourth argument.
1167 This list must be of the same form as the list returned by @code{ifactors(zn_order(g))}
1168 using the default option @code{factors_only : false}.
1169 However, compared to the running time of the logarithm algorithm
1170 providing the list of factors has only a quite small effect.
1172 The algorithm uses a Pohlig-Hellman-reduction and Pollard's Rho-method for
1173 discrete logarithms. The running time of @code{zn_log} primarily depends on the
1174 bitlength of the greatest prime factor of @code{zn_order(g)}.
1176 See also @mrefcomma{zn_primroot} @mrefcomma{zn_order} @mrefcomma{ifactors} @mrefdot{totient}
1178 Examples:
1180 @code{zn_log (a, g, n)} solves the congruence @code{g^x = a mod n}.
1182 @example
1183 (%i1) n : 22$
1184 (%i2) g : zn_primroot(n);
1185 (%o2) 7
1186 (%i3) ord_7 : zn_order(7, n);
1187 (%o3) 10
1188 (%i4) powers_7 : makelist(power_mod(g, x, n), x, 0, ord_7 - 1);
1189 (%o4) [1, 7, 5, 13, 3, 21, 15, 17, 9, 19]
1190 (%i5) zn_log(9, g, n);
1191 (%o5) 8
1192 (%i6) map(lambda([x], zn_log(x, g, n)), powers_7);
1193 (%o6) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
1194 (%i7) ord_5 : zn_order(5, n);
1195 (%o7) 5
1196 (%i8) powers_5 : makelist(power_mod(5,x,n), x, 0, ord_5 - 1);
1197 (%o8) [1, 5, 3, 15, 9]
1198 (%i9) zn_log(9, 5, n);
1199 (%o9) 4
1200 @end example
1202 The optional fourth argument must be of the same form as the list returned by
1203 @code{ifactors(zn_order(g))}.
1204 The running time primarily depends on the bitlength of the totient's greatest prime factor.
1206 @example
1207 (%i1) (p : 2^127-1, primep(p));
1208 (%o1) true
1209 (%i2) ifs : ifactors(p - 1)$
1210 (%i3) g : zn_primroot(p, ifs);
1211 (%o3) 43
1212 (%i4) a : power_mod(g, 4711, p)$
1213 (%i5) zn_log(a, g, p, ifs);
1214 (%o5) 4711
1215 (%i6) f_max : last(ifs);
1216 (%o6) [77158673929, 1]
1217 (%i7) ord_5 : zn_order(5,p,ifs)$
1218 (%i8) (p - 1)/ord_5;
1219 (%o8) 73
1220 (%i9) ifs_5 : ifactors(ord_5)$
1221 (%i10) a : power_mod(5, 4711, p)$
1222 (%i11) zn_log(a, 5, p, ifs_5);
1223 (%o11) 4711
1224 @end example
1226 @opencatbox{Categories:}
1227 @category{Number theory}
1228 @closecatbox
1229 @end deffn
1231 @c -----------------------------------------------------------------------------
1232 @anchor{zn_mult_table}
1233 @deffn {Function} zn_mult_table @
1234 @fname{zn_mult_table} (@var{n}) @
1235 @fname{zn_mult_table} (@var{n}, @var{gcd})
1237 Without the optional argument @var{gcd} @code{zn_mult_table(n)} shows a
1238 multiplication table of all elements in (Z/@var{n}Z)* which are all elements
1239 coprime to @var{n}.
1241 The optional second argument @var{gcd} allows to select a specific
1242 subset of (Z/@var{n}Z). If @var{gcd} is an integer, a multiplication table of
1243 all residues @code{x} with @code{gcd(x,n) = }@var{gcd} are returned.
1244 Additionally row and column headings are added for better readability.
1245 If necessary, these can be easily removed by @code{submatrix(1, table, 1)}.
1247 If @var{gcd} is set to @code{all}, the table is printed for all non-zero
1248 elements in (Z/@var{n}Z).
1250 The second example shows an alternative way to create a multiplication table
1251 for subgroups.
1253 See also @mrefcomma{zn_add_table} @mrefdot{zn_power_table}
1255 Examples:
1257 The default table shows all elements in (Z/@var{n}Z)* and allows to demonstrate
1258 and study basic properties of modular multiplication groups.
1259 E.g. the principal diagonal contains all quadratic residues,
1260 each row and column contains every element, the tables are symmetric, etc..
1262 If @var{gcd} is set to @code{all}, the table is printed for all non-zero
1263 elements in (Z/@var{n}Z).
1265 @example
1266 (%i1) zn_mult_table(8);
1267 [ 1 3 5 7 ]
1268 [ ]
1269 [ 3 1 7 5 ]
1270 (%o1) [ ]
1271 [ 5 7 1 3 ]
1272 [ ]
1273 [ 7 5 3 1 ]
1274 (%i2) zn_mult_table(8, all);
1275 [ 1 2 3 4 5 6 7 ]
1276 [ ]
1277 [ 2 4 6 0 2 4 6 ]
1278 [ ]
1279 [ 3 6 1 4 7 2 5 ]
1280 [ ]
1281 (%o2) [ 4 0 4 0 4 0 4 ]
1282 [ ]
1283 [ 5 2 7 4 1 6 3 ]
1284 [ ]
1285 [ 6 4 2 0 6 4 2 ]
1286 [ ]
1287 [ 7 6 5 4 3 2 1 ]
1288 @end example
1290 If @var{gcd} is an integer, row and column headings are added for better readability.
1292 If the subset chosen by @var{gcd} is a group there is another way to create
1293 a multiplication table. An isomorphic mapping from a group with @code{1} as
1294 identity builds a table which is easy to read. The mapping is accomplished via CRT.
1296 In the second version of @code{T36_4} the identity, here @code{28}, is placed in
1297 the top left corner, just like in table @code{T9}.
1299 @example
1300 (%i1) T36_4: zn_mult_table(36,4);
1301 [ * 4 8 16 20 28 32 ]
1302 [ ]
1303 [ 4 16 32 28 8 4 20 ]
1304 [ ]
1305 [ 8 32 28 20 16 8 4 ]
1306 [ ]
1307 (%o1) [ 16 28 20 4 32 16 8 ]
1308 [ ]
1309 [ 20 8 16 32 4 20 28 ]
1310 [ ]
1311 [ 28 4 8 16 20 28 32 ]
1312 [ ]
1313 [ 32 20 4 8 28 32 16 ]
1314 (%i2) T9: zn_mult_table(36/4);
1315 [ 1 2 4 5 7 8 ]
1316 [ ]
1317 [ 2 4 8 1 5 7 ]
1318 [ ]
1319 [ 4 8 7 2 1 5 ]
1320 (%o2) [ ]
1321 [ 5 1 2 7 8 4 ]
1322 [ ]
1323 [ 7 5 1 8 4 2 ]
1324 [ ]
1325 [ 8 7 5 4 2 1 ]
1326 (%i3) T36_4: matrixmap(lambda([x], chinese([0,x],[4,9])), T9);
1327 [ 28 20 4 32 16 8 ]
1328 [ ]
1329 [ 20 4 8 28 32 16 ]
1330 [ ]
1331 [ 4 8 16 20 28 32 ]
1332 (%o3) [ ]
1333 [ 32 28 20 16 8 4 ]
1334 [ ]
1335 [ 16 32 28 8 4 20 ]
1336 [ ]
1337 [ 8 16 32 4 20 28 ]
1338 @end example
1340 @opencatbox{Categories:}
1341 @category{Number theory}
1342 @closecatbox
1343 @end deffn
1345 @c -----------------------------------------------------------------------------
1346 @anchor{zn_nth_root}
1347 @deffn {Function} zn_nth_root @
1348 @fname{zn_nth_root} (@var{x}, @var{n}, @var{m}) @
1349 @fname{zn_nth_root} (@var{x}, @var{n}, @var{m}, [[@var{p1}, @var{e1}], @dots{}, [@var{pk}, @var{ek}]])
1351 Returns a list with all @var{n}-th roots of @var{x} from the multiplication
1352 subgroup of (Z/@var{m}Z) which contains @var{x}, or @code{false}, if @var{x}
1353 is no @var{n}-th power modulo @var{m} or not contained in any multiplication
1354 subgroup of (Z/@var{m}Z).
1356 @var{x} is an element of a multiplication subgroup modulo @var{m}, if the
1357 greatest common divisor @code{g = gcd(x,m)} is coprime to @code{m/g}.
1359 @code{zn_nth_root} is based on an algorithm by Adleman, Manders and Miller
1360 and on theorems about modulo multiplication groups by Daniel Shanks.
1362 The algorithm needs a prime factorization of the modulus @var{m}.
1363 So in case the factorization of @var{m} is known, the list of factors
1364 can be passed as the fourth argument. This optional argument
1365 must be of the same form as the list returned by @code{ifactors(m)}
1366 using the default option @code{factors_only: false}.
1369 Examples:
1371 A power table of the multiplication group modulo @code{14}
1372 followed by a list of lists containing all @var{n}-th roots of @code{1}
1373 with @var{n} from @code{1} to @code{6}.
1375 @example
1376 (%i1) zn_power_table(14);
1377 [ 1 1 1 1 1 1 ]
1378 [ ]
1379 [ 3 9 13 11 5 1 ]
1380 [ ]
1381 [ 5 11 13 9 3 1 ]
1382 (%o1) [ ]
1383 [ 9 11 1 9 11 1 ]
1384 [ ]
1385 [ 11 9 1 11 9 1 ]
1386 [ ]
1387 [ 13 1 13 1 13 1 ]
1388 (%i2) makelist(zn_nth_root(1,n,14), n,1,6);
1389 (%o2) [[1], [1, 13], [1, 9, 11], [1, 13], [1], [1, 3, 5, 9, 11, 13]]
1390 @end example
1392 In the following example @var{x} is not coprime to @var{m},
1393 but is a member of a multiplication subgroup of (Z/@var{m}Z)
1394 and any @var{n}-th root is a member of the same subgroup.
1396 The residue class @code{3} is no member of any multiplication subgroup of (Z/63Z)
1397 and is therefore not returned as a third root of @code{27}.
1399 Here @code{zn_power_table} shows all residues @code{x} in (Z/63Z)
1400 with @code{gcd(x,63) = 9}. This subgroup is isomorphic to (Z/7Z)*
1401 and its identity @code{36} is computed via CRT.
1403 @example
1404 (%i1) m: 7*9$
1406 (%i2) zn_power_table(m,9);
1407 [ 9 18 36 9 18 36 ]
1408 [ ]
1409 [ 18 9 36 18 9 36 ]
1410 [ ]
1411 [ 27 36 27 36 27 36 ]
1412 (%o2) [ ]
1413 [ 36 36 36 36 36 36 ]
1414 [ ]
1415 [ 45 9 27 18 54 36 ]
1416 [ ]
1417 [ 54 18 27 9 45 36 ]
1418 (%i3) zn_nth_root(27,3,m);
1419 (%o3) [27, 45, 54]
1420 (%i4) id7:1$ id63_9: chinese([id7,0],[7,9]);
1421 (%o5) 36
1422 @end example
1424 In the following RSA-like example, where the modulus @code{N} is squarefree,
1425 i.e. it splits into
1426 exclusively first power factors, every @code{x} from @code{0} to @code{N-1}
1427 is contained in a multiplication subgroup.
1429 The process of decryption needs the @code{e}-th root.
1430 @code{e} is coprime to @code{totient(N)} and therefore the @code{e}-th root is unique.
1431 In this case @code{zn_nth_root} effectively performs CRT-RSA.
1432 (Please note that @code{flatten} removes braces but no solutions.)
1434 @example
1435 (%i1) [p,q,e]: [5,7,17]$ N: p*q$
1437 (%i3) xs: makelist(x,x,0,N-1)$
1439 (%i4) ys: map(lambda([x],power_mod(x,e,N)),xs)$
1441 (%i5) zs: flatten(map(lambda([y], zn_nth_root(y,e,N)), ys))$
1443 (%i6) is(zs = xs);
1444 (%o6) true
1445 @end example
1447 In the following example the factorization of the modulus is known
1448 and passed as the fourth argument.
1450 @example
1451 (%i1) p: 2^107-1$ q: 2^127-1$ N: p*q$
1453 (%i4) ibase: obase: 16$
1455 (%i5) msg: 11223344556677889900aabbccddeeff$
1457 (%i6) enc: power_mod(msg, 10001, N);
1458 (%o6) 1a8db7892ae588bdc2be25dd5107a425001fe9c82161abc673241c8b383
1459 (%i7) zn_nth_root(enc, 10001, N, [[p,1],[q,1]]);
1460 (%o7) [11223344556677889900aabbccddeeff]
1461 @end example
1463 @opencatbox{Categories:}
1464 @category{Number theory}
1465 @closecatbox
1466 @end deffn
1468 @c -----------------------------------------------------------------------------
1469 @anchor{zn_order}
1470 @deffn {Function} zn_order @
1471 @fname{zn_order} (@var{x}, @var{n}) @
1472 @fname{zn_order} (@var{x}, @var{n}, [[@var{p1}, @var{e1}], @dots{}, [@var{pk}, @var{ek}]])
1474 Returns the order of @var{x} if it is a unit of the finite group (Z/@var{n}Z)*
1475 or returns @code{false}. @var{x} is a unit modulo @var{n} if it is coprime to @var{n}.
1477 The applied algorithm needs a prime factorization of @code{totient(n)}. This factorization
1478 might be time consuming in some cases and it can be useful to factor first
1479 and then to pass the list of factors to @code{zn_log} as the third argument.
1480 The list must be of the same form as the list returned by @code{ifactors(totient(n))}
1481 using the default option @code{factors_only : false}.
1483 See also @mrefcomma{zn_primroot} @mrefcomma{ifactors} @mrefdot{totient}
1485 Examples:
1487 @code{zn_order} computes the order of the unit @var{x} in (Z/@var{n}Z)*.
1489 @example
1490 (%i1) n: 22$
1491 (%i2) g: zn_primroot(n);
1492 (%o2) 7
1493 (%i3) units_22: sublist(makelist(i,i,1,21), lambda([x], gcd(x,n)=1));
1494 (%o3) [1, 3, 5, 7, 9, 13, 15, 17, 19, 21]
1495 (%i4) (ord_7: zn_order(7, n)) = totient(n);
1496 (%o4) 10 = 10
1497 (%i5) powers_7: makelist(power_mod(g,i,n), i,0,ord_7 - 1);
1498 (%o5) [1, 7, 5, 13, 3, 21, 15, 17, 9, 19]
1499 (%i6) map(lambda([x], zn_order(x, n)), powers_7);
1500 (%o6) [1, 10, 5, 10, 5, 2, 5, 10, 5, 10]
1501 (%i7) map(lambda([x], ord_7/gcd(x,ord_7)), makelist(i,i,0,ord_7-1));
1502 (%o7) [1, 10, 5, 10, 5, 2, 5, 10, 5, 10]
1503 (%i8) totient(totient(n));
1504 (%o8) 4
1505 @end example
1507 The optional third argument must be of the same form as the list returned by
1508 @code{ifactors(totient(n))}.
1510 @example
1511 (%i1) (p : 2^142 + 217, primep(p));
1512 (%o1) true
1513 (%i2) ifs: ifactors( totient(p) )$
1514 (%i3) g: zn_primroot(p, ifs);
1515 (%o3) 3
1516 (%i4) is( (ord_3 : zn_order(g, p, ifs)) = totient(p) );
1517 (%o4) true
1518 (%i5) map(lambda([x], ord_3/zn_order(x,p,ifs)), makelist(i,i,2,15));
1519 (%o5) [22, 1, 44, 10, 5, 2, 22, 2, 8, 2, 1, 1, 20, 1]
1520 @end example
1522 @opencatbox{Categories:}
1523 @category{Number theory}
1524 @closecatbox
1525 @end deffn
1527 @c -----------------------------------------------------------------------------
1528 @anchor{zn_power_table}
1529 @deffn {Function} zn_power_table @
1530 @fname{zn_power_table} (@var{n}) @
1531 @fname{zn_power_table} (@var{n}, @var{gcd}) @
1532 @fname{zn_power_table} (@var{n}, @var{gcd}, @var{max_exp})
1534 Without any optional argument @code{zn_power_table(n)}
1535 shows a power table of all elements in (Z/@var{n}Z)*
1536 which are all residue classes coprime to @var{n}.
1537 The exponent loops from @code{1} to the greatest characteristic factor of
1538 @code{totient(n)} (also known as Carmichael function or Carmichael lambda)
1539 and the table ends with a column of ones on the right side.
1541 The optional second argument @var{gcd} allows to select powers of a specific
1542 subset of (Z/@var{n}Z). If @var{gcd} is an integer, powers of all residue
1543 classes @code{x} with @code{gcd(x,n) = }@var{gcd} are returned,
1544 i.e. the default value for @var{gcd} is @code{1}.
1545 If @var{gcd} is set to @code{all}, the table contains powers of all elements
1546 in (Z/@var{n}Z).
1548 If the optional third argument @var{max_exp} is given, the exponent loops from
1549 @code{1} to @var{max_exp}.
1551 See also @mrefcomma{zn_add_table} @mrefdot{zn_mult_table}
1553 Examples:
1555 The default which is @var{gcd}@code{ = 1} allows to demonstrate and study basic
1556 theorems of e.g. Fermat and Euler.
1558 The argument @var{gcd} allows to select subsets of (Z/@var{n}Z) and to study
1559 multiplication subgroups and isomorphisms.
1560 E.g. the groups @code{G10} and @code{G10_2} are under multiplication both
1561 isomorphic to @code{G5}. @code{1} is the identity in @code{G5}.
1562 So are @code{1} resp. @code{6} the identities in @code{G10} resp. @code{G10_2}.
1563 There are corresponding mappings for primitive roots, n-th roots, etc..
1565 @example
1566 (%i1) zn_power_table(10);
1567 [ 1 1 1 1 ]
1568 [ ]
1569 [ 3 9 7 1 ]
1570 (%o1) [ ]
1571 [ 7 9 3 1 ]
1572 [ ]
1573 [ 9 1 9 1 ]
1574 (%i2) zn_power_table(10,2);
1575 [ 2 4 8 6 ]
1576 [ ]
1577 [ 4 6 4 6 ]
1578 (%o2) [ ]
1579 [ 6 6 6 6 ]
1580 [ ]
1581 [ 8 4 2 6 ]
1582 (%i3) zn_power_table(10,5);
1583 (%o3) [ 5 5 5 5 ]
1584 (%i4) zn_power_table(10,10);
1585 (%o4) [ 0 0 0 0 ]
1586 (%i5) G5: [1,2,3,4];
1587 (%o6) [1, 2, 3, 4]
1588 (%i6) G10_2: map(lambda([x], chinese([0,x],[2,5])), G5);
1589 (%o6) [6, 2, 8, 4]
1590 (%i7) G10: map(lambda([x], power_mod(3, zn_log(x,2,5), 10)), G5);
1591 (%o7) [1, 3, 7, 9]
1592 @end example
1594 If @var{gcd} is set to @code{all}, the table contains powers of all elements
1595 in (Z/@var{n}Z).
1597 The third argument @var{max_exp} allows to set the highest exponent.
1598 The following table shows a very small example of RSA.
1600 @example
1601 (%i1) N:2*5$ phi:totient(N)$ e:7$ d:inv_mod(e,phi)$
1603 (%i5) zn_power_table(N, all, e*d);
1604 [ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
1605 [ ]
1606 [ 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ]
1607 [ ]
1608 [ 2 4 8 6 2 4 8 6 2 4 8 6 2 4 8 6 2 4 8 6 2 ]
1609 [ ]
1610 [ 3 9 7 1 3 9 7 1 3 9 7 1 3 9 7 1 3 9 7 1 3 ]
1611 [ ]
1612 [ 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 ]
1613 (%o5) [ ]
1614 [ 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 ]
1615 [ ]
1616 [ 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 ]
1617 [ ]
1618 [ 7 9 3 1 7 9 3 1 7 9 3 1 7 9 3 1 7 9 3 1 7 ]
1619 [ ]
1620 [ 8 4 2 6 8 4 2 6 8 4 2 6 8 4 2 6 8 4 2 6 8 ]
1621 [ ]
1622 [ 9 1 9 1 9 1 9 1 9 1 9 1 9 1 9 1 9 1 9 1 9 ]
1623 @end example
1625 @opencatbox{Categories:}
1626 @category{Number theory}
1627 @closecatbox
1628 @end deffn
1630 @c -----------------------------------------------------------------------------
1631 @anchor{zn_primroot}
1632 @deffn {Function} zn_primroot @
1633 @fname{zn_primroot} (@var{n}) @
1634 @fname{zn_primroot} (@var{n}, [[@var{p1}, @var{e1}], @dots{}, [@var{pk}, @var{ek}]])
1636 If the multiplicative group (Z/@var{n}Z)* is cyclic, @code{zn_primroot} computes the
1637 smallest primitive root modulo @var{n}. (Z/@var{n}Z)* is cyclic if @var{n} is equal to
1638 @code{2}, @code{4}, @code{p^k} or @code{2*p^k}, where @code{p} is prime and
1639 greater than @code{2} and @code{k} is a natural number. @code{zn_primroot}
1640 performs an according pretest if the option variable @mref{zn_primroot_pretest}
1641 (default: @code{false}) is set to @code{true}. In any case the computation is limited
1642 by the upper bound @mrefdot{zn_primroot_limit}
1644 If (Z/@var{n}Z)* is not cyclic or if there is no primitive root up to
1645 @code{zn_primroot_limit}, @code{zn_primroot} returns @code{false}.
1647 The applied algorithm needs a prime factorization of @code{totient(n)}. This factorization
1648 might be time consuming in some cases and it can be useful to factor first
1649 and then to pass the list of factors to @code{zn_log} as an additional argument.
1650 The list must be of the same form as the list returned by @code{ifactors(totient(n))}
1651 using the default option @code{factors_only : false}.
1653 See also @mrefcomma{zn_primroot_p} @mrefcomma{zn_order} @mrefcomma{ifactors} @mrefdot{totient}
1655 Examples:
1657 @code{zn_primroot} computes the smallest primitive root modulo @var{n} or returns
1658 @code{false}.
1660 @example
1661 (%i1) n : 14$
1662 (%i2) g : zn_primroot(n);
1663 (%o2) 3
1664 (%i3) zn_order(g, n) = totient(n);
1665 (%o3) 6 = 6
1666 (%i4) n : 15$
1667 (%i5) zn_primroot(n);
1668 (%o5) false
1669 @end example
1671 The optional second argument must be of the same form as the list returned by
1672 @code{ifactors(totient(n))}.
1674 @example
1675 (%i1) (p : 2^142 + 217, primep(p));
1676 (%o1) true
1677 (%i2) ifs : ifactors( totient(p) )$
1678 (%i3) g : zn_primroot(p, ifs);
1679 (%o3) 3
1680 (%i4) [time(%o2), time(%o3)];
1681 (%o4) [[15.556972], [0.004]]
1682 (%i5) is(zn_order(g, p, ifs) = p - 1);
1683 (%o5) true
1684 (%i6) n : 2^142 + 216$
1685 (%i7) ifs : ifactors(totient(n))$
1686 (%i8) zn_primroot(n, ifs),
1687 zn_primroot_limit : 200, zn_primroot_verbose : true;
1688 `zn_primroot' stopped at zn_primroot_limit = 200
1689 (%o8) false
1690 @end example
1692 @opencatbox{Categories:}
1693 @category{Number theory}
1694 @closecatbox
1695 @end deffn
1697 @c -----------------------------------------------------------------------------
1698 @anchor{zn_primroot_limit}
1699 @defvr {Option variable} zn_primroot_limit
1700 Default value: @code{1000}
1702 If @mref{zn_primroot} cannot find a primitive root, it stops at this upper bound.
1703 If the option variable @mref{zn_primroot_verbose} (default: @code{false}) is
1704 set to @code{true}, a message will be printed when @code{zn_primroot_limit} is reached.
1706 @opencatbox{Categories:}
1707 @category{Number theory}
1708 @closecatbox
1709 @end defvr
1711 @c -----------------------------------------------------------------------------
1712 @anchor{zn_primroot_p}
1713 @deffn {Function} zn_primroot_p @
1714 @fname{zn_primroot_p} (@var{x}, @var{n}) @
1715 @fname{zn_primroot_p} (@var{x}, @var{n}, [[@var{p1}, @var{e1}], @dots{}, [@var{pk}, @var{ek}]])
1717 Checks whether @var{x} is a primitive root in the multiplicative group (Z/@var{n}Z)*.
1719 The applied algorithm needs a prime factorization of @code{totient(n)}. This factorization
1720 might be time consuming and in case @code{zn_primroot_p} will be consecutively
1721 applied to a list of candidates it can be useful to factor first and then to
1722 pass the list of factors to @code{zn_log} as a third argument.
1723 The list must be of the same form as the list returned by @code{ifactors(totient(n))}
1724 using the default option @code{factors_only : false}.
1726 See also @mrefcomma{zn_primroot} @mrefcomma{zn_order} @mrefcomma{ifactors} @mrefdot{totient}
1728 Examples:
1730 @code{zn_primroot_p} as a predicate function.
1732 @example
1733 (%i1) n : 14$
1734 (%i2) units_14 : sublist(makelist(i,i,1,13), lambda([i], gcd(i, n) = 1));
1735 (%o2) [1, 3, 5, 9, 11, 13]
1736 (%i3) zn_primroot_p(13, n);
1737 (%o3) false
1738 (%i4) sublist(units_14, lambda([x], zn_primroot_p(x, n)));
1739 (%o4) [3, 5]
1740 (%i5) map(lambda([x], zn_order(x, n)), units_14);
1741 (%o5) [1, 6, 6, 3, 3, 2]
1742 @end example
1744 The optional third argument must be of the same form as the list returned by
1745 @code{ifactors(totient(n))}.
1747 @example
1748 (%i1) (p: 2^142 + 217, primep(p));
1749 (%o1) true
1750 (%i2) ifs: ifactors( totient(p) )$
1751 (%i3) sublist(makelist(i,i,1,50), lambda([x], zn_primroot_p(x,p,ifs)));
1752 (%o3) [3, 12, 13, 15, 21, 24, 26, 27, 29, 33, 38, 42, 48]
1753 (%i4) [time(%o2), time(%o3)];
1754 (%o4) [[7.748484], [0.036002]]
1755 @end example
1757 @opencatbox{Categories:}
1758 @category{Predicate functions}
1759 @category{Number theory}
1760 @closecatbox
1761 @end deffn
1763 @c -----------------------------------------------------------------------------
1764 @anchor{zn_primroot_pretest}
1765 @defvr {Option variable} zn_primroot_pretest
1766 Default value: @code{false}
1768 The multiplicative group (Z/@var{n}Z)* is cyclic if @var{n} is equal to
1769 @code{2}, @code{4}, @code{p^k} or @code{2*p^k}, where @code{p} is prime and
1770 greater than @code{2} and @code{k} is a natural number.
1772 @code{zn_primroot_pretest} controls whether @mref{zn_primroot} will check
1773 if one of these cases occur before it computes the smallest primitive root.
1774 Only if @code{zn_primroot_pretest} is set to @code{true} this pretest will be
1775 performed.
1777 @opencatbox{Categories:}
1778 @category{Number theory}
1779 @closecatbox
1780 @end defvr
1782 @c -----------------------------------------------------------------------------
1783 @anchor{zn_primroot_verbose}
1784 @defvr {Option variable} zn_primroot_verbose
1785 Default value: @code{false}
1787 Controls whether @mref{zn_primroot} prints a message when reaching
1788 @mrefdot{zn_primroot_limit}
1790 @opencatbox{Categories:}
1791 @category{Number theory}
1792 @closecatbox
1793 @end defvr