1 /*
2     Copyright (C) 2020 Daniel Schultz
3 
4     This file is part of FLINT.
5 
6     FLINT is free software: you can redistribute it and/or modify it under
7     the terms of the GNU Lesser General Public License (LGPL) as published
8     by the Free Software Foundation; either version 2.1 of the License, or
9     (at your option) any later version.  See <https://www.gnu.org/licenses/>.
10 */
11 
12 #include "n_poly.h"
13 #include "nmod_vec.h"
14 
n_fq_poly_is_canonical(const n_poly_t A,const fq_nmod_ctx_t ctx)15 int n_fq_poly_is_canonical(const n_poly_t A, const fq_nmod_ctx_t ctx)
16 {
17     slong d = fq_nmod_ctx_degree(ctx);
18     slong i;
19 
20     if (A->length < 0)
21         return 0;
22 
23     if (d*A->length > A->alloc)
24         return 0;
25 
26     for (i = 0; i < A->length; i++)
27     {
28         if (!n_fq_is_canonical(A->coeffs + d*i, ctx))
29             return 0;
30         if (i + 1 == A->length && _n_fq_is_zero(A->coeffs + d*i, d))
31             return 0;
32     }
33     return 1;
34 }
35 
n_fq_poly_init2(n_fq_poly_t A,slong alloc,const fq_nmod_ctx_t ctx)36 void n_fq_poly_init2(
37     n_fq_poly_t A,
38     slong alloc,
39     const fq_nmod_ctx_t ctx)
40 {
41     slong d = fq_nmod_ctx_degree(ctx);
42     if (alloc > 0)
43     {
44         A->alloc = d*alloc;
45         A->coeffs = flint_malloc(A->alloc*sizeof(mp_limb_t));
46     }
47     else
48     {
49         A->alloc = 0;
50         A->coeffs = NULL;
51     }
52     A->length = 0;
53 }
54 
n_fq_poly_print_pretty(const n_fq_poly_t A,const char * x,const fq_nmod_ctx_t ctx)55 void n_fq_poly_print_pretty(
56     const n_fq_poly_t A,
57     const char * x,
58     const fq_nmod_ctx_t ctx)
59 {
60     slong d = fq_nmod_ctx_degree(ctx);
61     slong i;
62     int first;
63 
64     first = 1;
65     for (i = A->length - 1; i >= 0; i--)
66     {
67         if (i + 1 != A->length && _n_fq_is_zero(A->coeffs + d*i, d))
68             continue;
69 
70         if (!first)
71             flint_printf(" + ");
72 
73         first = 0;
74 
75         flint_printf("(");
76         n_fq_print_pretty(A->coeffs + d*i, ctx);
77         flint_printf(")*%s^%wd", x, i);
78     }
79 
80     if (first)
81         flint_printf("0");
82 }
83 
84 
n_fq_poly_randtest(n_poly_t A,flint_rand_t state,slong len,const fq_nmod_ctx_t ctx)85 void n_fq_poly_randtest(
86     n_poly_t A,
87     flint_rand_t state,
88     slong len,
89     const fq_nmod_ctx_t ctx)
90 {
91     slong d = fq_nmod_ctx_degree(ctx);
92     slong i;
93 
94     if (len < 1)
95     {
96         A->length = 0;
97         return;
98     }
99 
100     n_poly_fit_length(A, d*len);
101 
102     for (i = 0; i < d*len; i++)
103         A->coeffs[i] = n_randint(state, ctx->mod.n);
104 
105     A->length = len;
106     _n_fq_poly_normalise(A, d);
107 }
108 
109 
_n_fq_poly_one(n_poly_t A,slong d)110 void _n_fq_poly_one(n_poly_t A, slong d)
111 {
112     n_poly_fit_length(A, d);
113     A->length = 1;
114     _n_fq_one(A->coeffs + d*0, d);
115 }
116 
117 
n_fq_poly_is_one(n_poly_t A,const fq_nmod_ctx_t ctx)118 int n_fq_poly_is_one(n_poly_t A, const fq_nmod_ctx_t ctx)
119 {
120     return A->length == 1 && _n_fq_is_one(A->coeffs, fq_nmod_ctx_degree(ctx));
121 }
122 
123 
n_fq_poly_get_coeff_n_fq(mp_limb_t * c,const n_fq_poly_t A,slong e,const fq_nmod_ctx_t ctx)124 void n_fq_poly_get_coeff_n_fq(
125     mp_limb_t * c,
126     const n_fq_poly_t A,
127     slong e,
128     const fq_nmod_ctx_t ctx)
129 {
130     slong d = fq_nmod_ctx_degree(ctx);
131 
132     if (e >= A->length)
133         _n_fq_zero(c, d);
134     else
135         _n_fq_set(c, A->coeffs + d*e, d);
136 }
137 
n_fq_poly_get_coeff_fq_nmod(fq_nmod_t c,const n_fq_poly_t A,slong e,const fq_nmod_ctx_t ctx)138 void n_fq_poly_get_coeff_fq_nmod(
139     fq_nmod_t c,
140     const n_fq_poly_t A,
141     slong e,
142     const fq_nmod_ctx_t ctx)
143 {
144     slong d = fq_nmod_ctx_degree(ctx);
145 
146     if (e >= A->length)
147         fq_nmod_zero(c, ctx);
148     else
149         n_fq_get_fq_nmod(c, A->coeffs + d*e, ctx);
150 }
151 
n_fq_poly_set_coeff_n_fq(n_fq_poly_t A,slong j,const mp_limb_t * c,const fq_nmod_ctx_t ctx)152 void n_fq_poly_set_coeff_n_fq(
153     n_fq_poly_t A,
154     slong j,
155     const mp_limb_t * c,
156     const fq_nmod_ctx_t ctx)
157 {
158     slong d = fq_nmod_ctx_degree(ctx);
159 
160     FLINT_ASSERT(n_fq_is_canonical(c, ctx));
161     FLINT_ASSERT(n_fq_poly_is_canonical(A, ctx));
162 
163     n_poly_fit_length(A, d*(j + 1));
164 
165     if (j + 1 <= A->length)
166     {
167         _n_fq_set(A->coeffs + d*j, c, d);
168         if (j + 1 == A->length)
169             _n_fq_poly_normalise(A, d);
170     }
171     else if (!_n_fq_is_zero(c, d)) /* extend polynomial */
172     {
173         flint_mpn_zero(A->coeffs + d*A->length, d*(j - A->length));
174         _n_fq_set(A->coeffs + d*j, c, d);
175         A->length = j + 1;
176     }
177 
178     FLINT_ASSERT(n_fq_poly_is_canonical(A, ctx));
179 }
180 
n_fq_poly_set_coeff_fq_nmod(n_poly_t A,slong j,const fq_nmod_t c,const fq_nmod_ctx_t ctx)181 void n_fq_poly_set_coeff_fq_nmod(
182     n_poly_t A,
183     slong j,
184     const fq_nmod_t c,
185     const fq_nmod_ctx_t ctx)
186 {
187     slong d = fq_nmod_ctx_degree(ctx);
188 
189     FLINT_ASSERT(n_fq_poly_is_canonical(A, ctx));
190 
191     n_poly_fit_length(A, d*(j + 1));
192 
193     if (j + 1 <= A->length)
194     {
195         n_fq_set_fq_nmod(A->coeffs + d*j, c, ctx);
196         if (j + 1 == A->length)
197             _n_fq_poly_normalise(A, d);
198     }
199     else if (!fq_nmod_is_zero(c, ctx))
200     {
201         flint_mpn_zero(A->coeffs + d*A->length, d*(j - A->length));
202         n_fq_set_fq_nmod(A->coeffs + d*j, c, ctx);
203         A->length = j + 1;
204     }
205 
206     FLINT_ASSERT(n_fq_poly_is_canonical(A, ctx));
207 }
208 
209 
n_fq_poly_scalar_mul_ui(n_poly_t A,const n_poly_t B,ulong c,const fq_nmod_ctx_t ctx)210 void n_fq_poly_scalar_mul_ui(
211     n_poly_t A,
212     const n_poly_t B,
213     ulong c,
214     const fq_nmod_ctx_t ctx)
215 {
216     slong d = fq_nmod_ctx_degree(ctx);
217 
218     if (c >= ctx->mod.n)
219     {
220         NMOD_RED(c, c, ctx->mod);
221     }
222 
223     if (B->length < 1 || c == 0)
224     {
225         n_poly_zero(A);
226         return;
227     }
228 
229     n_poly_fit_length(A, d*B->length);
230     _nmod_vec_scalar_mul_nmod(A->coeffs, B->coeffs, d*B->length, c, ctx->mod);
231     A->length = B->length;
232     _n_fq_poly_normalise(A, d);
233 }
234 
235 
n_fq_poly_equal(const n_poly_t A,const n_poly_t B,const fq_nmod_ctx_t ctx)236 int n_fq_poly_equal(
237     const n_poly_t A,
238     const n_poly_t B,
239     const fq_nmod_ctx_t ctx)
240 {
241     slong d = fq_nmod_ctx_degree(ctx);
242     slong i;
243 
244     if (A->length != B->length)
245         return 0;
246 
247     for (i = 0; i < d*B->length; i++)
248         if (A->coeffs[i] != B->coeffs[i])
249             return 0;
250 
251     return 1;
252 }
253 
n_fq_poly_set(n_poly_t A,const n_poly_t B,const fq_nmod_ctx_t ctx)254 void n_fq_poly_set(
255     n_poly_t A,
256     const n_poly_t B,
257     const fq_nmod_ctx_t ctx)
258 {
259     slong d = fq_nmod_ctx_degree(ctx);
260 
261     if (A == B)
262         return;
263 
264     n_poly_fit_length(A, d*B->length);
265     _nmod_vec_set(A->coeffs, B->coeffs, d*B->length);
266     A->length = B->length;
267 }
268 
n_fq_poly_set_n_fq(n_poly_t A,const mp_limb_t * c,const fq_nmod_ctx_t ctx)269 void n_fq_poly_set_n_fq(
270     n_poly_t A,
271     const mp_limb_t * c,
272     const fq_nmod_ctx_t ctx)
273 {
274     slong d = fq_nmod_ctx_degree(ctx);
275 
276     n_poly_fit_length(A, d);
277     _nmod_vec_set(A->coeffs, c, d);
278     A->length = 1;
279     _n_fq_poly_normalise(A, d);
280 }
281 
n_fq_poly_set_fq_nmod(n_poly_t A,const fq_nmod_t c,const fq_nmod_ctx_t ctx)282 void n_fq_poly_set_fq_nmod(
283     n_poly_t A,
284     const fq_nmod_t c,
285     const fq_nmod_ctx_t ctx)
286 {
287     slong d = fq_nmod_ctx_degree(ctx);
288 
289     n_poly_fit_length(A, d);
290     n_fq_set_fq_nmod(A->coeffs, c, ctx);
291     A->length = 1;
292     _n_fq_poly_normalise(A, d);
293 }
294 
n_fq_poly_shift_right(n_poly_t A,const n_poly_t B,slong n,const fq_nmod_ctx_t ctx)295 void n_fq_poly_shift_right(
296     n_poly_t A,
297     const n_poly_t B,
298     slong n,
299     const fq_nmod_ctx_t ctx)
300 {
301     slong d = fq_nmod_ctx_degree(ctx);
302 
303     if (n < 1)
304     {
305         n_fq_poly_set(A, B, ctx);
306         return;
307     }
308     else if (B->length <= n)
309     {
310         n_poly_zero(A);
311         return;
312     }
313     else
314     {
315         n_poly_fit_length(A, d*(B->length - n));
316         flint_mpn_copyi(A->coeffs, B->coeffs + d*n, d*(B->length - n));
317         A->length = B->length - n;
318     }
319 }
320 
n_fq_poly_shift_left(n_poly_t A,const n_poly_t B,slong n,const fq_nmod_ctx_t ctx)321 void n_fq_poly_shift_left(
322     n_poly_t A,
323     const n_poly_t B,
324     slong n,
325     const fq_nmod_ctx_t ctx)
326 {
327     slong d = fq_nmod_ctx_degree(ctx);
328 
329     if (n < 1)
330     {
331         n_fq_poly_set(A, B, ctx);
332     }
333     else if (n_fq_poly_is_zero(B))
334     {
335         n_fq_poly_zero(A);
336     }
337     else
338     {
339         n_poly_fit_length(A, d*(B->length + n));
340         flint_mpn_copyd(A->coeffs + d*n, B->coeffs, d*B->length);
341         flint_mpn_zero(A->coeffs, d*n);
342         A->length = B->length + n;
343     }
344 }
345 
n_fq_poly_truncate(n_poly_t A,slong len,const fq_nmod_ctx_t ctx)346 void n_fq_poly_truncate(n_poly_t A, slong len, const fq_nmod_ctx_t ctx)
347 {
348     slong d = fq_nmod_ctx_degree(ctx);
349     if (A->length > len)
350     {
351         A->length = len;
352         _n_fq_poly_normalise(A, d);
353     }
354 }
355 
356 
n_fq_poly_evaluate_fq_nmod(fq_nmod_t e,const n_poly_t A,const fq_nmod_t c,const fq_nmod_ctx_t ctx)357 void n_fq_poly_evaluate_fq_nmod(
358     fq_nmod_t e,
359     const n_poly_t A,
360     const fq_nmod_t c,
361     const fq_nmod_ctx_t ctx)
362 {
363     fq_nmod_poly_t AA;
364     fq_nmod_poly_init(AA, ctx);
365     n_fq_poly_get_fq_nmod_poly(AA, A, ctx);
366     fq_nmod_poly_evaluate_fq_nmod(e, AA, c, ctx);
367     fq_nmod_poly_clear(AA, ctx);
368 }
369 
370 
n_fq_poly_evaluate_n_fq(mp_limb_t * e,const n_poly_t A,const mp_limb_t * c,const fq_nmod_ctx_t ctx)371 void n_fq_poly_evaluate_n_fq(
372     mp_limb_t * e,
373     const n_poly_t A,
374     const mp_limb_t * c,
375     const fq_nmod_ctx_t ctx)
376 {
377     slong d = fq_nmod_ctx_degree(ctx);
378     slong i;
379     mp_limb_t * u = FLINT_ARRAY_ALLOC(d, mp_limb_t);
380     mp_limb_t * t = FLINT_ARRAY_ALLOC(d, mp_limb_t);
381 
382     _n_fq_zero(t, d);
383     for (i = 0; i < A->length; i++)
384     {
385         n_fq_pow_ui(u, c, i, ctx);
386         n_fq_mul(u, u, A->coeffs + d*i, ctx);
387         n_fq_add(t, t, u, ctx);
388     }
389 
390     _n_fq_set(e, t, d);
391 
392     flint_free(u);
393     flint_free(t);
394 }
395 
n_fq_poly_eval_pow(mp_limb_t * ev,const n_fq_poly_t P,n_fq_poly_t alphapow,const fq_nmod_ctx_t ctx)396 void n_fq_poly_eval_pow(
397     mp_limb_t * ev,
398     const n_fq_poly_t P,
399     n_fq_poly_t alphapow,
400     const fq_nmod_ctx_t ctx)
401 {
402     slong d = fq_nmod_ctx_degree(ctx);
403     const mp_limb_t * Pcoeffs = P->coeffs;
404     slong i, Plen = P->length;
405     mp_limb_t * alpha_powers = alphapow->coeffs;
406     mp_limb_t * t;
407     slong k;
408     TMP_INIT;
409 
410     TMP_START;
411     t = TMP_ALLOC(d*FLINT_MAX(N_FQ_MUL_ITCH, N_FQ_LAZY_ITCH)*sizeof(mp_limb_t));
412 
413     if (Plen > alphapow->length)
414     {
415         slong oldlength = alphapow->length;
416         FLINT_ASSERT(2 <= oldlength);
417         n_poly_fit_length(alphapow, d*Plen);
418         alphapow->length = Plen;
419         alpha_powers = alphapow->coeffs;
420         for (k = oldlength; k < Plen; k++)
421         {
422             _n_fq_mul(alpha_powers + d*k, alpha_powers + d*(k - 1),
423                                           alpha_powers + d*1, ctx, t);
424         }
425     }
426 
427     _nmod_vec_zero(t, 6*d);
428 
429     switch (_n_fq_dot_lazy_size(Plen, ctx))
430     {
431 #define lazycase(n)                                                           \
432 case n:                                                                       \
433     for (i = 0; i < Plen; i++)                                                \
434         _n_fq_madd2_lazy##n(t, Pcoeffs + d*i, alpha_powers + d*i, d);         \
435     _n_fq_reduce2_lazy##n(t, d, ctx->mod);                                    \
436     break;
437 
438     lazycase(1)
439     lazycase(2)
440     lazycase(3)
441 
442     default:
443         for (i = 0; i < Plen; i++)
444             _n_fq_madd2(t, Pcoeffs + d*i, alpha_powers + d*i, ctx, t + 2*d);
445         break;
446     }
447 
448     _n_fq_reduce2(ev, t, ctx, t + 2*d);
449 
450     TMP_END;
451 }
452 
453 
n_fq_poly_set_fq_nmod_poly(n_poly_t A,const fq_nmod_poly_t B,const fq_nmod_ctx_t ctx)454 void n_fq_poly_set_fq_nmod_poly(
455     n_poly_t A,
456     const fq_nmod_poly_t B,
457     const fq_nmod_ctx_t ctx)
458 {
459     slong i;
460     slong d = fq_nmod_ctx_degree(ctx);
461 
462     n_poly_fit_length(A, d*B->length);
463 
464     for (i = 0; i < B->length; i++)
465         n_fq_set_fq_nmod(A->coeffs + d*i, B->coeffs + i, ctx);
466 
467     A->length = B->length;
468 }
469 
n_fq_poly_get_fq_nmod_poly(fq_nmod_poly_t A,const n_poly_t B,const fq_nmod_ctx_t ctx)470 void n_fq_poly_get_fq_nmod_poly(
471     fq_nmod_poly_t A,
472     const n_poly_t B,
473     const fq_nmod_ctx_t ctx)
474 {
475     slong d = fq_nmod_ctx_degree(ctx);
476     slong i;
477 
478     FLINT_ASSERT(B->alloc >= d*B->length);
479 
480     fq_nmod_poly_fit_length(A, B->length, ctx);
481 
482     for (i = 0; i < B->length; i++)
483         n_fq_get_fq_nmod(A->coeffs + i, B->coeffs + d*i, ctx);
484 
485     A->length = B->length;
486 }
487 
488 
n_fq_poly_scalar_mul_n_fq(n_poly_t A,const n_poly_t B,const mp_limb_t * c,const fq_nmod_ctx_t ctx)489 void n_fq_poly_scalar_mul_n_fq(
490     n_poly_t A,
491     const n_poly_t B,
492     const mp_limb_t * c,
493     const fq_nmod_ctx_t ctx)
494 {
495     slong i, d = fq_nmod_ctx_degree(ctx);
496     n_poly_fit_length(A, d*B->length);
497     for (i = 0; i < B->length; i++)
498         n_fq_mul(A->coeffs + d*i, B->coeffs + d*i, c, ctx);
499     A->length = B->length;
500     _n_fq_poly_normalise(A, d);
501 }
502 
n_fq_poly_make_monic(n_poly_t A,const n_poly_t B,const fq_nmod_ctx_t ctx)503 void n_fq_poly_make_monic(
504     n_poly_t A,
505     const n_poly_t B,
506     const fq_nmod_ctx_t ctx)
507 {
508     slong d = fq_nmod_ctx_degree(ctx);
509     slong itch = FLINT_MAX(N_FQ_MUL_ITCH, N_FQ_INV_ITCH);
510     mp_limb_t * tmp, * inv;
511     slong i, Blen = B->length;
512 
513     if (Blen < 1)
514     {
515         n_poly_zero(A);
516         return;
517     }
518 
519     n_poly_fit_length(A, d*Blen);
520 
521     tmp = FLINT_ARRAY_ALLOC(d*(itch + 1), mp_limb_t);
522     inv = tmp + d*itch;
523 
524     _n_fq_inv(inv, B->coeffs + d*(Blen - 1), ctx, tmp);
525 
526     for (i = 0; i < Blen - 1; i++)
527         _n_fq_mul(A->coeffs + d*i, B->coeffs + d*i, inv, ctx, tmp);
528 
529     _n_fq_one(A->coeffs + d*(Blen - 1), d);
530     A->length = Blen;
531 
532     flint_free(tmp);
533 }
534 
535 
536 
n_fq_poly_scalar_addmul_n_fq(n_fq_poly_t A,const n_fq_poly_t B,const n_fq_poly_t C,const mp_limb_t * s,const fq_nmod_ctx_t ctx)537 void n_fq_poly_scalar_addmul_n_fq(
538     n_fq_poly_t A,
539     const n_fq_poly_t B,
540     const n_fq_poly_t C,
541     const mp_limb_t * s,
542     const fq_nmod_ctx_t ctx)
543 {
544     slong d = fq_nmod_ctx_degree(ctx);
545     slong i;
546     mp_limb_t * Acoeffs;
547     mp_limb_t * Bcoeffs;
548     mp_limb_t * Ccoeffs;
549     slong Blen = B->length;
550     slong Clen = C->length;
551     mp_limb_t * t;
552     TMP_INIT;
553 
554     n_poly_fit_length(A, d*FLINT_MAX(Blen, Clen));
555 
556     Bcoeffs = B->coeffs;
557     Ccoeffs = C->coeffs;
558 
559     TMP_START;
560     t = TMP_ALLOC(d*N_FQ_MUL_ITCH*sizeof(mp_limb_t));
561 
562     if (Blen > Clen)
563     {
564         n_poly_fit_length(A, d*Blen);
565         Acoeffs = A->coeffs;
566         for (i = 0; i < Clen; i++)
567             _n_fq_addmul(Acoeffs + d*i, Bcoeffs + d*i, Ccoeffs + d*i, s, ctx, t);
568         if (A != B)
569             _nmod_vec_set(Acoeffs + d*Clen, Bcoeffs + d*Clen, d*(Blen - Clen));
570         A->length = Blen;
571     }
572     else if (Blen < Clen)
573     {
574         n_poly_fit_length(A, d*Clen);
575         Acoeffs = A->coeffs;
576         for (i = 0; i < Blen; i++)
577             _n_fq_addmul(Acoeffs + d*i, Bcoeffs + d*i, Ccoeffs + d*i, s, ctx, t);
578         for ( ; i < Clen; i++)
579             _n_fq_mul(Acoeffs + d*i, Ccoeffs + d*i, s, ctx, t);
580         A->length = Clen;
581     }
582     else
583     {
584         n_poly_fit_length(A, d*Blen);
585         Acoeffs = A->coeffs;
586         for (i = 0; i < Blen; i++)
587             _n_fq_addmul(Acoeffs + d*i, Bcoeffs + d*i, Ccoeffs + d*i, s, ctx, t);
588         A->length = Blen;
589         _n_fq_poly_normalise(A, d);
590     }
591 
592     TMP_END;
593 }
594 
595 
596 /* multiply A by (x^k - c) */
n_fq_poly_shift_left_scalar_submul(n_poly_t A,slong k,const mp_limb_t * c,const fq_nmod_ctx_t ctx)597 void n_fq_poly_shift_left_scalar_submul(
598     n_poly_t A,
599     slong k,
600     const mp_limb_t * c,
601     const fq_nmod_ctx_t ctx)
602 {
603     slong d = fq_nmod_ctx_degree(ctx);
604     mp_limb_t * Acoeffs;
605     slong i;
606     slong Alen = A->length;
607     mp_limb_t * u = FLINT_ARRAY_ALLOC(d, mp_limb_t);
608 
609     n_poly_fit_length(A, d*(Alen + k));
610 
611     Acoeffs = A->coeffs;
612 
613     flint_mpn_copyd(Acoeffs + d*k, Acoeffs, d*Alen);
614     flint_mpn_zero(Acoeffs, d*k);
615 
616     for (i = 0; i < A->length; i++)
617     {
618         n_fq_mul(u, c, Acoeffs + d*(i + k), ctx);
619         _n_fq_sub(Acoeffs + d*i, Acoeffs + d*i, u, d, fq_nmod_ctx_mod(ctx));
620     }
621 
622     A->length = Alen + k;
623 
624     flint_free(u);
625 }
626 
627 
n_fq_poly_remove(n_poly_t f,const n_poly_t g,const fq_nmod_ctx_t ctx)628 ulong n_fq_poly_remove(
629     n_poly_t f,
630     const n_poly_t g,
631     const fq_nmod_ctx_t ctx)
632 {
633     n_poly_t q, r;
634     ulong i = 0;
635 
636     n_poly_init(q);
637     n_poly_init(r);
638 
639     while (1)
640     {
641         if (f->length < g->length)
642             break;
643         n_fq_poly_divrem(q, r, f, g, ctx);
644         if (r->length == 0)
645             n_poly_swap(q, f);
646         else
647             break;
648         i++;
649     }
650 
651     n_poly_clear(q);
652     n_poly_clear(r);
653 
654     return i;
655 }
656