Annotation of OpenXM_contrib/gmp/tests/mpz/t-gcd.c, Revision 1.1.1.1
1.1 ohara 1: /* Test mpz_gcd, mpz_gcdext, mpz_mul, mpz_tdiv_r, mpz_add, mpz_cmp,
2: mpz_cmp_ui, mpz_init_set, mpz_set, mpz_clear.
3:
4: Copyright 1991, 1993, 1994, 1996, 1997, 2000, 2001, 2002 Free Software
5: Foundation, Inc.
6:
7: This file is part of the GNU MP Library.
8:
9: The GNU MP Library is free software; you can redistribute it and/or modify
10: it under the terms of the GNU Lesser General Public License as published by
11: the Free Software Foundation; either version 2.1 of the License, or (at your
12: option) any later version.
13:
14: The GNU MP Library is distributed in the hope that it will be useful, but
15: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17: License for more details.
18:
19: You should have received a copy of the GNU Lesser General Public License
20: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
21: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
22: MA 02111-1307, USA. */
23:
24: #include <stdio.h>
25: #include <stdlib.h>
26:
27: #include "gmp.h"
28: #include "gmp-impl.h"
29: #include "tests.h"
30:
31: void dump_abort _PROTO ((int, mpz_t, mpz_t));
32: void debug_mp _PROTO ((mpz_t, int));
33:
34: int
35: main (int argc, char **argv)
36: {
37: mpz_t op1, op2, x;
38: mpz_t gcd, gcd2, s, t;
39: mpz_t temp1, temp2;
40: mp_size_t op1_size, op2_size, x_size;
41: int i;
42: int reps = 2000;
43: gmp_randstate_ptr rands;
44: mpz_t bs;
45: unsigned long bsi, size_range;
46:
47: tests_start ();
48: rands = RANDS;
49:
50: mpz_init (bs);
51:
52: if (argc == 2)
53: reps = atoi (argv[1]);
54:
55: mpz_init (op1);
56: mpz_init (op2);
57: mpz_init (x);
58: mpz_init (gcd);
59: mpz_init (gcd2);
60: mpz_init (temp1);
61: mpz_init (temp2);
62: mpz_init (s);
63: mpz_init (t);
64:
65: for (i = 0; i < reps; i++)
66: {
67: mpz_urandomb (bs, rands, 32);
68: size_range = mpz_get_ui (bs) % 12 + 2; /* 0..8191 bit operands */
69:
70: mpz_urandomb (bs, rands, size_range);
71: op1_size = mpz_get_ui (bs);
72: mpz_rrandomb (op1, rands, op1_size);
73:
74: mpz_urandomb (bs, rands, size_range);
75: op2_size = mpz_get_ui (bs);
76: mpz_rrandomb (op2, rands, op2_size);
77:
78: mpz_urandomb (bs, rands, size_range);
79: x_size = mpz_get_ui (bs);
80: mpz_rrandomb (x, rands, x_size);
81:
82: mpz_urandomb (bs, rands, 2);
83: bsi = mpz_get_ui (bs);
84: if ((bsi & 1) != 0)
85: mpz_neg (op1, op1);
86: if ((bsi & 2) != 0)
87: mpz_neg (op2, op2);
88:
89: /* printf ("%ld %ld\n", SIZ (op1), SIZ (op2)); */
90:
91: mpz_mul (op1, op1, x);
92: mpz_mul (op2, op2, x);
93:
94: mpz_gcd (gcd, op1, op2);
95: /* We know GCD will be at least X, since we multiplied both operands
96: with it. */
97: if (mpz_cmp (gcd, x) < 0 && mpz_sgn (op1) != 0 && mpz_sgn (op2) != 0)
98: dump_abort (i, op1, op2);
99:
100: if (mpz_fits_ulong_p (op2))
101: {
102: mpz_gcd_ui (gcd2, op1, mpz_get_ui (op2));
103: if (mpz_cmp (gcd, gcd2))
104: dump_abort (i, op1, op2);
105: }
106:
107: mpz_gcdext (gcd2, s, t, op1, op2);
108: if (mpz_cmp (gcd, gcd2))
109: dump_abort (i, op1, op2);
110:
111: mpz_gcdext (gcd2, s, NULL, op1, op2);
112: if (mpz_cmp (gcd, gcd2))
113: dump_abort (i, op1, op2);
114:
115: mpz_mul (temp1, s, op1);
116: mpz_mul (temp2, t, op2);
117: mpz_add (gcd2, temp1, temp2);
118: if (mpz_cmp (gcd, gcd2))
119: dump_abort (i, op1, op2);
120: }
121:
122: mpz_clear (bs);
123: mpz_clear (op1);
124: mpz_clear (op2);
125: mpz_clear (x);
126: mpz_clear (gcd);
127: mpz_clear (gcd2);
128: mpz_clear (temp1);
129: mpz_clear (temp2);
130: mpz_clear (s);
131: mpz_clear (t);
132:
133: tests_end ();
134: exit (0);
135: }
136:
137: void
138: dump_abort (int testnr, mpz_t op1, mpz_t op2)
139: {
140: fprintf (stderr, "ERROR in test %d\n", testnr);
141: fprintf (stderr, "op1 = "); debug_mp (op1, -16);
142: fprintf (stderr, "op2 = "); debug_mp (op2, -16);
143: abort();
144: }
145:
146: void
147: debug_mp (mpz_t x, int base)
148: {
149: mpz_out_str (stderr, base, x); fputc ('\n', stderr);
150: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>