1 /* Copyright (c) 1984, 1986, 1987, 1988, 1989 AT&T */
2 /* All Rights Reserved */
3
4
5 /*
6 * Copyright (c) 1980 Regents of the University of California.
7 * All rights reserved. The Berkeley software License Agreement
8 * specifies the terms and conditions for redistribution.
9 */
10 /* Portions Copyright(c) 1988, Sun Microsystems Inc. */
11 /* All Rights Reserved */
12
13 /*
14 * Copyright (c) 1997, by Sun Microsystems, Inc.
15 * All rights reserved.
16 */
17
18 #ident "%Z%%M% %I% %E% SMI" /* SVr4.0 1.1 */
19
20 /* LINTLIBRARY */
21
22 #include <mp.h>
23 #include "libmp.h"
24 #include <sys/types.h>
25
26 void
mp_gcd(MINT * a,MINT * b,MINT * c)27 mp_gcd(MINT *a, MINT *b, MINT *c)
28 {
29 MINT x, y, z, w;
30
31 x.len = y.len = z.len = w.len = 0;
32 _mp_move(a, &x);
33 _mp_move(b, &y);
34 while (y.len != 0) {
35 mp_mdiv(&x, &y, &w, &z);
36 _mp_move(&y, &x);
37 _mp_move(&z, &y);
38 }
39 _mp_move(&x, c);
40 _mp_xfree(&x);
41 _mp_xfree(&y);
42 _mp_xfree(&z);
43 _mp_xfree(&w);
44 }
45
46 void
mp_invert(MINT * x1,MINT * x0,MINT * c)47 mp_invert(MINT *x1, MINT *x0, MINT *c)
48 {
49 MINT u2, u3;
50 MINT v2, v3;
51 MINT zero;
52 MINT q, r;
53 MINT t;
54 MINT x0_prime;
55 static MINT *one = NULL;
56
57 /*
58 * Minimize calls to allocators. Don't use pointers for local
59 * variables, for the one "initialized" multiple precision
60 * variable, do it just once.
61 */
62 if (one == NULL)
63 one = mp_itom(1);
64
65 zero.len = q.len = r.len = t.len = 0;
66
67 x0_prime.len = u2.len = u3.len = 0;
68 _mp_move(x0, &u3);
69 _mp_move(x0, &x0_prime);
70
71 v2.len = v3.len = 0;
72 _mp_move(one, &v2);
73 _mp_move(x1, &v3);
74
75 while (mp_mcmp(&v3, &zero) != 0) {
76 /* invariant: x0*u1 + x1*u2 = u3 */
77 /* invariant: x0*v1 + x2*v2 = v3 */
78 /* invariant: x(n+1) = x(n-1) % x(n) */
79 mp_mdiv(&u3, &v3, &q, &r);
80 _mp_move(&v3, &u3);
81 _mp_move(&r, &v3);
82
83 mp_mult(&q, &v2, &t);
84 mp_msub(&u2, &t, &t);
85 _mp_move(&v2, &u2);
86 _mp_move(&t, &v2);
87 }
88 /* now x0*u1 + x1*u2 == 1, therefore, (u2*x1) % x0 == 1 */
89 _mp_move(&u2, c);
90 if (mp_mcmp(c, &zero) < 0) {
91 mp_madd(&x0_prime, c, c);
92 }
93 _mp_xfree(&zero);
94 _mp_xfree(&v2);
95 _mp_xfree(&v3);
96 _mp_xfree(&u2);
97 _mp_xfree(&u3);
98 _mp_xfree(&q);
99 _mp_xfree(&r);
100 _mp_xfree(&t);
101 }
102