ceb3f0c0 |
1 | /* -*-c-*- |
2 | * |
bc985cef |
3 | * $Id: gfreduce.c,v 1.3 2004/03/23 15:19:32 mdw Exp $ |
ceb3f0c0 |
4 | * |
5 | * Efficient reduction modulo sparse binary polynomials |
6 | * |
7 | * (c) 2004 Straylight/Edgeware |
8 | */ |
9 | |
10 | /*----- Licensing notice --------------------------------------------------* |
11 | * |
12 | * This file is part of Catacomb. |
13 | * |
14 | * Catacomb is free software; you can redistribute it and/or modify |
15 | * it under the terms of the GNU Library General Public License as |
16 | * published by the Free Software Foundation; either version 2 of the |
17 | * License, or (at your option) any later version. |
18 | * |
19 | * Catacomb is distributed in the hope that it will be useful, |
20 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
21 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
22 | * GNU Library General Public License for more details. |
23 | * |
24 | * You should have received a copy of the GNU Library General Public |
25 | * License along with Catacomb; if not, write to the Free |
26 | * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, |
27 | * MA 02111-1307, USA. |
28 | */ |
29 | |
30 | /*----- Revision history --------------------------------------------------* |
31 | * |
32 | * $Log: gfreduce.c,v $ |
bc985cef |
33 | * Revision 1.3 2004/03/23 15:19:32 mdw |
34 | * Test elliptic curves more thoroughly. |
35 | * |
c3caa2fa |
36 | * Revision 1.2 2004/03/21 22:52:06 mdw |
37 | * Merge and close elliptic curve branch. |
38 | * |
ceb3f0c0 |
39 | * Revision 1.1.2.1 2004/03/21 22:39:46 mdw |
40 | * Elliptic curves on binary fields work. |
41 | * |
42 | */ |
43 | |
44 | /*----- Header files ------------------------------------------------------*/ |
45 | |
46 | #include <mLib/alloc.h> |
47 | #include <mLib/darray.h> |
48 | #include <mLib/macros.h> |
49 | |
50 | #include "gf.h" |
51 | #include "gfreduce.h" |
52 | #include "gfreduce-exp.h" |
53 | #include "fibrand.h" |
54 | #include "mprand.h" |
55 | |
56 | /*----- Data structures ---------------------------------------------------*/ |
57 | |
58 | DA_DECL(instr_v, gfreduce_instr); |
59 | |
60 | /*----- Main code ---------------------------------------------------------*/ |
61 | |
62 | /* --- What's going on here? --- * |
63 | * |
64 | * Let's face it, @gfx_div@ sucks. It works (I hope), but it's not in any |
65 | * sense fast. Here, we do efficient reduction modulo sparse polynomials. |
66 | * |
67 | * Suppose we have a polynomial @X@ we're trying to reduce mod @P@. If we |
68 | * take the topmost nonzero word of @X@, call it @w@, then we can eliminate |
69 | * it by subtracting off @w P x^{k}@ for an appropriate value of @k@. The |
70 | * trick is in observing that if @P@ is sparse we can do this multiplication |
71 | * and subtraction efficiently, just by XORing appropriate shifts of @w@ into |
72 | * @X@. |
73 | * |
74 | * The first tricky bit is in working out when to stop. I'll use eight-bit |
75 | * words to demonstrate what I'm talking about. |
76 | * |
77 | * xxxxxx xxxxxxxx xxxxxxxx xxxxxxxx xxxxxxxx xxxxxxxx |
78 | * 001ppppp pppppppp pppppppp pppppppp |
79 | * |<rp>| |
80 | * |<------------ bp ------------->| |
81 | * |<------------ nw --------------->| |
82 | * |
83 | * The trick of taking whole words off of @X@ stops working when there are |
84 | * only @nw@ words left. Then we have to mask off the bottom bits of @w@ |
85 | * before continuing. |
86 | */ |
87 | |
88 | /* --- @gfreduce_create@ --- * |
89 | * |
90 | * Arguments: @gfreduce *r@ = structure to fill in |
91 | * @mp *x@ = a (hopefully sparse) polynomial |
92 | * |
93 | * Returns: --- |
94 | * |
95 | * Use: Initializes a context structure for reduction. |
96 | */ |
97 | |
98 | void gfreduce_create(gfreduce *r, mp *p) |
99 | { |
100 | instr_v iv = DA_INIT; |
101 | unsigned long d, dw; |
102 | mpscan sc; |
103 | unsigned long i; |
104 | gfreduce_instr *ip; |
105 | unsigned f = 0; |
106 | size_t w, ww, wi, wl, ll; |
107 | |
108 | /* --- Sort out the easy stuff --- */ |
109 | |
110 | d = mp_bits(p); assert(d); d--; |
111 | r->lim = d/MPW_BITS; |
112 | dw = d%MPW_BITS; |
113 | if (!dw) |
114 | r->mask = 0; |
115 | else { |
116 | r->mask = MPW(((mpw)-1) << dw); |
117 | r->lim++; |
118 | } |
119 | r->p = mp_copy(p); |
120 | |
121 | /* --- Stash a new instruction --- */ |
122 | |
123 | #define INSTR(op_, arg_) do { \ |
124 | DA_ENSURE(&iv, 1); \ |
125 | DA(&iv)[DA_LEN(&iv)].op = (op_); \ |
126 | DA(&iv)[DA_LEN(&iv)].arg = (arg_); \ |
127 | DA_EXTEND(&iv, 1); \ |
128 | } while (0) |
129 | |
130 | #define f_lsr 1u |
131 | |
132 | w = (d + MPW_BITS - 1)/MPW_BITS; |
133 | INSTR(GFRI_LOAD, w); |
134 | wi = DA_LEN(&iv); |
135 | f = 0; |
136 | ll = 0; |
137 | for (i = 0, mp_scan(&sc, p); mp_step(&sc) && i < d; i++) { |
138 | if (!mp_bit(&sc)) |
139 | continue; |
140 | ww = (d - i + MPW_BITS - 1)/MPW_BITS; |
141 | if (ww != w) { |
142 | wl = DA_LEN(&iv); |
143 | INSTR(GFRI_STORE, w); |
144 | if (!ll) |
145 | ll = DA_LEN(&iv); |
146 | if (!(f & f_lsr)) |
147 | INSTR(GFRI_LOAD, ww); |
148 | else { |
149 | INSTR(GFRI_LOAD, w - 1); |
150 | for (; wi < wl; wi++) { |
151 | ip = &DA(&iv)[wi]; |
152 | assert(ip->op == GFRI_LSL); |
153 | if (ip->arg) |
154 | INSTR(GFRI_LSR, MPW_BITS - ip->arg); |
155 | } |
156 | if (w - 1 != ww) { |
157 | INSTR(GFRI_STORE, w - 1); |
158 | INSTR(GFRI_LOAD, ww); |
159 | } |
160 | f &= ~f_lsr; |
161 | } |
162 | w = ww; |
163 | wi = DA_LEN(&iv); |
164 | } |
165 | INSTR(GFRI_LSL, (i - d)%MPW_BITS); |
166 | if ((i - d)%MPW_BITS) |
167 | f |= f_lsr; |
168 | } |
169 | wl = DA_LEN(&iv); |
170 | INSTR(GFRI_STORE, w); |
171 | if (!ll) |
172 | ll = DA_LEN(&iv); |
173 | if (f & f_lsr) { |
174 | INSTR(GFRI_LOAD, w - 1); |
175 | for (; wi < wl; wi++) { |
176 | ip = &DA(&iv)[wi]; |
177 | assert(ip->op == GFRI_LSL); |
178 | if (ip->arg) |
179 | INSTR(GFRI_LSR, MPW_BITS - ip->arg); |
180 | } |
181 | INSTR(GFRI_STORE, w - 1); |
182 | } |
183 | |
184 | #undef INSTR |
185 | |
186 | r->in = DA_LEN(&iv); |
187 | r->iv = xmalloc(r->in * sizeof(gfreduce_instr)); |
188 | r->liv = r->iv + ll; |
189 | memcpy(r->iv, DA(&iv), r->in * sizeof(gfreduce_instr)); |
190 | DA_DESTROY(&iv); |
191 | } |
192 | |
193 | /* --- @gfreduce_destroy@ --- * |
194 | * |
195 | * Arguments: @gfreduce *r@ = structure to free |
196 | * |
197 | * Returns: --- |
198 | * |
199 | * Use: Reclaims the resources from a reduction context. |
200 | */ |
201 | |
202 | void gfreduce_destroy(gfreduce *r) |
203 | { |
204 | mp_drop(r->p); |
205 | xfree(r->iv); |
206 | } |
207 | |
208 | /* --- @gfreduce_dump@ --- * |
209 | * |
210 | * Arguments: @gfreduce *r@ = structure to dump |
211 | * @FILE *fp@ = file to dump on |
212 | * |
213 | * Returns: --- |
214 | * |
215 | * Use: Dumps a reduction context. |
216 | */ |
217 | |
218 | void gfreduce_dump(gfreduce *r, FILE *fp) |
219 | { |
220 | size_t i; |
221 | |
222 | fprintf(fp, "poly = "); mp_writefile(r->p, fp, 16); |
223 | fprintf(fp, "\n lim = %lu; mask = %lx\n", |
224 | (unsigned long)r->lim, (unsigned long)r->mask); |
225 | for (i = 0; i < r->in; i++) { |
226 | static const char *opname[] = { "load", "lsl", "lsr", "store" }; |
227 | assert(r->iv[i].op < N(opname)); |
228 | fprintf(fp, " %s %lu\n", |
229 | opname[r->iv[i].op], |
230 | (unsigned long)r->iv[i].arg); |
231 | } |
232 | } |
233 | |
234 | /* --- @gfreduce_do@ --- * |
235 | * |
236 | * Arguments: @gfreduce *r@ = reduction context |
237 | * @mp *d@ = destination |
238 | * @mp *x@ = source |
239 | * |
240 | * Returns: Destination, @x@ reduced modulo the reduction poly. |
241 | */ |
242 | |
243 | static void run(const gfreduce_instr *i, const gfreduce_instr *il, |
244 | mpw *v, mpw z) |
245 | { |
246 | mpw w = 0; |
247 | |
248 | for (; i < il; i++) { |
249 | switch (i->op) { |
250 | case GFRI_LOAD: w = *(v - i->arg); break; |
251 | case GFRI_LSL: w ^= z << i->arg; break; |
252 | case GFRI_LSR: w ^= z >> i->arg; break; |
253 | case GFRI_STORE: *(v - i->arg) = MPW(w); break; |
254 | default: abort(); |
255 | } |
256 | } |
257 | } |
258 | |
259 | mp *gfreduce_do(gfreduce *r, mp *d, mp *x) |
260 | { |
261 | mpw *v, *vl; |
262 | const gfreduce_instr *il; |
263 | mpw z; |
264 | |
265 | /* --- Try to reuse the source's space --- */ |
266 | |
267 | MP_COPY(x); |
268 | if (d) MP_DROP(d); |
269 | MP_DEST(x, MP_LEN(x), x->f); |
270 | |
271 | /* --- Do the reduction --- */ |
272 | |
273 | il = r->iv + r->in; |
274 | if (MP_LEN(x) >= r->lim) { |
275 | v = x->v + r->lim; |
276 | vl = x->vl; |
277 | while (vl-- > v) { |
278 | while (*vl) { |
279 | z = *vl; |
280 | *vl = 0; |
281 | run(r->iv, il, vl, z); |
282 | } |
283 | } |
284 | if (r->mask) { |
285 | while (*vl & r->mask) { |
286 | z = *vl & r->mask; |
287 | *vl &= ~r->mask; |
288 | run(r->iv, il, vl, z); |
289 | } |
290 | } |
291 | } |
292 | |
293 | /* --- Done --- */ |
294 | |
295 | MP_SHRINK(x); |
296 | return (x); |
297 | } |
298 | |
299 | /* --- @gfreduce_sqrt@ --- * |
300 | * |
301 | * Arguments: @gfreduce *r@ = pointer to reduction context |
302 | * @mp *d@ = destination |
303 | * @mp *x@ = some polynomial |
304 | * |
305 | * Returns: The square root of @x@ modulo @r->p@, or null. |
306 | */ |
307 | |
308 | mp *gfreduce_sqrt(gfreduce *r, mp *d, mp *x) |
309 | { |
310 | mp *y = MP_COPY(x); |
311 | mp *z, *spare = MP_NEW; |
312 | unsigned long m = mp_bits(r->p) - 1; |
313 | unsigned long i; |
314 | |
315 | for (i = 0; i < m - 1; i++) { |
316 | mp *t = gf_sqr(spare, y); |
317 | spare = y; |
318 | y = gfreduce_do(r, t, t); |
319 | } |
320 | z = gf_sqr(spare, y); |
321 | z = gfreduce_do(r, z, z); |
322 | if (!MP_EQ(x, z)) { |
323 | mp_drop(y); |
324 | y = 0; |
325 | } |
326 | mp_drop(z); |
327 | mp_drop(d); |
328 | return (y); |
329 | } |
330 | |
331 | /* --- @gfreduce_trace@ --- * |
332 | * |
333 | * Arguments: @gfreduce *r@ = pointer to reduction context |
334 | * @mp *x@ = some polynomial |
335 | * |
336 | * Returns: The trace of @x@. (%$\Tr(x)=x + x^2 + \cdots + x^{2^{m-1}}$% |
337 | * if %$x \in \gf{2^m}$%). |
338 | */ |
339 | |
340 | int gfreduce_trace(gfreduce *r, mp *x) |
341 | { |
342 | mp *y = MP_COPY(x); |
343 | mp *spare = MP_NEW; |
344 | unsigned long m = mp_bits(r->p) - 1; |
345 | unsigned long i; |
346 | int rc; |
347 | |
348 | for (i = 0; i < m - 1; i++) { |
349 | mp *t = gf_sqr(spare, y); |
350 | spare = y; |
351 | y = gfreduce_do(r, t, t); |
352 | y = gf_add(y, y, x); |
353 | } |
354 | rc = !MP_ISZERO(y); |
355 | mp_drop(spare); |
356 | mp_drop(y); |
357 | return (rc); |
358 | } |
359 | |
360 | /* --- @gfreduce_halftrace@ --- * |
361 | * |
362 | * Arguments: @gfreduce *r@ = pointer to reduction context |
363 | * @mp *d@ = destination |
364 | * @mp *x@ = some polynomial |
365 | * |
366 | * Returns: The half-trace of @x@. |
367 | * (%$\HfTr(x)= x + x^{2^2} + \cdots + x^{2^{m-1}}$% |
368 | * if %$x \in \gf{2^m}$% with %$m$% odd). |
369 | */ |
370 | |
371 | mp *gfreduce_halftrace(gfreduce *r, mp *d, mp *x) |
372 | { |
373 | mp *y = MP_COPY(x); |
374 | mp *spare = MP_NEW; |
375 | unsigned long m = mp_bits(r->p) - 1; |
376 | unsigned long i; |
377 | |
378 | mp_drop(d); |
379 | for (i = 0; i < m - 1; i += 2) { |
380 | mp *t = gf_sqr(spare, y); |
381 | spare = y; |
382 | y = gfreduce_do(r, t, t); |
383 | t = gf_sqr(spare, y); |
384 | spare = y; |
385 | y = gfreduce_do(r, t, t); |
386 | y = gf_add(y, y, x); |
387 | } |
388 | mp_drop(spare); |
389 | return (y); |
390 | } |
391 | |
392 | /* --- @gfreduce_quadsolve@ --- * |
393 | * |
394 | * Arguments: @gfreduce *r@ = pointer to reduction context |
395 | * @mp *d@ = destination |
396 | * @mp *x@ = some polynomial |
397 | * |
398 | * Returns: A polynomial @y@ such that %$y^2 + y = x$%, or null. |
399 | */ |
400 | |
401 | mp *gfreduce_quadsolve(gfreduce *r, mp *d, mp *x) |
402 | { |
403 | unsigned long m = mp_bits(r->p) - 1; |
404 | mp *t; |
405 | |
406 | MP_COPY(x); |
407 | if (m & 1) |
408 | d = gfreduce_halftrace(r, d, x); |
409 | else { |
410 | mp *z, *w, *rho = MP_NEW; |
411 | mp *spare = MP_NEW; |
412 | grand *fr = fibrand_create(0); |
413 | unsigned long i; |
414 | |
415 | for (;;) { |
416 | rho = mprand(rho, m, fr, 0); |
417 | z = MP_ZERO; |
418 | w = MP_COPY(rho); |
419 | for (i = 0; i < m - 1; i++) { |
420 | t = gf_sqr(spare, z); spare = z; z = gfreduce_do(r, t, t); |
421 | t = gf_sqr(spare, w); spare = w; w = gfreduce_do(r, t, t); |
422 | t = gf_mul(spare, w, x); t = gfreduce_do(r, t, t); spare = t; |
423 | z = gf_add(z, z, t); |
424 | w = gf_add(w, w, rho); |
425 | } |
426 | if (!MP_ISZERO(w)) |
427 | break; |
428 | MP_DROP(z); |
429 | MP_DROP(w); |
430 | } |
431 | if (d) MP_DROP(d); |
432 | MP_DROP(w); |
433 | MP_DROP(spare); |
434 | MP_DROP(rho); |
435 | fr->ops->destroy(fr); |
436 | d = z; |
437 | } |
438 | |
439 | t = gf_sqr(MP_NEW, d); t = gfreduce_do(r, t, t); t = gf_add(t, t, d); |
440 | if (!MP_EQ(t, x)) { |
441 | MP_DROP(d); |
442 | d = 0; |
443 | } |
444 | MP_DROP(t); |
445 | MP_DROP(x); |
bc985cef |
446 | if (d) d->v[0] &= ~(mpw)1; |
ceb3f0c0 |
447 | return (d); |
448 | } |
449 | |
450 | /* --- @gfreduce_exp@ --- * |
451 | * |
452 | * Arguments: @gfreduce *gr@ = pointer to reduction context |
453 | * @mp *d@ = fake destination |
454 | * @mp *a@ = base |
455 | * @mp *e@ = exponent |
456 | * |
457 | * Returns: Result, %$a^e \bmod m$%. |
458 | */ |
459 | |
460 | mp *gfreduce_exp(gfreduce *gr, mp *d, mp *a, mp *e) |
461 | { |
462 | mp *x = MP_ONE; |
463 | mp *spare = (e->f & MP_BURN) ? MP_NEWSEC : MP_NEW; |
464 | |
465 | MP_SHRINK(e); |
466 | if (!MP_LEN(e)) |
467 | ; |
468 | else if (MP_LEN(e) < EXP_THRESH) |
469 | EXP_SIMPLE(x, a, e); |
470 | else |
471 | EXP_WINDOW(x, a, e); |
472 | mp_drop(d); |
473 | mp_drop(spare); |
474 | return (x); |
475 | } |
476 | |
477 | /*----- Test rig ----------------------------------------------------------*/ |
478 | |
479 | #ifdef TEST_RIG |
480 | |
481 | #define MP(x) mp_readstring(MP_NEW, #x, 0, 0) |
482 | |
483 | static int vreduce(dstr *v) |
484 | { |
485 | mp *d = *(mp **)v[0].buf; |
486 | mp *n = *(mp **)v[1].buf; |
487 | mp *r = *(mp **)v[2].buf; |
488 | mp *c; |
489 | int ok = 1; |
490 | gfreduce rr; |
491 | |
492 | gfreduce_create(&rr, d); |
493 | c = gfreduce_do(&rr, MP_NEW, n); |
494 | if (!MP_EQ(c, r)) { |
495 | fprintf(stderr, "\n*** reduction failed\n*** "); |
496 | gfreduce_dump(&rr, stderr); |
497 | fprintf(stderr, "\n*** n = "); mp_writefile(n, stderr, 16); |
498 | fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16); |
499 | fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16); |
500 | fprintf(stderr, "\n"); |
501 | ok = 0; |
502 | } |
503 | gfreduce_destroy(&rr); |
504 | mp_drop(n); mp_drop(d); mp_drop(r); mp_drop(c); |
505 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
506 | return (ok); |
507 | } |
508 | |
509 | static int vmodexp(dstr *v) |
510 | { |
511 | mp *p = *(mp **)v[0].buf; |
512 | mp *g = *(mp **)v[1].buf; |
513 | mp *x = *(mp **)v[2].buf; |
514 | mp *r = *(mp **)v[3].buf; |
515 | mp *c; |
516 | int ok = 1; |
517 | gfreduce rr; |
518 | |
519 | gfreduce_create(&rr, p); |
520 | c = gfreduce_exp(&rr, MP_NEW, g, x); |
521 | if (!MP_EQ(c, r)) { |
522 | fprintf(stderr, "\n*** modexp failed\n*** "); |
523 | fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16); |
524 | fprintf(stderr, "\n*** g = "); mp_writefile(g, stderr, 16); |
525 | fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16); |
526 | fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16); |
527 | fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16); |
528 | fprintf(stderr, "\n"); |
529 | ok = 0; |
530 | } |
531 | gfreduce_destroy(&rr); |
532 | mp_drop(p); mp_drop(g); mp_drop(r); mp_drop(x); mp_drop(c); |
533 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
534 | return (ok); |
535 | } |
536 | |
537 | static int vsqrt(dstr *v) |
538 | { |
539 | mp *p = *(mp **)v[0].buf; |
540 | mp *x = *(mp **)v[1].buf; |
541 | mp *r = *(mp **)v[2].buf; |
542 | mp *c; |
543 | int ok = 1; |
544 | gfreduce rr; |
545 | |
546 | gfreduce_create(&rr, p); |
547 | c = gfreduce_sqrt(&rr, MP_NEW, x); |
548 | if (!MP_EQ(c, r)) { |
549 | fprintf(stderr, "\n*** sqrt failed\n*** "); |
550 | fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16); |
551 | fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16); |
552 | fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16); |
553 | fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16); |
554 | fprintf(stderr, "\n"); |
555 | ok = 0; |
556 | } |
557 | gfreduce_destroy(&rr); |
558 | mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c); |
559 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
560 | return (ok); |
561 | } |
562 | |
563 | static int vtr(dstr *v) |
564 | { |
565 | mp *p = *(mp **)v[0].buf; |
566 | mp *x = *(mp **)v[1].buf; |
567 | int r = *(int *)v[2].buf, c; |
568 | int ok = 1; |
569 | gfreduce rr; |
570 | |
571 | gfreduce_create(&rr, p); |
572 | c = gfreduce_trace(&rr, x); |
573 | if (c != r) { |
574 | fprintf(stderr, "\n*** trace failed\n*** "); |
575 | fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16); |
576 | fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16); |
577 | fprintf(stderr, "\n*** c = %d", c); |
578 | fprintf(stderr, "\n*** r = %d", r); |
579 | fprintf(stderr, "\n"); |
580 | ok = 0; |
581 | } |
582 | gfreduce_destroy(&rr); |
583 | mp_drop(p); mp_drop(x); |
584 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
585 | return (ok); |
586 | } |
587 | |
588 | static int vhftr(dstr *v) |
589 | { |
590 | mp *p = *(mp **)v[0].buf; |
591 | mp *x = *(mp **)v[1].buf; |
592 | mp *r = *(mp **)v[2].buf; |
593 | mp *c; |
594 | int ok = 1; |
595 | gfreduce rr; |
596 | |
597 | gfreduce_create(&rr, p); |
598 | c = gfreduce_halftrace(&rr, MP_NEW, x); |
599 | if (!MP_EQ(c, r)) { |
600 | fprintf(stderr, "\n*** halftrace failed\n*** "); |
601 | fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16); |
602 | fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16); |
603 | fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16); |
604 | fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16); |
605 | fprintf(stderr, "\n"); |
606 | ok = 0; |
607 | } |
608 | gfreduce_destroy(&rr); |
609 | mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c); |
610 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
611 | return (ok); |
612 | } |
613 | |
614 | static int vquad(dstr *v) |
615 | { |
616 | mp *p = *(mp **)v[0].buf; |
617 | mp *x = *(mp **)v[1].buf; |
618 | mp *r = *(mp **)v[2].buf; |
619 | mp *c; |
620 | int ok = 1; |
621 | gfreduce rr; |
622 | |
623 | gfreduce_create(&rr, p); |
624 | c = gfreduce_quadsolve(&rr, MP_NEW, x); |
625 | if (!MP_EQ(c, r)) { |
626 | fprintf(stderr, "\n*** quadsolve failed\n*** "); |
627 | fprintf(stderr, "\n*** p = "); mp_writefile(p, stderr, 16); |
628 | fprintf(stderr, "\n*** x = "); mp_writefile(x, stderr, 16); |
629 | fprintf(stderr, "\n*** c = "); mp_writefile(c, stderr, 16); |
630 | fprintf(stderr, "\n*** r = "); mp_writefile(r, stderr, 16); |
631 | fprintf(stderr, "\n"); |
632 | ok = 0; |
633 | } |
634 | gfreduce_destroy(&rr); |
635 | mp_drop(p); mp_drop(r); mp_drop(x); mp_drop(c); |
636 | assert(mparena_count(MPARENA_GLOBAL) == 0); |
637 | return (ok); |
638 | } |
639 | |
640 | static test_chunk defs[] = { |
641 | { "reduce", vreduce, { &type_mp, &type_mp, &type_mp, 0 } }, |
642 | { "modexp", vmodexp, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } }, |
643 | { "sqrt", vsqrt, { &type_mp, &type_mp, &type_mp, 0 } }, |
644 | { "trace", vtr, { &type_mp, &type_mp, &type_int, 0 } }, |
645 | { "halftrace", vhftr, { &type_mp, &type_mp, &type_mp, 0 } }, |
646 | { "quadsolve", vquad, { &type_mp, &type_mp, &type_mp, 0 } }, |
647 | { 0, 0, { 0 } } |
648 | }; |
649 | |
650 | int main(int argc, char *argv[]) |
651 | { |
652 | test_run(argc, argv, defs, SRCDIR"/tests/gfreduce"); |
653 | return (0); |
654 | } |
655 | |
656 | #endif |
657 | |
658 | /*----- That's all, folks -------------------------------------------------*/ |