My Project
Loading...
Searching...
No Matches
facMul.cc
Go to the documentation of this file.
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file facMul.cc
5 *
6 * This file implements functions for fast multiplication and division with
7 * remainder.
8 *
9 * Nomenclature rules: kronSub* -> plain Kronecker substitution
10 * reverseSubst* -> reverse Kronecker substitution
11 * kronSubRecipro* -> reciprocal Kronecker substitution as
12 * described in D. Harvey "Faster
13 * polynomial multiplication via
14 * multipoint Kronecker substitution"
15 *
16 * @author Martin Lee
17 *
18 **/
19/*****************************************************************************/
20
21#include "debug.h"
22
23#include "config.h"
24
25
26#include "canonicalform.h"
27#include "facMul.h"
28#include "cf_util.h"
29#include "cf_iter.h"
30#include "cf_algorithm.h"
32
33#ifdef HAVE_NTL
34#include <NTL/lzz_pEX.h>
35#include "NTLconvert.h"
36#endif
37
38#ifdef HAVE_FLINT
39#include "FLINTconvert.h"
40#endif
41
42// univariate polys
43#if defined(HAVE_NTL) || defined(HAVE_FLINT)
44
45#ifdef HAVE_FLINT
46void kronSubQa (fmpz_poly_t result, const CanonicalForm& A, int d)
47{
48 int degAy= degree (A);
49 fmpz_poly_init2 (result, d*(degAy + 1));
50 _fmpz_poly_set_length (result, d*(degAy + 1));
52 for (CFIterator i= A; i.hasTerms(); i++)
53 {
54 if (i.coeff().inBaseDomain())
55 convertCF2initFmpz (fmpz_poly_get_coeff_ptr (result, i.exp()*d), i.coeff());
56 else
57 for (j= i.coeff(); j.hasTerms(); j++)
58 convertCF2initFmpz (fmpz_poly_get_coeff_ptr (result, i.exp()*d+j.exp()),
59 j.coeff());
60 }
61 _fmpz_poly_normalise(result);
62}
63
64
66reverseSubstQa (const fmpz_poly_t F, int d, const Variable& x,
67 const Variable& alpha, const CanonicalForm& den)
68{
70 int i= 0;
71 int degf= fmpz_poly_degree (F);
72 int k= 0;
73 int degfSubK;
74 int repLength;
75 fmpq_poly_t buf;
76 fmpq_poly_t mipo;
78 while (degf >= k)
79 {
80 degfSubK= degf - k;
81 if (degfSubK >= d)
82 repLength= d;
83 else
84 repLength= degfSubK + 1;
85
86 fmpq_poly_init2 (buf, repLength);
87 _fmpq_poly_set_length (buf, repLength);
88 _fmpz_vec_set (buf->coeffs, F->coeffs + k, repLength);
89 _fmpq_poly_normalise (buf);
90 fmpq_poly_rem (buf, buf, mipo);
91
93 fmpq_poly_clear (buf);
94 i++;
95 k= d*i;
96 }
97 fmpq_poly_clear (mipo);
98 result /= den;
99 return result;
100}
101
104 const Variable& alpha)
105{
106 CanonicalForm A= F;
108
109 CanonicalForm denA= bCommonDen (A);
110 CanonicalForm denB= bCommonDen (B);
111
112 A *= denA;
113 B *= denB;
114 int degAa= degree (A, alpha);
115 int degBa= degree (B, alpha);
116 int d= degAa + 1 + degBa;
117
118 fmpz_poly_t FLINTA,FLINTB;
119 kronSubQa (FLINTA, A, d);
120 kronSubQa (FLINTB, B, d);
121
122 fmpz_poly_mul (FLINTA, FLINTA, FLINTB);
123
124 denA *= denB;
125 A= reverseSubstQa (FLINTA, d, F.mvar(), alpha, denA);
126
127 fmpz_poly_clear (FLINTA);
128 fmpz_poly_clear (FLINTB);
129 return A;
130}
131
134{
135 CanonicalForm A= F;
137
138 CanonicalForm denA= bCommonDen (A);
139 CanonicalForm denB= bCommonDen (B);
140
141 A *= denA;
142 B *= denB;
143 fmpz_poly_t FLINTA,FLINTB;
144 convertFacCF2Fmpz_poly_t (FLINTA, A);
145 convertFacCF2Fmpz_poly_t (FLINTB, B);
146 fmpz_poly_mul (FLINTA, FLINTA, FLINTB);
147 denA *= denB;
148 A= convertFmpz_poly_t2FacCF (FLINTA, F.mvar());
149 A /= denA;
150 fmpz_poly_clear (FLINTA);
151 fmpz_poly_clear (FLINTB);
152
153 return A;
154}
155
156/*CanonicalForm
157mulFLINTQ2 (const CanonicalForm& F, const CanonicalForm& G)
158{
159 CanonicalForm A= F;
160 CanonicalForm B= G;
161
162 fmpq_poly_t FLINTA,FLINTB;
163 convertFacCF2Fmpq_poly_t (FLINTA, A);
164 convertFacCF2Fmpq_poly_t (FLINTB, B);
165
166 fmpq_poly_mul (FLINTA, FLINTA, FLINTB);
167 A= convertFmpq_poly_t2FacCF (FLINTA, F.mvar());
168 fmpq_poly_clear (FLINTA);
169 fmpq_poly_clear (FLINTB);
170 return A;
171}*/
172
175{
176 CanonicalForm A= F;
178
179 fmpq_poly_t FLINTA,FLINTB;
180 convertFacCF2Fmpq_poly_t (FLINTA, A);
181 convertFacCF2Fmpq_poly_t (FLINTB, B);
182
183 fmpq_poly_div (FLINTA, FLINTA, FLINTB);
184 A= convertFmpq_poly_t2FacCF (FLINTA, F.mvar());
185
186 fmpq_poly_clear (FLINTA);
187 fmpq_poly_clear (FLINTB);
188 return A;
189}
190
193{
194 CanonicalForm A= F;
196
197 fmpq_poly_t FLINTA,FLINTB;
198 convertFacCF2Fmpq_poly_t (FLINTA, A);
199 convertFacCF2Fmpq_poly_t (FLINTB, B);
200
201 fmpq_poly_rem (FLINTA, FLINTA, FLINTB);
202 A= convertFmpq_poly_t2FacCF (FLINTA, F.mvar());
203
204 fmpq_poly_clear (FLINTA);
205 fmpq_poly_clear (FLINTB);
206 return A;
207}
208
211 const Variable& alpha, int m)
212{
213 CanonicalForm A= F;
215
216 CanonicalForm denA= bCommonDen (A);
217 CanonicalForm denB= bCommonDen (B);
218
219 A *= denA;
220 B *= denB;
221
222 int degAa= degree (A, alpha);
223 int degBa= degree (B, alpha);
224 int d= degAa + 1 + degBa;
225
226 fmpz_poly_t FLINTA,FLINTB;
227 kronSubQa (FLINTA, A, d);
228 kronSubQa (FLINTB, B, d);
229
230 int k= d*m;
231 fmpz_poly_mullow (FLINTA, FLINTA, FLINTB, k);
232
233 denA *= denB;
234 A= reverseSubstQa (FLINTA, d, F.mvar(), alpha, denA);
235 fmpz_poly_clear (FLINTA);
236 fmpz_poly_clear (FLINTB);
237 return A;
238}
239
242{
243 if (F.inCoeffDomain() && G.inCoeffDomain())
244 return F*G;
245 if (F.inCoeffDomain())
246 return mod (F*G, power (G.mvar(), m));
247 if (G.inCoeffDomain())
248 return mod (F*G, power (F.mvar(), m));
251 return mulFLINTQaTrunc (F, G, alpha, m);
252
253 CanonicalForm A= F;
255
256 CanonicalForm denA= bCommonDen (A);
257 CanonicalForm denB= bCommonDen (B);
258
259 A *= denA;
260 B *= denB;
261 fmpz_poly_t FLINTA,FLINTB;
262 convertFacCF2Fmpz_poly_t (FLINTA, A);
263 convertFacCF2Fmpz_poly_t (FLINTB, B);
264 fmpz_poly_mullow (FLINTA, FLINTA, FLINTB, m);
265 denA *= denB;
266 A= convertFmpz_poly_t2FacCF (FLINTA, F.mvar());
267 A /= denA;
268 fmpz_poly_clear (FLINTA);
269 fmpz_poly_clear (FLINTB);
270
271 return A;
272}
273
275{
276 if (d == 0)
277 return F;
278 if (F.inCoeffDomain())
279 return F*power (x,d);
281 CFIterator i= F;
282 while (d - i.exp() < 0)
283 i++;
284
285 for (; i.hasTerms() && (d - i.exp() >= 0); i++)
286 result += i.coeff()*power (x, d - i.exp());
287 return result;
288}
289
291newtonInverse (const CanonicalForm& F, const int n, const Variable& x)
292{
293 int l= ilog2(n);
294
296 if (F.inCoeffDomain())
297 g= F;
298 else
299 g= F [0];
300
301 if (!F.inCoeffDomain())
302 ASSERT (F.mvar() == x, "main variable of F and x differ");
303 ASSERT (!g.isZero(), "expected a unit");
304
305 if (!g.isOne())
306 g = 1/g;
308 int exp= 0;
309 if (n & 1)
310 {
311 result= g;
312 exp= 1;
313 }
315
316 for (int i= 1; i <= l; i++)
317 {
318 h= mulNTL (g, mod (F, power (x, (1 << i))));
319 h= mod (h, power (x, (1 << i)) - 1);
320 h= div (h, power (x, (1 << (i - 1))));
321 g -= power (x, (1 << (i - 1)))*
322 mulFLINTQTrunc (g, h, 1 << (i-1));
323
324 if (n & (1 << i))
325 {
326 if (exp)
327 {
328 h= mulNTL (result, mod (F, power (x, exp + (1 << i))));
329 h= mod (h, power (x, exp + (1 << i)) - 1);
330 h= div (h, power (x, exp));
331 result -= power(x, exp)*mulFLINTQTrunc (g, h, 1 << i);
332 exp += (1 << i);
333 }
334 else
335 {
336 exp= (1 << i);
337 result= g;
338 }
339 }
340 }
341
342 return result;
343}
344
345void
348{
349 ASSERT (F.level() == G.level(), "F and G have different level");
350 CanonicalForm A= F;
352 Variable x= A.mvar();
353 int degA= degree (A);
354 int degB= degree (B);
355 int m= degA - degB;
356
357 if (m < 0)
358 {
359 R= A;
360 Q= 0;
361 return;
362 }
363
364 if (degB <= 1)
365 divrem (A, B, Q, R);
366 else
367 {
368 R= uniReverse (A, degA, x);
369
370 CanonicalForm revB= uniReverse (B, degB, x);
371 revB= newtonInverse (revB, m + 1, x);
372 Q= mulFLINTQTrunc (R, revB, m + 1);
373 Q= uniReverse (Q, m, x);
374
375 R= A - mulNTL (Q, B);
376 }
377}
378
379void
381{
382 ASSERT (F.level() == G.level(), "F and G have different level");
383 CanonicalForm A= F;
385 Variable x= A.mvar();
386 int degA= degree (A);
387 int degB= degree (B);
388 int m= degA - degB;
389
390 if (m < 0)
391 {
392 Q= 0;
393 return;
394 }
395
396 if (degB <= 1)
397 Q= div (A, B);
398 else
399 {
400 CanonicalForm R= uniReverse (A, degA, x);
401 CanonicalForm revB= uniReverse (B, degB, x);
402 revB= newtonInverse (revB, m + 1, x);
403 Q= mulFLINTQTrunc (R, revB, m + 1);
404 Q= uniReverse (Q, m, x);
405 }
406}
407
408#endif
409
411mulNTL (const CanonicalForm& F, const CanonicalForm& G, const modpk& b)
412{
414 return F*G;
415 if (getCharacteristic() == 0)
416 {
418 if ((!F.inCoeffDomain() && !G.inCoeffDomain()) &&
420 {
421 if (b.getp() != 0)
422 {
424 bool is_rat= isOn (SW_RATIONAL);
425 if (!is_rat)
426 On (SW_RATIONAL);
427 mipo *=bCommonDen (mipo);
428 if (!is_rat)
430#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
431 fmpz_t FLINTp;
432 fmpz_mod_poly_t FLINTmipo;
433 fq_ctx_t fq_con;
434 fq_poly_t FLINTF, FLINTG;
435
436 fmpz_init (FLINTp);
437
438 convertCF2initFmpz (FLINTp, b.getpk());
439
440 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, mipo, FLINTp);
441
442 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
443 fmpz_mod_ctx_t fmpz_ctx;
444 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
445 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
446 #else
447 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
448 #endif
449
450 convertFacCF2Fq_poly_t (FLINTF, F, fq_con);
452
453 fq_poly_mul (FLINTF, FLINTF, FLINTG, fq_con);
454
456 alpha, fq_con);
457
458 fmpz_clear (FLINTp);
459 fq_poly_clear (FLINTF, fq_con);
460 fq_poly_clear (FLINTG, fq_con);
461 fq_ctx_clear (fq_con);
462 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
463 fmpz_mod_poly_clear (FLINTmipo,fmpz_ctx);
464 fmpz_mod_ctx_clear(fmpz_ctx);
465 #else
466 fmpz_mod_poly_clear(FLINTmipo);
467 #endif
468 return b (result);
469#endif
470#ifdef HAVE_NTL
471 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
472 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (mipo));
473 ZZ_pE::init (NTLmipo);
474 ZZ_pEX NTLg= convertFacCF2NTLZZ_pEX (G, NTLmipo);
475 ZZ_pEX NTLf= convertFacCF2NTLZZ_pEX (F, NTLmipo);
476 mul (NTLf, NTLf, NTLg);
477
478 return b (convertNTLZZ_pEX2CF (NTLf, F.mvar(), alpha));
479#endif
480 }
481#ifdef HAVE_FLINT
483 return result;
484#else
485 return F*G;
486#endif
487 }
488 else if (!F.inCoeffDomain() && !G.inCoeffDomain())
489 {
490#ifdef HAVE_FLINT
491 if (b.getp() != 0)
492 {
493 fmpz_t FLINTpk;
494 fmpz_init (FLINTpk);
495 convertCF2initFmpz (FLINTpk, b.getpk());
496 fmpz_mod_poly_t FLINTF, FLINTG;
497 convertFacCF2Fmpz_mod_poly_t (FLINTF, F, FLINTpk);
498 convertFacCF2Fmpz_mod_poly_t (FLINTG, G, FLINTpk);
499 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
500 fmpz_mod_ctx_t fmpz_ctx;
501 fmpz_mod_ctx_init(fmpz_ctx,FLINTpk);
502 fmpz_mod_poly_mul (FLINTF, FLINTF, FLINTG, fmpz_ctx);
503 #else
504 fmpz_mod_poly_mul (FLINTF, FLINTF, FLINTG);
505 #endif
507 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
508 fmpz_mod_poly_clear (FLINTG,fmpz_ctx);
509 fmpz_mod_poly_clear (FLINTF,fmpz_ctx);
510 fmpz_mod_ctx_clear(fmpz_ctx);
511 #else
512 fmpz_mod_poly_clear (FLINTG);
513 fmpz_mod_poly_clear (FLINTF);
514 #endif
515 fmpz_clear (FLINTpk);
516 return result;
517 }
518 return mulFLINTQ (F, G);
519#endif
520#ifdef HAVE_NTL
521 if (b.getp() != 0)
522 {
523 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
524 ZZX ZZf= convertFacCF2NTLZZX (F);
525 ZZX ZZg= convertFacCF2NTLZZX (G);
526 ZZ_pX NTLf= to_ZZ_pX (ZZf);
527 ZZ_pX NTLg= to_ZZ_pX (ZZg);
528 mul (NTLf, NTLf, NTLg);
529 return b (convertNTLZZX2CF (to_ZZX (NTLf), F.mvar()));
530 }
531 return F*G;
532#endif
533 }
534 if (b.getp() != 0)
535 {
536 if (!F.inBaseDomain() && !G.inBaseDomain())
537 {
539 {
540#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
541 fmpz_t FLINTp;
542 fmpz_mod_poly_t FLINTmipo;
543 fq_ctx_t fq_con;
544
545 fmpz_init (FLINTp);
546 convertCF2initFmpz (FLINTp, b.getpk());
547
549 bool rat=isOn(SW_RATIONAL);
552 mipo *= cd;
553 if (!rat) Off(SW_RATIONAL);
554 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, mipo, FLINTp);
555
556 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
557 fmpz_mod_ctx_t fmpz_ctx;
558 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
559 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
560 #else
561 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
562 #endif
563
565
566 if (F.inCoeffDomain() && !G.inCoeffDomain())
567 {
568 fq_poly_t FLINTG;
569 fmpz_poly_t FLINTF;
570 convertFacCF2Fmpz_poly_t (FLINTF, F);
572
573 fq_poly_scalar_mul_fq (FLINTG, FLINTG, FLINTF, fq_con);
574
575 result= convertFq_poly_t2FacCF (FLINTG, G.mvar(), alpha, fq_con);
576 fmpz_poly_clear (FLINTF);
577 fq_poly_clear (FLINTG, fq_con);
578 }
579 else if (!F.inCoeffDomain() && G.inCoeffDomain())
580 {
581 fq_poly_t FLINTF;
582 fmpz_poly_t FLINTG;
583
584 convertFacCF2Fmpz_poly_t (FLINTG, G);
585 convertFacCF2Fq_poly_t (FLINTF, F, fq_con);
586
587 fq_poly_scalar_mul_fq (FLINTF, FLINTF, FLINTG, fq_con);
588
590 fmpz_poly_clear (FLINTG);
591 fq_poly_clear (FLINTF, fq_con);
592 }
593 else
594 {
595 fq_t FLINTF, FLINTG;
596
597 convertFacCF2Fq_t (FLINTF, F, fq_con);
598 convertFacCF2Fq_t (FLINTG, G, fq_con);
599
600 fq_mul (FLINTF, FLINTF, FLINTG, fq_con);
601
602 result= convertFq_t2FacCF (FLINTF, alpha);
603 fq_clear (FLINTF, fq_con);
604 fq_clear (FLINTG, fq_con);
605 }
606
607 fmpz_clear (FLINTp);
608 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
609 fmpz_mod_poly_clear (FLINTmipo,fmpz_ctx);
610 fmpz_mod_ctx_clear(fmpz_ctx);
611 #else
612 fmpz_mod_poly_clear (FLINTmipo);
613 #endif
614 fq_ctx_clear (fq_con);
615
616 return b (result);
617#endif
618#ifdef HAVE_NTL
619 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
620 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (getMipo (alpha)));
621 ZZ_pE::init (NTLmipo);
622
623 if (F.inCoeffDomain() && !G.inCoeffDomain())
624 {
625 ZZ_pEX NTLg= convertFacCF2NTLZZ_pEX (G, NTLmipo);
626 ZZ_pX NTLf= convertFacCF2NTLZZpX (F);
627 mul (NTLg, to_ZZ_pE (NTLf), NTLg);
628 return b (convertNTLZZ_pEX2CF (NTLg, G.mvar(), alpha));
629 }
630 else if (!F.inCoeffDomain() && G.inCoeffDomain())
631 {
632 ZZ_pX NTLg= convertFacCF2NTLZZpX (G);
633 ZZ_pEX NTLf= convertFacCF2NTLZZ_pEX (F, NTLmipo);
634 mul (NTLf, NTLf, to_ZZ_pE (NTLg));
635 return b (convertNTLZZ_pEX2CF (NTLf, F.mvar(), alpha));
636 }
637 else
638 {
639 ZZ_pX NTLg= convertFacCF2NTLZZpX (G);
640 ZZ_pX NTLf= convertFacCF2NTLZZpX (F);
641 ZZ_pE result;
642 mul (result, to_ZZ_pE (NTLg), to_ZZ_pE (NTLf));
643 return b (convertNTLZZpX2CF (rep (result), alpha));
644 }
645#endif
646 }
647 }
648 return b (F*G);
649 }
650 return F*G;
651 }
652 else if (F.inCoeffDomain() || G.inCoeffDomain())
653 return F*G;
654 ASSERT (F.isUnivariate() && G.isUnivariate(), "expected univariate polys");
655 ASSERT (F.level() == G.level(), "expected polys of same level");
656#ifdef HAVE_NTL
657#if (!defined(HAVE_FLINT) || __FLINT_RELEASE < 20400)
659 {
661 zz_p::init (getCharacteristic());
662 }
663#endif
664#endif
668 {
669 if (!getReduce (alpha))
670 {
671 result= 0;
672 for (CFIterator i= F; i.hasTerms(); i++)
673 result += i.coeff()*G*power (F.mvar(),i.exp());
674 return result;
675 }
676#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
677 nmod_poly_t FLINTmipo;
678 fq_nmod_ctx_t fq_con;
679
680 nmod_poly_init (FLINTmipo, getCharacteristic());
682
683 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
684
685 fq_nmod_poly_t FLINTF, FLINTG;
688
689 fq_nmod_poly_mul (FLINTF, FLINTF, FLINTG, fq_con);
690
692
693 fq_nmod_poly_clear (FLINTF, fq_con);
694 fq_nmod_poly_clear (FLINTG, fq_con);
695 nmod_poly_clear (FLINTmipo);
697 return result;
698#elif defined(AHVE_NTL)
699 zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
700 zz_pE::init (NTLMipo);
701 zz_pEX NTLF= convertFacCF2NTLzz_pEX (F, NTLMipo);
702 zz_pEX NTLG= convertFacCF2NTLzz_pEX (G, NTLMipo);
703 mul (NTLF, NTLF, NTLG);
705 return result;
706#endif
707 }
708 else
709 {
710#ifdef HAVE_FLINT
711 nmod_poly_t FLINTF, FLINTG;
712 convertFacCF2nmod_poly_t (FLINTF, F);
713 convertFacCF2nmod_poly_t (FLINTG, G);
714 nmod_poly_mul (FLINTF, FLINTF, FLINTG);
715 result= convertnmod_poly_t2FacCF (FLINTF, F.mvar());
716 nmod_poly_clear (FLINTF);
717 nmod_poly_clear (FLINTG);
718 return result;
719#endif
720#ifdef HAVE_NTL
721 zz_pX NTLF= convertFacCF2NTLzzpX (F);
722 zz_pX NTLG= convertFacCF2NTLzzpX (G);
723 mul (NTLF, NTLF, NTLG);
724 return convertNTLzzpX2CF(NTLF, F.mvar());
725#endif
726 }
727 return F*G;
728}
729
731modNTL (const CanonicalForm& F, const CanonicalForm& G, const modpk& b)
732{
734 return mod (F, G);
735 if (F.inCoeffDomain() && G.isUnivariate() && !G.inCoeffDomain())
736 {
737 if (b.getp() != 0)
738 return b(F);
739 return F;
740 }
741 else if (F.inCoeffDomain() && G.inCoeffDomain())
742 {
743 if (b.getp() != 0)
744 return b(F%G);
745 return mod (F, G);
746 }
747 else if (F.isUnivariate() && G.inCoeffDomain())
748 {
749 if (b.getp() != 0)
750 return b(F%G);
751 return mod (F,G);
752 }
753
754 if (getCharacteristic() == 0)
755 {
757 if (!hasFirstAlgVar (F, alpha) && !hasFirstAlgVar (G, alpha))
758 {
759#ifdef HAVE_FLINT
760 if (b.getp() != 0)
761 {
762 fmpz_t FLINTpk;
763 fmpz_init (FLINTpk);
764 convertCF2initFmpz (FLINTpk, b.getpk());
765 fmpz_mod_poly_t FLINTF, FLINTG;
766 convertFacCF2Fmpz_mod_poly_t (FLINTF, F, FLINTpk);
767 convertFacCF2Fmpz_mod_poly_t (FLINTG, G, FLINTpk);
768 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
769 fmpz_mod_ctx_t fmpz_ctx;
770 fmpz_mod_ctx_init(fmpz_ctx,FLINTpk);
771 fmpz_mod_poly_divrem (FLINTG, FLINTF, FLINTF, FLINTG, fmpz_ctx);
772 #else
773 fmpz_mod_poly_divrem (FLINTG, FLINTF, FLINTF, FLINTG);
774 #endif
776 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
777 fmpz_mod_poly_clear (FLINTG, fmpz_ctx);
778 fmpz_mod_poly_clear (FLINTF, fmpz_ctx);
779 fmpz_mod_ctx_clear(fmpz_ctx);
780 #else
781 fmpz_mod_poly_clear (FLINTG);
782 fmpz_mod_poly_clear (FLINTF);
783 #endif
784 fmpz_clear (FLINTpk);
785 return result;
786 }
787 return modFLINTQ (F, G);
788#else
789 if (b.getp() != 0)
790 {
791 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
792 ZZX ZZf= convertFacCF2NTLZZX (F);
793 ZZX ZZg= convertFacCF2NTLZZX (G);
794 ZZ_pX NTLf= to_ZZ_pX (ZZf);
795 ZZ_pX NTLg= to_ZZ_pX (ZZg);
796 rem (NTLf, NTLf, NTLg);
797 return b (convertNTLZZX2CF (to_ZZX (NTLf), F.mvar()));
798 }
799 return mod (F, G);
800#endif
801 }
802 else
803 {
804 if (b.getp() != 0)
805 {
806#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
807 fmpz_t FLINTp;
808 fmpz_mod_poly_t FLINTmipo;
809 fq_ctx_t fq_con;
810 fq_poly_t FLINTF, FLINTG;
811
812 fmpz_init (FLINTp);
813
814 convertCF2initFmpz (FLINTp, b.getpk());
815
817 bool rat=isOn(SW_RATIONAL);
820 mipo *= cd;
821 if (!rat) Off(SW_RATIONAL);
822 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, mipo, FLINTp);
823
824 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
825 fmpz_mod_ctx_t fmpz_ctx;
826 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
827 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
828 #else
829 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
830 #endif
831
832 convertFacCF2Fq_poly_t (FLINTF, F, fq_con);
834
835 fq_poly_rem (FLINTF, FLINTF, FLINTG, fq_con);
836
838 alpha, fq_con);
839
840 fmpz_clear (FLINTp);
841 fq_poly_clear (FLINTF, fq_con);
842 fq_poly_clear (FLINTG, fq_con);
843 fq_ctx_clear (fq_con);
844 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
845 fmpz_mod_poly_clear (FLINTmipo, fmpz_ctx);
846 fmpz_mod_ctx_clear(fmpz_ctx);
847 #else
848 fmpz_mod_poly_clear (FLINTmipo);
849 #endif
850
851 return b(result);
852#else
853 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
854 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (getMipo (alpha)));
855 ZZ_pE::init (NTLmipo);
856 ZZ_pEX NTLg= convertFacCF2NTLZZ_pEX (G, NTLmipo);
857 ZZ_pEX NTLf= convertFacCF2NTLZZ_pEX (F, NTLmipo);
858 rem (NTLf, NTLf, NTLg);
859 return b (convertNTLZZ_pEX2CF (NTLf, F.mvar(), alpha));
860#endif
861 }
862#ifdef HAVE_FLINT
864 newtonDivrem (F, G, Q, R);
865 return R;
866#else
867 return mod (F,G);
868#endif
869 }
870 }
871
872 ASSERT (F.isUnivariate() && G.isUnivariate(), "expected univariate polys");
873 ASSERT (F.level() == G.level(), "expected polys of same level");
874#if (!defined(HAVE_FLINT) || __FLINT_RELEASE < 20400)
876 {
878 zz_p::init (getCharacteristic());
879 }
880#endif
884 {
885#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
886 nmod_poly_t FLINTmipo;
887 fq_nmod_ctx_t fq_con;
888
889 nmod_poly_init (FLINTmipo, getCharacteristic());
891
892 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
893
894 fq_nmod_poly_t FLINTF, FLINTG;
897
898 fq_nmod_poly_rem (FLINTF, FLINTF, FLINTG, fq_con);
899
901
902 fq_nmod_poly_clear (FLINTF, fq_con);
903 fq_nmod_poly_clear (FLINTG, fq_con);
904 nmod_poly_clear (FLINTmipo);
906#else
907 zz_pX NTLMipo= convertFacCF2NTLzzpX(getMipo (alpha));
908 zz_pE::init (NTLMipo);
909 zz_pEX NTLF= convertFacCF2NTLzz_pEX (F, NTLMipo);
910 zz_pEX NTLG= convertFacCF2NTLzz_pEX (G, NTLMipo);
911 rem (NTLF, NTLF, NTLG);
913#endif
914 }
915 else
916 {
917#ifdef HAVE_FLINT
918 nmod_poly_t FLINTF, FLINTG;
919 convertFacCF2nmod_poly_t (FLINTF, F);
920 convertFacCF2nmod_poly_t (FLINTG, G);
921 nmod_poly_divrem (FLINTG, FLINTF, FLINTF, FLINTG);
922 result= convertnmod_poly_t2FacCF (FLINTF, F.mvar());
923 nmod_poly_clear (FLINTF);
924 nmod_poly_clear (FLINTG);
925#else
926 zz_pX NTLF= convertFacCF2NTLzzpX (F);
927 zz_pX NTLG= convertFacCF2NTLzzpX (G);
928 rem (NTLF, NTLF, NTLG);
929 result= convertNTLzzpX2CF(NTLF, F.mvar());
930#endif
931 }
932 return result;
933}
934
936divNTL (const CanonicalForm& F, const CanonicalForm& G, const modpk& b)
937{
939 return div (F, G);
940 if (F.inCoeffDomain() && G.isUnivariate() && !G.inCoeffDomain())
941 {
942 return 0;
943 }
944 else if (F.inCoeffDomain() && G.inCoeffDomain())
945 {
946 if (b.getp() != 0)
947 {
948 if (!F.inBaseDomain() || !G.inBaseDomain())
949 {
953#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
954 fmpz_t FLINTp;
955 fmpz_mod_poly_t FLINTmipo;
956 fq_ctx_t fq_con;
957 fq_t FLINTF, FLINTG;
958
959 fmpz_init (FLINTp);
960 convertCF2initFmpz (FLINTp, b.getpk());
961
962 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, getMipo (alpha), FLINTp);
963
964 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
965 fmpz_mod_ctx_t fmpz_ctx;
966 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
967 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
968 #else
969 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
970 #endif
971
972 convertFacCF2Fq_t (FLINTF, F, fq_con);
973 convertFacCF2Fq_t (FLINTG, G, fq_con);
974
975 fq_inv (FLINTG, FLINTG, fq_con);
976 fq_mul (FLINTF, FLINTF, FLINTG, fq_con);
977
979
980 fmpz_clear (FLINTp);
981 fq_clear (FLINTF, fq_con);
982 fq_clear (FLINTG, fq_con);
983 fq_ctx_clear (fq_con);
984 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
985 fmpz_mod_poly_clear (FLINTmipo, fmpz_ctx);
986 fmpz_mod_ctx_clear(fmpz_ctx);
987 #else
988 fmpz_mod_poly_clear (FLINTmipo);
989 #endif
990 return b (result);
991#else
992 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
993 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (getMipo (alpha)));
994 ZZ_pE::init (NTLmipo);
995 ZZ_pX NTLg= convertFacCF2NTLZZpX (G);
996 ZZ_pX NTLf= convertFacCF2NTLZZpX (F);
997 ZZ_pE result;
998 div (result, to_ZZ_pE (NTLf), to_ZZ_pE (NTLg));
999 return b (convertNTLZZpX2CF (rep (result), alpha));
1000#endif
1001 }
1002 return b(div (F,G));
1003 }
1004 return div (F, G);
1005 }
1006 else if (F.isUnivariate() && G.inCoeffDomain())
1007 {
1008 if (b.getp() != 0)
1009 {
1010 if (!G.inBaseDomain())
1011 {
1014#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
1015 fmpz_t FLINTp;
1016 fmpz_mod_poly_t FLINTmipo;
1017 fq_ctx_t fq_con;
1018 fq_poly_t FLINTF;
1019 fq_t FLINTG;
1020
1021 fmpz_init (FLINTp);
1022 convertCF2initFmpz (FLINTp, b.getpk());
1023
1024 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, getMipo (alpha), FLINTp);
1025
1026 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1027 fmpz_mod_ctx_t fmpz_ctx;
1028 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
1029 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
1030 #else
1031 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
1032 #endif
1033
1034 convertFacCF2Fq_poly_t (FLINTF, F, fq_con);
1035 convertFacCF2Fq_t (FLINTG, G, fq_con);
1036
1037 fq_inv (FLINTG, FLINTG, fq_con);
1038 fq_poly_scalar_mul_fq (FLINTF, FLINTF, FLINTG, fq_con);
1039
1041 alpha, fq_con);
1042
1043 fmpz_clear (FLINTp);
1044 fq_poly_clear (FLINTF, fq_con);
1045 fq_clear (FLINTG, fq_con);
1046 fq_ctx_clear (fq_con);
1047 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1048 fmpz_mod_poly_clear (FLINTmipo, fmpz_ctx);
1049 fmpz_mod_ctx_clear(fmpz_ctx);
1050 #else
1051 fmpz_mod_poly_clear (FLINTmipo);
1052 #endif
1053 return b (result);
1054#else
1055 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
1056 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (getMipo (alpha)));
1057 ZZ_pE::init (NTLmipo);
1058 ZZ_pX NTLg= convertFacCF2NTLZZpX (G);
1059 ZZ_pEX NTLf= convertFacCF2NTLZZ_pEX (F, NTLmipo);
1060 div (NTLf, NTLf, to_ZZ_pE (NTLg));
1061 return b (convertNTLZZ_pEX2CF (NTLf, F.mvar(), alpha));
1062#endif
1063 }
1064 return b(div (F,G));
1065 }
1066 return div (F, G);
1067 }
1068
1069 if (getCharacteristic() == 0)
1070 {
1071
1073 if (!hasFirstAlgVar (F, alpha) && !hasFirstAlgVar (G, alpha))
1074 {
1075#ifdef HAVE_FLINT
1076 if (b.getp() != 0)
1077 {
1078 fmpz_t FLINTpk;
1079 fmpz_init (FLINTpk);
1080 convertCF2initFmpz (FLINTpk, b.getpk());
1081 fmpz_mod_poly_t FLINTF, FLINTG;
1082 convertFacCF2Fmpz_mod_poly_t (FLINTF, F, FLINTpk);
1083 convertFacCF2Fmpz_mod_poly_t (FLINTG, G, FLINTpk);
1084 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1085 fmpz_mod_ctx_t fmpz_ctx;
1086 fmpz_mod_ctx_init(fmpz_ctx,FLINTpk);
1087 fmpz_mod_poly_divrem (FLINTF, FLINTG, FLINTF, FLINTG, fmpz_ctx);
1088 #else
1089 fmpz_mod_poly_divrem (FLINTF, FLINTG, FLINTF, FLINTG);
1090 #endif
1092 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1093 fmpz_mod_poly_clear (FLINTG, fmpz_ctx);
1094 fmpz_mod_poly_clear (FLINTF, fmpz_ctx);
1095 fmpz_mod_ctx_clear(fmpz_ctx);
1096 #else
1097 fmpz_mod_poly_clear (FLINTG);
1098 fmpz_mod_poly_clear (FLINTF);
1099 #endif
1100 fmpz_clear (FLINTpk);
1101 return result;
1102 }
1103 return divFLINTQ (F,G);
1104#else
1105 if (b.getp() != 0)
1106 {
1107 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
1108 ZZX ZZf= convertFacCF2NTLZZX (F);
1109 ZZX ZZg= convertFacCF2NTLZZX (G);
1110 ZZ_pX NTLf= to_ZZ_pX (ZZf);
1111 ZZ_pX NTLg= to_ZZ_pX (ZZg);
1112 div (NTLf, NTLf, NTLg);
1113 return b (convertNTLZZX2CF (to_ZZX (NTLf), F.mvar()));
1114 }
1115 return div (F, G);
1116#endif
1117 }
1118 else
1119 {
1120 if (b.getp() != 0)
1121 {
1122#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
1123 fmpz_t FLINTp;
1124 fmpz_mod_poly_t FLINTmipo;
1125 fq_ctx_t fq_con;
1126 fq_poly_t FLINTF, FLINTG;
1127
1128 fmpz_init (FLINTp);
1129 convertCF2initFmpz (FLINTp, b.getpk());
1130
1131 convertFacCF2Fmpz_mod_poly_t (FLINTmipo, getMipo (alpha), FLINTp);
1132
1133 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1134 fmpz_mod_ctx_t fmpz_ctx;
1135 fmpz_mod_ctx_init(fmpz_ctx,FLINTp);
1136 fq_ctx_init_modulus (fq_con, FLINTmipo, fmpz_ctx, "Z");
1137 #else
1138 fq_ctx_init_modulus (fq_con, FLINTmipo, "Z");
1139 #endif
1140
1141 convertFacCF2Fq_poly_t (FLINTF, F, fq_con);
1142 convertFacCF2Fq_poly_t (FLINTG, G, fq_con);
1143
1144 fq_poly_divrem (FLINTF, FLINTG, FLINTF, FLINTG, fq_con);
1145
1147 alpha, fq_con);
1148
1149 fmpz_clear (FLINTp);
1150 fq_poly_clear (FLINTF, fq_con);
1151 fq_poly_clear (FLINTG, fq_con);
1152 fq_ctx_clear (fq_con);
1153 #if (HAVE_FLINT && __FLINT_RELEASE >= 20700)
1154 fmpz_mod_poly_clear (FLINTmipo, fmpz_ctx);
1155 fmpz_mod_ctx_clear(fmpz_ctx);
1156 #else
1157 fmpz_mod_poly_clear (FLINTmipo);
1158 #endif
1159 return b (result);
1160#else
1161 ZZ_p::init (convertFacCF2NTLZZ (b.getpk()));
1162 ZZ_pX NTLmipo= to_ZZ_pX (convertFacCF2NTLZZX (getMipo (alpha)));
1163 ZZ_pE::init (NTLmipo);
1164 ZZ_pEX NTLg= convertFacCF2NTLZZ_pEX (G, NTLmipo);
1165 ZZ_pEX NTLf= convertFacCF2NTLZZ_pEX (F, NTLmipo);
1166 div (NTLf, NTLf, NTLg);
1167 return b (convertNTLZZ_pEX2CF (NTLf, F.mvar(), alpha));
1168#endif
1169 }
1170#ifdef HAVE_FLINT
1172 newtonDiv (F, G, Q);
1173 return Q;
1174#else
1175 return div (F,G);
1176#endif
1177 }
1178 }
1179
1180 ASSERT (F.isUnivariate() && G.isUnivariate(), "expected univariate polys");
1181 ASSERT (F.level() == G.level(), "expected polys of same level");
1182#if (!defined(HAVE_FLINT) || __FLINT_RELEASE < 20400)
1184 {
1186 zz_p::init (getCharacteristic());
1187 }
1188#endif
1191 if (hasFirstAlgVar (F, alpha) || hasFirstAlgVar (G, alpha))
1192 {
1193#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
1194 nmod_poly_t FLINTmipo;
1195 fq_nmod_ctx_t fq_con;
1196
1197 nmod_poly_init (FLINTmipo, getCharacteristic());
1198 convertFacCF2nmod_poly_t (FLINTmipo, getMipo (alpha));
1199
1200 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
1201
1202 fq_nmod_poly_t FLINTF, FLINTG;
1205
1206 fq_nmod_poly_divrem (FLINTF, FLINTG, FLINTF, FLINTG, fq_con);
1207
1209
1210 fq_nmod_poly_clear (FLINTF, fq_con);
1211 fq_nmod_poly_clear (FLINTG, fq_con);
1212 nmod_poly_clear (FLINTmipo);
1214#else
1215 zz_pX NTLMipo= convertFacCF2NTLzzpX(getMipo (alpha));
1216 zz_pE::init (NTLMipo);
1217 zz_pEX NTLF= convertFacCF2NTLzz_pEX (F, NTLMipo);
1218 zz_pEX NTLG= convertFacCF2NTLzz_pEX (G, NTLMipo);
1219 div (NTLF, NTLF, NTLG);
1220 result= convertNTLzz_pEX2CF(NTLF, F.mvar(), alpha);
1221#endif
1222 }
1223 else
1224 {
1225#ifdef HAVE_FLINT
1226 nmod_poly_t FLINTF, FLINTG;
1227 convertFacCF2nmod_poly_t (FLINTF, F);
1228 convertFacCF2nmod_poly_t (FLINTG, G);
1229 nmod_poly_div (FLINTF, FLINTF, FLINTG);
1230 result= convertnmod_poly_t2FacCF (FLINTF, F.mvar());
1231 nmod_poly_clear (FLINTF);
1232 nmod_poly_clear (FLINTG);
1233#else
1234 zz_pX NTLF= convertFacCF2NTLzzpX (F);
1235 zz_pX NTLG= convertFacCF2NTLzzpX (G);
1236 div (NTLF, NTLF, NTLG);
1237 result= convertNTLzzpX2CF(NTLF, F.mvar());
1238#endif
1239 }
1240 return result;
1241}
1242
1243// end univariate polys
1244//*************************
1245// bivariate polys
1246
1247#ifdef HAVE_FLINT
1248void kronSubFp (nmod_poly_t result, const CanonicalForm& A, int d)
1249{
1250 int degAy= degree (A);
1251 nmod_poly_init2 (result, getCharacteristic(), d*(degAy + 1));
1252 result->length= d*(degAy + 1);
1253 flint_mpn_zero (result->coeffs, d*(degAy+1));
1254
1255 nmod_poly_t buf;
1256
1257 int k;
1258 for (CFIterator i= A; i.hasTerms(); i++)
1259 {
1260 convertFacCF2nmod_poly_t (buf, i.coeff());
1261 k= i.exp()*d;
1262 flint_mpn_copyi (result->coeffs+k, buf->coeffs, nmod_poly_length(buf));
1263
1265 }
1266 _nmod_poly_normalise (result);
1267}
1268
1269#if ( __FLINT_RELEASE >= 20400)
1270void
1271kronSubFq (fq_nmod_poly_t result, const CanonicalForm& A, int d,
1272 const fq_nmod_ctx_t fq_con)
1273{
1274 int degAy= degree (A);
1275 fq_nmod_poly_init2 (result, d*(degAy + 1), fq_con);
1276 _fq_nmod_poly_set_length (result, d*(degAy + 1), fq_con);
1277 _fq_nmod_vec_zero (result->coeffs, d*(degAy + 1), fq_con);
1278
1279 fq_nmod_poly_t buf1;
1280
1281 nmod_poly_t buf2;
1282
1283 int k;
1284
1285 for (CFIterator i= A; i.hasTerms(); i++)
1286 {
1287 if (i.coeff().inCoeffDomain())
1288 {
1289 convertFacCF2nmod_poly_t (buf2, i.coeff());
1290 fq_nmod_poly_init2 (buf1, 1, fq_con);
1291 fq_nmod_poly_set_coeff (buf1, 0, buf2, fq_con);
1293 }
1294 else
1296
1297 k= i.exp()*d;
1298 _fq_nmod_vec_set (result->coeffs+k, buf1->coeffs,
1299 fq_nmod_poly_length (buf1, fq_con), fq_con);
1300
1302 }
1303
1304 _fq_nmod_poly_normalise (result, fq_con);
1305}
1306#endif
1307
1308/*void kronSubQa (fmpq_poly_t result, const CanonicalForm& A, int d1, int d2)
1309{
1310 int degAy= degree (A);
1311 fmpq_poly_init2 (result, d1*(degAy + 1));
1312
1313 fmpq_poly_t buf;
1314 fmpq_t coeff;
1315
1316 int k, l, bufRepLength;
1317 CFIterator j;
1318 for (CFIterator i= A; i.hasTerms(); i++)
1319 {
1320 if (i.coeff().inCoeffDomain())
1321 {
1322 k= d1*i.exp();
1323 convertFacCF2Fmpq_poly_t (buf, i.coeff());
1324 bufRepLength= (int) fmpq_poly_length(buf);
1325 for (l= 0; l < bufRepLength; l++)
1326 {
1327 fmpq_poly_get_coeff_fmpq (coeff, buf, l);
1328 fmpq_poly_set_coeff_fmpq (result, l + k, coeff);
1329 }
1330 fmpq_poly_clear (buf);
1331 }
1332 else
1333 {
1334 for (j= i.coeff(); j.hasTerms(); j++)
1335 {
1336 k= d1*i.exp();
1337 k += d2*j.exp();
1338 convertFacCF2Fmpq_poly_t (buf, j.coeff());
1339 bufRepLength= (int) fmpq_poly_length(buf);
1340 for (l= 0; l < bufRepLength; l++)
1341 {
1342 fmpq_poly_get_coeff_fmpq (coeff, buf, l);
1343 fmpq_poly_set_coeff_fmpq (result, k + l, coeff);
1344 }
1345 fmpq_poly_clear (buf);
1346 }
1347 }
1348 }
1349 fmpq_clear (coeff);
1350 _fmpq_poly_normalise (result);
1351}*/
1352
1353void kronSubQa (fmpz_poly_t result, const CanonicalForm& A, int d1, int d2)
1354{
1355 int degAy= degree (A);
1356 fmpz_poly_init2 (result, d1*(degAy + 1));
1357 _fmpz_poly_set_length (result, d1*(degAy + 1));
1358
1359 fmpz_poly_t buf;
1360
1361 int k;
1362 CFIterator j;
1363 for (CFIterator i= A; i.hasTerms(); i++)
1364 {
1365 if (i.coeff().inCoeffDomain())
1366 {
1367 k= d1*i.exp();
1368 convertFacCF2Fmpz_poly_t (buf, i.coeff());
1369 _fmpz_vec_set (result->coeffs + k, buf->coeffs, buf->length);
1370 fmpz_poly_clear (buf);
1371 }
1372 else
1373 {
1374 for (j= i.coeff(); j.hasTerms(); j++)
1375 {
1376 k= d1*i.exp();
1377 k += d2*j.exp();
1378 convertFacCF2Fmpz_poly_t (buf, j.coeff());
1379 _fmpz_vec_set (result->coeffs + k, buf->coeffs, buf->length);
1380 fmpz_poly_clear (buf);
1381 }
1382 }
1383 }
1384 _fmpz_poly_normalise (result);
1385}
1386
1387void
1388kronSubReciproFp (nmod_poly_t subA1, nmod_poly_t subA2, const CanonicalForm& A,
1389 int d)
1390{
1391 int degAy= degree (A);
1392 mp_limb_t ninv= n_preinvert_limb (getCharacteristic());
1393 nmod_poly_init2_preinv (subA1, getCharacteristic(), ninv, d*(degAy + 2));
1394 nmod_poly_init2_preinv (subA2, getCharacteristic(), ninv, d*(degAy + 2));
1395
1396 nmod_poly_t buf;
1397
1398 int k, kk, j, bufRepLength;
1399 for (CFIterator i= A; i.hasTerms(); i++)
1400 {
1401 convertFacCF2nmod_poly_t (buf, i.coeff());
1402
1403 k= i.exp()*d;
1404 kk= (degAy - i.exp())*d;
1405 bufRepLength= (int) nmod_poly_length (buf);
1406 for (j= 0; j < bufRepLength; j++)
1407 {
1408 nmod_poly_set_coeff_ui (subA1, j + k,
1409 n_addmod (nmod_poly_get_coeff_ui (subA1, j+k),
1410 nmod_poly_get_coeff_ui (buf, j),
1412 )
1413 );
1414 nmod_poly_set_coeff_ui (subA2, j + kk,
1415 n_addmod (nmod_poly_get_coeff_ui (subA2, j + kk),
1416 nmod_poly_get_coeff_ui (buf, j),
1418 )
1419 );
1420 }
1422 }
1423 _nmod_poly_normalise (subA1);
1424 _nmod_poly_normalise (subA2);
1425}
1426
1427#if ( __FLINT_RELEASE >= 20400)
1428void
1429kronSubReciproFq (fq_nmod_poly_t subA1, fq_nmod_poly_t subA2,
1430 const CanonicalForm& A, int d, const fq_nmod_ctx_t fq_con)
1431{
1432 int degAy= degree (A);
1433 fq_nmod_poly_init2 (subA1, d*(degAy + 2), fq_con);
1434 fq_nmod_poly_init2 (subA2, d*(degAy + 2), fq_con);
1435
1436 _fq_nmod_poly_set_length (subA1, d*(degAy + 2), fq_con);
1437 _fq_nmod_vec_zero (subA1->coeffs, d*(degAy + 2), fq_con);
1438
1439 _fq_nmod_poly_set_length (subA2, d*(degAy + 2), fq_con);
1440 _fq_nmod_vec_zero (subA2->coeffs, d*(degAy + 2), fq_con);
1441
1442 fq_nmod_poly_t buf1;
1443
1444 nmod_poly_t buf2;
1445
1446 int k, kk;
1447 for (CFIterator i= A; i.hasTerms(); i++)
1448 {
1449 if (i.coeff().inCoeffDomain())
1450 {
1451 convertFacCF2nmod_poly_t (buf2, i.coeff());
1452 fq_nmod_poly_init2 (buf1, 1, fq_con);
1453 fq_nmod_poly_set_coeff (buf1, 0, buf2, fq_con);
1455 }
1456 else
1458
1459 k= i.exp()*d;
1460 kk= (degAy - i.exp())*d;
1461 _fq_nmod_vec_add (subA1->coeffs+k, subA1->coeffs+k, buf1->coeffs,
1462 fq_nmod_poly_length(buf1, fq_con), fq_con);
1463 _fq_nmod_vec_add (subA2->coeffs+kk, subA2->coeffs+kk, buf1->coeffs,
1464 fq_nmod_poly_length(buf1, fq_con), fq_con);
1465
1467 }
1468 _fq_nmod_poly_normalise (subA1, fq_con);
1469 _fq_nmod_poly_normalise (subA2, fq_con);
1470}
1471#endif
1472
1473void
1474kronSubReciproQ (fmpz_poly_t subA1, fmpz_poly_t subA2, const CanonicalForm& A,
1475 int d)
1476{
1477 int degAy= degree (A);
1478 fmpz_poly_init2 (subA1, d*(degAy + 2));
1479 fmpz_poly_init2 (subA2, d*(degAy + 2));
1480
1481 fmpz_poly_t buf;
1482
1483 int k, kk;
1484 for (CFIterator i= A; i.hasTerms(); i++)
1485 {
1486 convertFacCF2Fmpz_poly_t (buf, i.coeff());
1487
1488 k= i.exp()*d;
1489 kk= (degAy - i.exp())*d;
1490 _fmpz_vec_add (subA1->coeffs+k, subA1->coeffs + k, buf->coeffs, buf->length);
1491 _fmpz_vec_add (subA2->coeffs+kk, subA2->coeffs + kk, buf->coeffs, buf->length);
1492 fmpz_poly_clear (buf);
1493 }
1494
1495 _fmpz_poly_normalise (subA1);
1496 _fmpz_poly_normalise (subA2);
1497}
1498
1499CanonicalForm reverseSubstQ (const fmpz_poly_t F, int d)
1500{
1501 Variable y= Variable (2);
1502 Variable x= Variable (1);
1503
1504 fmpz_poly_t buf;
1506 int i= 0;
1507 int degf= fmpz_poly_degree(F);
1508 int k= 0;
1509 int degfSubK, repLength;
1510 while (degf >= k)
1511 {
1512 degfSubK= degf - k;
1513 if (degfSubK >= d)
1514 repLength= d;
1515 else
1516 repLength= degfSubK + 1;
1517
1518 fmpz_poly_init2 (buf, repLength);
1519 _fmpz_poly_set_length (buf, repLength);
1520 _fmpz_vec_set (buf->coeffs, F->coeffs+k, repLength);
1521 _fmpz_poly_normalise (buf);
1522
1524 i++;
1525 k= d*i;
1526 fmpz_poly_clear (buf);
1527 }
1528
1529 return result;
1530}
1531
1532/*CanonicalForm
1533reverseSubstQa (const fmpq_poly_t F, int d1, int d2, const Variable& alpha,
1534 const fmpq_poly_t mipo)
1535{
1536 Variable y= Variable (2);
1537 Variable x= Variable (1);
1538
1539 fmpq_poly_t f;
1540 fmpq_poly_init (f);
1541 fmpq_poly_set (f, F);
1542
1543 fmpq_poly_t buf;
1544 CanonicalForm result= 0, result2;
1545 int i= 0;
1546 int degf= fmpq_poly_degree(f);
1547 int k= 0;
1548 int degfSubK;
1549 int repLength;
1550 fmpq_t coeff;
1551 while (degf >= k)
1552 {
1553 degfSubK= degf - k;
1554 if (degfSubK >= d1)
1555 repLength= d1;
1556 else
1557 repLength= degfSubK + 1;
1558
1559 fmpq_init (coeff);
1560 int j= 0;
1561 int l;
1562 result2= 0;
1563 while (j*d2 < repLength)
1564 {
1565 fmpq_poly_init2 (buf, d2);
1566 for (l= 0; l < d2; l++)
1567 {
1568 fmpq_poly_get_coeff_fmpq (coeff, f, k + j*d2 + l);
1569 fmpq_poly_set_coeff_fmpq (buf, l, coeff);
1570 }
1571 _fmpq_poly_normalise (buf);
1572 fmpq_poly_rem (buf, buf, mipo);
1573 result2 += convertFmpq_poly_t2FacCF (buf, alpha)*power (x, j);
1574 j++;
1575 fmpq_poly_clear (buf);
1576 }
1577 if (repLength - j*d2 != 0 && j*d2 - repLength < d2)
1578 {
1579 j--;
1580 repLength -= j*d2;
1581 fmpq_poly_init2 (buf, repLength);
1582 j++;
1583 for (l= 0; l < repLength; l++)
1584 {
1585 fmpq_poly_get_coeff_fmpq (coeff, f, k + j*d2 + l);
1586 fmpq_poly_set_coeff_fmpq (buf, l, coeff);
1587 }
1588 _fmpq_poly_normalise (buf);
1589 fmpq_poly_rem (buf, buf, mipo);
1590 result2 += convertFmpq_poly_t2FacCF (buf, alpha)*power (x, j);
1591 fmpq_poly_clear (buf);
1592 }
1593 fmpq_clear (coeff);
1594
1595 result += result2*power (y, i);
1596 i++;
1597 k= d1*i;
1598 }
1599
1600 fmpq_poly_clear (f);
1601 return result;
1602}*/
1603
1605reverseSubstQa (const fmpz_poly_t F, int d1, int d2, const Variable& alpha,
1606 const fmpq_poly_t mipo)
1607{
1608 Variable y= Variable (2);
1609 Variable x= Variable (1);
1610
1611 fmpq_poly_t buf;
1612 CanonicalForm result= 0, result2;
1613 int i= 0;
1614 int degf= fmpz_poly_degree(F);
1615 int k= 0;
1616 int degfSubK;
1617 int repLength;
1618 while (degf >= k)
1619 {
1620 degfSubK= degf - k;
1621 if (degfSubK >= d1)
1622 repLength= d1;
1623 else
1624 repLength= degfSubK + 1;
1625
1626 int j= 0;
1627 result2= 0;
1628 while (j*d2 < repLength)
1629 {
1630 fmpq_poly_init2 (buf, d2);
1631 _fmpq_poly_set_length (buf, d2);
1632 _fmpz_vec_set (buf->coeffs, F->coeffs + k + j*d2, d2);
1633 _fmpq_poly_normalise (buf);
1634 fmpq_poly_rem (buf, buf, mipo);
1635 result2 += convertFmpq_poly_t2FacCF (buf, alpha)*power (x, j);
1636 j++;
1637 fmpq_poly_clear (buf);
1638 }
1639 if (repLength - j*d2 != 0 && j*d2 - repLength < d2)
1640 {
1641 j--;
1642 repLength -= j*d2;
1643 fmpq_poly_init2 (buf, repLength);
1644 _fmpq_poly_set_length (buf, repLength);
1645 j++;
1646 _fmpz_vec_set (buf->coeffs, F->coeffs + k + j*d2, repLength);
1647 _fmpq_poly_normalise (buf);
1648 fmpq_poly_rem (buf, buf, mipo);
1649 result2 += convertFmpq_poly_t2FacCF (buf, alpha)*power (x, j);
1650 fmpq_poly_clear (buf);
1651 }
1652
1653 result += result2*power (y, i);
1654 i++;
1655 k= d1*i;
1656 }
1657
1658 return result;
1659}
1660
1662reverseSubstReciproFp (const nmod_poly_t F, const nmod_poly_t G, int d, int k)
1663{
1664 Variable y= Variable (2);
1665 Variable x= Variable (1);
1666
1667 nmod_poly_t f, g;
1668 mp_limb_t ninv= n_preinvert_limb (getCharacteristic());
1669 nmod_poly_init_preinv (f, getCharacteristic(), ninv);
1670 nmod_poly_init_preinv (g, getCharacteristic(), ninv);
1671 nmod_poly_set (f, F);
1672 nmod_poly_set (g, G);
1673 int degf= nmod_poly_degree(f);
1674 int degg= nmod_poly_degree(g);
1675
1676
1677 nmod_poly_t buf1,buf2, buf3;
1678
1679 if (nmod_poly_length (f) < (long) d*(k+1)) //zero padding
1680 nmod_poly_fit_length (f,(long)d*(k+1));
1681
1683 int i= 0;
1684 int lf= 0;
1685 int lg= d*k;
1686 int degfSubLf= degf;
1687 int deggSubLg= degg-lg;
1688 int repLengthBuf2, repLengthBuf1, ind, tmp;
1689 while (degf >= lf || lg >= 0)
1690 {
1691 if (degfSubLf >= d)
1692 repLengthBuf1= d;
1693 else if (degfSubLf < 0)
1694 repLengthBuf1= 0;
1695 else
1696 repLengthBuf1= degfSubLf + 1;
1697 nmod_poly_init2_preinv (buf1, getCharacteristic(), ninv, repLengthBuf1);
1698
1699 for (ind= 0; ind < repLengthBuf1; ind++)
1700 nmod_poly_set_coeff_ui (buf1, ind, nmod_poly_get_coeff_ui (f, ind+lf));
1701 _nmod_poly_normalise (buf1);
1702
1703 repLengthBuf1= nmod_poly_length (buf1);
1704
1705 if (deggSubLg >= d - 1)
1706 repLengthBuf2= d - 1;
1707 else if (deggSubLg < 0)
1708 repLengthBuf2= 0;
1709 else
1710 repLengthBuf2= deggSubLg + 1;
1711
1712 nmod_poly_init2_preinv (buf2, getCharacteristic(), ninv, repLengthBuf2);
1713 for (ind= 0; ind < repLengthBuf2; ind++)
1714 nmod_poly_set_coeff_ui (buf2, ind, nmod_poly_get_coeff_ui (g, ind + lg));
1715
1716 _nmod_poly_normalise (buf2);
1717 repLengthBuf2= nmod_poly_length (buf2);
1718
1719 nmod_poly_init2_preinv (buf3, getCharacteristic(), ninv, repLengthBuf2 + d);
1720 for (ind= 0; ind < repLengthBuf1; ind++)
1721 nmod_poly_set_coeff_ui (buf3, ind, nmod_poly_get_coeff_ui (buf1, ind));
1722 for (ind= repLengthBuf1; ind < d; ind++)
1723 nmod_poly_set_coeff_ui (buf3, ind, 0);
1724 for (ind= 0; ind < repLengthBuf2; ind++)
1725 nmod_poly_set_coeff_ui (buf3, ind+d, nmod_poly_get_coeff_ui (buf2, ind));
1726 _nmod_poly_normalise (buf3);
1727
1728 result += convertnmod_poly_t2FacCF (buf3, x)*power (y, i);
1729 i++;
1730
1731
1732 lf= i*d;
1733 degfSubLf= degf - lf;
1734
1735 lg= d*(k-i);
1736 deggSubLg= degg - lg;
1737
1738 if (lg >= 0 && deggSubLg > 0)
1739 {
1740 if (repLengthBuf2 > degfSubLf + 1)
1741 degfSubLf= repLengthBuf2 - 1;
1742 tmp= tmin (repLengthBuf1, deggSubLg + 1);
1743 for (ind= 0; ind < tmp; ind++)
1744 nmod_poly_set_coeff_ui (g, ind + lg,
1745 n_submod (nmod_poly_get_coeff_ui (g, ind + lg),
1746 nmod_poly_get_coeff_ui (buf1, ind),
1748 )
1749 );
1750 }
1751 if (lg < 0)
1752 {
1755 nmod_poly_clear (buf3);
1756 break;
1757 }
1758 if (degfSubLf >= 0)
1759 {
1760 for (ind= 0; ind < repLengthBuf2; ind++)
1761 nmod_poly_set_coeff_ui (f, ind + lf,
1762 n_submod (nmod_poly_get_coeff_ui (f, ind + lf),
1763 nmod_poly_get_coeff_ui (buf2, ind),
1765 )
1766 );
1767 }
1770 nmod_poly_clear (buf3);
1771 }
1772
1775
1776 return result;
1777}
1778
1779#if ( __FLINT_RELEASE >= 20400)
1781reverseSubstReciproFq (const fq_nmod_poly_t F, const fq_nmod_poly_t G, int d,
1782 int k, const Variable& alpha, const fq_nmod_ctx_t fq_con)
1783{
1784 Variable y= Variable (2);
1785 Variable x= Variable (1);
1786
1787 fq_nmod_poly_t f, g;
1788 int degf= fq_nmod_poly_degree(F, fq_con);
1789 int degg= fq_nmod_poly_degree(G, fq_con);
1790
1791 fq_nmod_poly_t buf1,buf2, buf3;
1792
1795 fq_nmod_poly_set (f, F, fq_con);
1796 fq_nmod_poly_set (g, G, fq_con);
1797 if (fq_nmod_poly_length (f, fq_con) < (long) d*(k + 1)) //zero padding
1798 fq_nmod_poly_fit_length (f, (long) d*(k + 1), fq_con);
1799
1801 int i= 0;
1802 int lf= 0;
1803 int lg= d*k;
1804 int degfSubLf= degf;
1805 int deggSubLg= degg-lg;
1806 int repLengthBuf2, repLengthBuf1, tmp;
1807 while (degf >= lf || lg >= 0)
1808 {
1809 if (degfSubLf >= d)
1810 repLengthBuf1= d;
1811 else if (degfSubLf < 0)
1812 repLengthBuf1= 0;
1813 else
1814 repLengthBuf1= degfSubLf + 1;
1815 fq_nmod_poly_init2 (buf1, repLengthBuf1, fq_con);
1816 _fq_nmod_poly_set_length (buf1, repLengthBuf1, fq_con);
1817
1818 _fq_nmod_vec_set (buf1->coeffs, f->coeffs + lf, repLengthBuf1, fq_con);
1819 _fq_nmod_poly_normalise (buf1, fq_con);
1820
1821 repLengthBuf1= fq_nmod_poly_length (buf1, fq_con);
1822
1823 if (deggSubLg >= d - 1)
1824 repLengthBuf2= d - 1;
1825 else if (deggSubLg < 0)
1826 repLengthBuf2= 0;
1827 else
1828 repLengthBuf2= deggSubLg + 1;
1829
1830 fq_nmod_poly_init2 (buf2, repLengthBuf2, fq_con);
1831 _fq_nmod_poly_set_length (buf2, repLengthBuf2, fq_con);
1832 _fq_nmod_vec_set (buf2->coeffs, g->coeffs + lg, repLengthBuf2, fq_con);
1833
1834 _fq_nmod_poly_normalise (buf2, fq_con);
1835 repLengthBuf2= fq_nmod_poly_length (buf2, fq_con);
1836
1837 fq_nmod_poly_init2 (buf3, repLengthBuf2 + d, fq_con);
1838 _fq_nmod_poly_set_length (buf3, repLengthBuf2 + d, fq_con);
1839 _fq_nmod_vec_set (buf3->coeffs, buf1->coeffs, repLengthBuf1, fq_con);
1840 _fq_nmod_vec_set (buf3->coeffs + d, buf2->coeffs, repLengthBuf2, fq_con);
1841
1842 _fq_nmod_poly_normalise (buf3, fq_con);
1843
1845 i++;
1846
1847
1848 lf= i*d;
1849 degfSubLf= degf - lf;
1850
1851 lg= d*(k - i);
1852 deggSubLg= degg - lg;
1853
1854 if (lg >= 0 && deggSubLg > 0)
1855 {
1856 if (repLengthBuf2 > degfSubLf + 1)
1857 degfSubLf= repLengthBuf2 - 1;
1858 tmp= tmin (repLengthBuf1, deggSubLg + 1);
1859 _fq_nmod_vec_sub (g->coeffs + lg, g->coeffs + lg, buf1-> coeffs,
1860 tmp, fq_con);
1861 }
1862 if (lg < 0)
1863 {
1866 fq_nmod_poly_clear (buf3, fq_con);
1867 break;
1868 }
1869 if (degfSubLf >= 0)
1870 _fq_nmod_vec_sub (f->coeffs + lf, f->coeffs + lf, buf2->coeffs,
1871 repLengthBuf2, fq_con);
1874 fq_nmod_poly_clear (buf3, fq_con);
1875 }
1876
1879
1880 return result;
1881}
1882#endif
1883
1885reverseSubstReciproQ (const fmpz_poly_t F, const fmpz_poly_t G, int d, int k)
1886{
1887 Variable y= Variable (2);
1888 Variable x= Variable (1);
1889
1890 fmpz_poly_t f, g;
1891 fmpz_poly_init (f);
1892 fmpz_poly_init (g);
1893 fmpz_poly_set (f, F);
1894 fmpz_poly_set (g, G);
1895 int degf= fmpz_poly_degree(f);
1896 int degg= fmpz_poly_degree(g);
1897
1898 fmpz_poly_t buf1,buf2, buf3;
1899
1900 if (fmpz_poly_length (f) < (long) d*(k+1)) //zero padding
1901 fmpz_poly_fit_length (f,(long)d*(k+1));
1902
1904 int i= 0;
1905 int lf= 0;
1906 int lg= d*k;
1907 int degfSubLf= degf;
1908 int deggSubLg= degg-lg;
1909 int repLengthBuf2, repLengthBuf1, ind, tmp;
1910 fmpz_t tmp1, tmp2;
1911 while (degf >= lf || lg >= 0)
1912 {
1913 if (degfSubLf >= d)
1914 repLengthBuf1= d;
1915 else if (degfSubLf < 0)
1916 repLengthBuf1= 0;
1917 else
1918 repLengthBuf1= degfSubLf + 1;
1919
1920 fmpz_poly_init2 (buf1, repLengthBuf1);
1921
1922 for (ind= 0; ind < repLengthBuf1; ind++)
1923 {
1924 fmpz_poly_get_coeff_fmpz (tmp1, f, ind + lf);
1925 fmpz_poly_set_coeff_fmpz (buf1, ind, tmp1);
1926 }
1927 _fmpz_poly_normalise (buf1);
1928
1929 repLengthBuf1= fmpz_poly_length (buf1);
1930
1931 if (deggSubLg >= d - 1)
1932 repLengthBuf2= d - 1;
1933 else if (deggSubLg < 0)
1934 repLengthBuf2= 0;
1935 else
1936 repLengthBuf2= deggSubLg + 1;
1937
1938 fmpz_poly_init2 (buf2, repLengthBuf2);
1939
1940 for (ind= 0; ind < repLengthBuf2; ind++)
1941 {
1942 fmpz_poly_get_coeff_fmpz (tmp1, g, ind + lg);
1943 fmpz_poly_set_coeff_fmpz (buf2, ind, tmp1);
1944 }
1945
1946 _fmpz_poly_normalise (buf2);
1947 repLengthBuf2= fmpz_poly_length (buf2);
1948
1949 fmpz_poly_init2 (buf3, repLengthBuf2 + d);
1950 for (ind= 0; ind < repLengthBuf1; ind++)
1951 {
1952 fmpz_poly_get_coeff_fmpz (tmp1, buf1, ind);
1953 fmpz_poly_set_coeff_fmpz (buf3, ind, tmp1);
1954 }
1955 for (ind= repLengthBuf1; ind < d; ind++)
1956 fmpz_poly_set_coeff_ui (buf3, ind, 0);
1957 for (ind= 0; ind < repLengthBuf2; ind++)
1958 {
1959 fmpz_poly_get_coeff_fmpz (tmp1, buf2, ind);
1960 fmpz_poly_set_coeff_fmpz (buf3, ind + d, tmp1);
1961 }
1962 _fmpz_poly_normalise (buf3);
1963
1964 result += convertFmpz_poly_t2FacCF (buf3, x)*power (y, i);
1965 i++;
1966
1967
1968 lf= i*d;
1969 degfSubLf= degf - lf;
1970
1971 lg= d*(k-i);
1972 deggSubLg= degg - lg;
1973
1974 if (lg >= 0 && deggSubLg > 0)
1975 {
1976 if (repLengthBuf2 > degfSubLf + 1)
1977 degfSubLf= repLengthBuf2 - 1;
1978 tmp= tmin (repLengthBuf1, deggSubLg + 1);
1979 for (ind= 0; ind < tmp; ind++)
1980 {
1981 fmpz_poly_get_coeff_fmpz (tmp1, g, ind + lg);
1982 fmpz_poly_get_coeff_fmpz (tmp2, buf1, ind);
1983 fmpz_sub (tmp1, tmp1, tmp2);
1984 fmpz_poly_set_coeff_fmpz (g, ind + lg, tmp1);
1985 }
1986 }
1987 if (lg < 0)
1988 {
1989 fmpz_poly_clear (buf1);
1990 fmpz_poly_clear (buf2);
1991 fmpz_poly_clear (buf3);
1992 break;
1993 }
1994 if (degfSubLf >= 0)
1995 {
1996 for (ind= 0; ind < repLengthBuf2; ind++)
1997 {
1998 fmpz_poly_get_coeff_fmpz (tmp1, f, ind + lf);
1999 fmpz_poly_get_coeff_fmpz (tmp2, buf2, ind);
2000 fmpz_sub (tmp1, tmp1, tmp2);
2001 fmpz_poly_set_coeff_fmpz (f, ind + lf, tmp1);
2002 }
2003 }
2004 fmpz_poly_clear (buf1);
2005 fmpz_poly_clear (buf2);
2006 fmpz_poly_clear (buf3);
2007 }
2008
2009 fmpz_poly_clear (f);
2010 fmpz_poly_clear (g);
2011 fmpz_clear (tmp1);
2012 fmpz_clear (tmp2);
2013
2014 return result;
2015}
2016
2017#if ( __FLINT_RELEASE >= 20400)
2019reverseSubstFq (const fq_nmod_poly_t F, int d, const Variable& alpha,
2020 const fq_nmod_ctx_t fq_con)
2021{
2022 Variable y= Variable (2);
2023 Variable x= Variable (1);
2024
2025 fq_nmod_poly_t buf;
2027 int i= 0;
2028 int degf= fq_nmod_poly_degree(F, fq_con);
2029 int k= 0;
2030 int degfSubK, repLength;
2031 while (degf >= k)
2032 {
2033 degfSubK= degf - k;
2034 if (degfSubK >= d)
2035 repLength= d;
2036 else
2037 repLength= degfSubK + 1;
2038
2039 fq_nmod_poly_init2 (buf, repLength, fq_con);
2040 _fq_nmod_poly_set_length (buf, repLength, fq_con);
2041 _fq_nmod_vec_set (buf->coeffs, F->coeffs+k, repLength, fq_con);
2042 _fq_nmod_poly_normalise (buf, fq_con);
2043
2045 i++;
2046 k= d*i;
2048 }
2049
2050 return result;
2051}
2052#endif
2053
2054CanonicalForm reverseSubstFp (const nmod_poly_t F, int d)
2055{
2056 Variable y= Variable (2);
2057 Variable x= Variable (1);
2058
2059 mp_limb_t ninv= n_preinvert_limb (getCharacteristic());
2060
2061 nmod_poly_t buf;
2063 int i= 0;
2064 int degf= nmod_poly_degree(F);
2065 int k= 0;
2066 int degfSubK, repLength, j;
2067 while (degf >= k)
2068 {
2069 degfSubK= degf - k;
2070 if (degfSubK >= d)
2071 repLength= d;
2072 else
2073 repLength= degfSubK + 1;
2074
2075 nmod_poly_init2_preinv (buf, getCharacteristic(), ninv, repLength);
2076 for (j= 0; j < repLength; j++)
2077 nmod_poly_set_coeff_ui (buf, j, nmod_poly_get_coeff_ui (F, j + k));
2078 _nmod_poly_normalise (buf);
2079
2081 i++;
2082 k= d*i;
2084 }
2085
2086 return result;
2087}
2088
2092{
2093 int d1= degree (F, 1) + degree (G, 1) + 1;
2094 d1 /= 2;
2095 d1 += 1;
2096
2097 nmod_poly_t F1, F2;
2098 kronSubReciproFp (F1, F2, F, d1);
2099
2100 nmod_poly_t G1, G2;
2101 kronSubReciproFp (G1, G2, G, d1);
2102
2103 int k= d1*degree (M);
2104 nmod_poly_mullow (F1, F1, G1, (long) k);
2105
2106 int degtailF= degree (tailcoeff (F), 1);;
2107 int degtailG= degree (tailcoeff (G), 1);
2108 int taildegF= taildegree (F);
2109 int taildegG= taildegree (G);
2110
2111 int b= nmod_poly_degree (F2) + nmod_poly_degree (G2) - k - degtailF - degtailG
2112 + d1*(2+taildegF + taildegG);
2113 nmod_poly_mulhigh (F2, F2, G2, b);
2114 nmod_poly_shift_right (F2, F2, b);
2115 int d2= tmax (nmod_poly_degree (F2)/d1, nmod_poly_degree (F1)/d1);
2116
2117
2119
2122 nmod_poly_clear (G1);
2123 nmod_poly_clear (G2);
2124 return result;
2125}
2126
2130{
2131 CanonicalForm A= F;
2132 CanonicalForm B= G;
2133
2134 int degAx= degree (A, 1);
2135 int degAy= degree (A, 2);
2136 int degBx= degree (B, 1);
2137 int degBy= degree (B, 2);
2138 int d1= degAx + 1 + degBx;
2139 int d2= tmax (degAy, degBy);
2140
2141 if (d1 > 128 && d2 > 160 && (degAy == degBy) && (2*degAy > degree (M)))
2142 return mulMod2FLINTFpReci (A, B, M);
2143
2144 nmod_poly_t FLINTA, FLINTB;
2145 kronSubFp (FLINTA, A, d1);
2146 kronSubFp (FLINTB, B, d1);
2147
2148 int k= d1*degree (M);
2149 nmod_poly_mullow (FLINTA, FLINTA, FLINTB, (long) k);
2150
2151 A= reverseSubstFp (FLINTA, d1);
2152
2153 nmod_poly_clear (FLINTA);
2154 nmod_poly_clear (FLINTB);
2155 return A;
2156}
2157
2158#if ( __FLINT_RELEASE >= 20400)
2161 CanonicalForm& M, const Variable& alpha,
2162 const fq_nmod_ctx_t fq_con)
2163{
2164 int d1= degree (F, 1) + degree (G, 1) + 1;
2165 d1 /= 2;
2166 d1 += 1;
2167
2168 fq_nmod_poly_t F1, F2;
2169 kronSubReciproFq (F1, F2, F, d1, fq_con);
2170
2171 fq_nmod_poly_t G1, G2;
2172 kronSubReciproFq (G1, G2, G, d1, fq_con);
2173
2174 int k= d1*degree (M);
2175 fq_nmod_poly_mullow (F1, F1, G1, (long) k, fq_con);
2176
2177 int degtailF= degree (tailcoeff (F), 1);
2178 int degtailG= degree (tailcoeff (G), 1);
2179 int taildegF= taildegree (F);
2180 int taildegG= taildegree (G);
2181
2182 int b= k + degtailF + degtailG - d1*(2+taildegF + taildegG);
2183
2184 fq_nmod_poly_reverse (F2, F2, fq_nmod_poly_length (F2, fq_con), fq_con);
2185 fq_nmod_poly_reverse (G2, G2, fq_nmod_poly_length (G2, fq_con), fq_con);
2186 fq_nmod_poly_mullow (F2, F2, G2, b+1, fq_con);
2187 fq_nmod_poly_reverse (F2, F2, b+1, fq_con);
2188
2189 int d2= tmax (fq_nmod_poly_degree (F2, fq_con)/d1,
2190 fq_nmod_poly_degree (F1, fq_con)/d1);
2191
2193
2198 return result;
2199}
2200
2203 CanonicalForm& M, const Variable& alpha,
2204 const fq_nmod_ctx_t fq_con)
2205{
2206 CanonicalForm A= F;
2207 CanonicalForm B= G;
2208
2209 int degAx= degree (A, 1);
2210 int degAy= degree (A, 2);
2211 int degBx= degree (B, 1);
2212 int degBy= degree (B, 2);
2213 int d1= degAx + 1 + degBx;
2214 int d2= tmax (degAy, degBy);
2215
2216 if (d1 > 128 && d2 > 160 && (degAy == degBy) && (2*degAy > degree (M)))
2217 return mulMod2FLINTFqReci (A, B, M, alpha, fq_con);
2218
2219 fq_nmod_poly_t FLINTA, FLINTB;
2220 kronSubFq (FLINTA, A, d1, fq_con);
2221 kronSubFq (FLINTB, B, d1, fq_con);
2222
2223 int k= d1*degree (M);
2224 fq_nmod_poly_mullow (FLINTA, FLINTA, FLINTB, (long) k, fq_con);
2225
2226 A= reverseSubstFq (FLINTA, d1, alpha, fq_con);
2227
2228 fq_nmod_poly_clear (FLINTA, fq_con);
2229 fq_nmod_poly_clear (FLINTB, fq_con);
2230 return A;
2231}
2232#endif
2233
2237{
2238 int d1= degree (F, 1) + degree (G, 1) + 1;
2239 d1 /= 2;
2240 d1 += 1;
2241
2242 fmpz_poly_t F1, F2;
2243 kronSubReciproQ (F1, F2, F, d1);
2244
2245 fmpz_poly_t G1, G2;
2246 kronSubReciproQ (G1, G2, G, d1);
2247
2248 int k= d1*degree (M);
2249 fmpz_poly_mullow (F1, F1, G1, (long) k);
2250
2251 int degtailF= degree (tailcoeff (F), 1);;
2252 int degtailG= degree (tailcoeff (G), 1);
2253 int taildegF= taildegree (F);
2254 int taildegG= taildegree (G);
2255
2256 int b= fmpz_poly_degree (F2) + fmpz_poly_degree (G2) - k - degtailF - degtailG
2257 + d1*(2+taildegF + taildegG);
2258 fmpz_poly_mulhigh_n (F2, F2, G2, b);
2259 fmpz_poly_shift_right (F2, F2, b);
2260 int d2= tmax (fmpz_poly_degree (F2)/d1, fmpz_poly_degree (F1)/d1);
2261
2263
2264 fmpz_poly_clear (F1);
2265 fmpz_poly_clear (F2);
2266 fmpz_poly_clear (G1);
2267 fmpz_poly_clear (G2);
2268 return result;
2269}
2270
2274{
2275 CanonicalForm A= F;
2276 CanonicalForm B= G;
2277
2278 int degAx= degree (A, 1);
2279 int degBx= degree (B, 1);
2280 int d1= degAx + 1 + degBx;
2281
2284 A *= f;
2285 B *= g;
2286
2287 fmpz_poly_t FLINTA, FLINTB;
2288 kronSubQa (FLINTA, A, d1);
2289 kronSubQa (FLINTB, B, d1);
2290 int k= d1*degree (M);
2291
2292 fmpz_poly_mullow (FLINTA, FLINTA, FLINTB, (long) k);
2293 A= reverseSubstQ (FLINTA, d1);
2294 fmpz_poly_clear (FLINTA);
2295 fmpz_poly_clear (FLINTB);
2296 return A/(f*g);
2297}
2298
2299/*CanonicalForm
2300mulMod2FLINTQa (const CanonicalForm& F, const CanonicalForm& G,
2301 const CanonicalForm& M)
2302{
2303 Variable a;
2304 if (!hasFirstAlgVar (F,a) && !hasFirstAlgVar (G, a))
2305 return mulMod2FLINTQ (F, G, M);
2306 CanonicalForm A= F;
2307
2308 int degFx= degree (F, 1);
2309 int degFa= degree (F, a);
2310 int degGx= degree (G, 1);
2311 int degGa= degree (G, a);
2312
2313 int d2= degFa+degGa+1;
2314 int d1= degFx + 1 + degGx;
2315 d1 *= d2;
2316
2317 fmpq_poly_t FLINTF, FLINTG;
2318 kronSubQa (FLINTF, F, d1, d2);
2319 kronSubQa (FLINTG, G, d1, d2);
2320
2321 fmpq_poly_mullow (FLINTF, FLINTF, FLINTG, d1*degree (M));
2322
2323 fmpq_poly_t mipo;
2324 convertFacCF2Fmpq_poly_t (mipo, getMipo (a));
2325 CanonicalForm result= reverseSubstQa (FLINTF, d1, d2, a, mipo);
2326 fmpq_poly_clear (FLINTF);
2327 fmpq_poly_clear (FLINTG);
2328 return result;
2329}*/
2330
2333 const CanonicalForm& M)
2334{
2335 Variable a;
2336 if (!hasFirstAlgVar (F,a) && !hasFirstAlgVar (G, a))
2337 return mulMod2FLINTQ (F, G, M);
2338 CanonicalForm A= F, B= G;
2339
2340 int degFx= degree (F, 1);
2341 int degFa= degree (F, a);
2342 int degGx= degree (G, 1);
2343 int degGa= degree (G, a);
2344
2345 int d2= degFa+degGa+1;
2346 int d1= degFx + 1 + degGx;
2347 d1 *= d2;
2348
2351 A *= f;
2352 B *= g;
2353
2354 fmpz_poly_t FLINTF, FLINTG;
2355 kronSubQa (FLINTF, A, d1, d2);
2356 kronSubQa (FLINTG, B, d1, d2);
2357
2358 fmpz_poly_mullow (FLINTF, FLINTF, FLINTG, d1*degree (M));
2359
2360 fmpq_poly_t mipo;
2362 A= reverseSubstQa (FLINTF, d1, d2, a, mipo);
2363 fmpz_poly_clear (FLINTF);
2364 fmpz_poly_clear (FLINTG);
2365 return A/(f*g);
2366}
2367
2368#endif
2369
2370#ifndef HAVE_FLINT
2371zz_pX kronSubFp (const CanonicalForm& A, int d)
2372{
2373 int degAy= degree (A);
2374 zz_pX result;
2375 result.rep.SetLength (d*(degAy + 1));
2376
2377 zz_p *resultp;
2378 resultp= result.rep.elts();
2379 zz_pX buf;
2380 zz_p *bufp;
2381 int j, k, bufRepLength;
2382
2383 for (CFIterator i= A; i.hasTerms(); i++)
2384 {
2385 if (i.coeff().inCoeffDomain())
2386 buf= convertFacCF2NTLzzpX (i.coeff());
2387 else
2388 buf= convertFacCF2NTLzzpX (i.coeff());
2389
2390 k= i.exp()*d;
2391 bufp= buf.rep.elts();
2392 bufRepLength= (int) buf.rep.length();
2393 for (j= 0; j < bufRepLength; j++)
2394 resultp [j + k]= bufp [j];
2395 }
2396 result.normalize();
2397
2398 return result;
2399}
2400#endif
2401
2402#if (!(HAVE_FLINT && __FLINT_RELEASE >= 20400))
2403zz_pEX kronSubFq (const CanonicalForm& A, int d, const Variable& alpha)
2404{
2405 int degAy= degree (A);
2406 zz_pEX result;
2407 result.rep.SetLength (d*(degAy + 1));
2408
2409 Variable v;
2410 zz_pE *resultp;
2411 resultp= result.rep.elts();
2412 zz_pEX buf1;
2413 zz_pE *buf1p;
2414 zz_pX buf2;
2415 zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
2416 int j, k, buf1RepLength;
2417
2418 for (CFIterator i= A; i.hasTerms(); i++)
2419 {
2420 if (i.coeff().inCoeffDomain())
2421 {
2422 buf2= convertFacCF2NTLzzpX (i.coeff());
2423 buf1= to_zz_pEX (to_zz_pE (buf2));
2424 }
2425 else
2426 buf1= convertFacCF2NTLzz_pEX (i.coeff(), NTLMipo);
2427
2428 k= i.exp()*d;
2429 buf1p= buf1.rep.elts();
2430 buf1RepLength= (int) buf1.rep.length();
2431 for (j= 0; j < buf1RepLength; j++)
2432 resultp [j + k]= buf1p [j];
2433 }
2434 result.normalize();
2435
2436 return result;
2437}
2438
2439void
2440kronSubReciproFq (zz_pEX& subA1, zz_pEX& subA2,const CanonicalForm& A, int d,
2441 const Variable& alpha)
2442{
2443 int degAy= degree (A);
2444 subA1.rep.SetLength ((long) d*(degAy + 2));
2445 subA2.rep.SetLength ((long) d*(degAy + 2));
2446
2447 Variable v;
2448 zz_pE *subA1p;
2449 zz_pE *subA2p;
2450 subA1p= subA1.rep.elts();
2451 subA2p= subA2.rep.elts();
2452 zz_pEX buf;
2453 zz_pE *bufp;
2454 zz_pX buf2;
2455 zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
2456 int j, k, kk, bufRepLength;
2457
2458 for (CFIterator i= A; i.hasTerms(); i++)
2459 {
2460 if (i.coeff().inCoeffDomain())
2461 {
2462 buf2= convertFacCF2NTLzzpX (i.coeff());
2463 buf= to_zz_pEX (to_zz_pE (buf2));
2464 }
2465 else
2466 buf= convertFacCF2NTLzz_pEX (i.coeff(), NTLMipo);
2467
2468 k= i.exp()*d;
2469 kk= (degAy - i.exp())*d;
2470 bufp= buf.rep.elts();
2471 bufRepLength= (int) buf.rep.length();
2472 for (j= 0; j < bufRepLength; j++)
2473 {
2474 subA1p [j + k] += bufp [j];
2475 subA2p [j + kk] += bufp [j];
2476 }
2477 }
2478 subA1.normalize();
2479 subA2.normalize();
2480}
2481#endif
2482
2483#ifndef HAVE_FLINT
2484void
2485kronSubReciproFp (zz_pX& subA1, zz_pX& subA2, const CanonicalForm& A, int d)
2486{
2487 int degAy= degree (A);
2488 subA1.rep.SetLength ((long) d*(degAy + 2));
2489 subA2.rep.SetLength ((long) d*(degAy + 2));
2490
2491 zz_p *subA1p;
2492 zz_p *subA2p;
2493 subA1p= subA1.rep.elts();
2494 subA2p= subA2.rep.elts();
2495 zz_pX buf;
2496 zz_p *bufp;
2497 int j, k, kk, bufRepLength;
2498
2499 for (CFIterator i= A; i.hasTerms(); i++)
2500 {
2501 buf= convertFacCF2NTLzzpX (i.coeff());
2502
2503 k= i.exp()*d;
2504 kk= (degAy - i.exp())*d;
2505 bufp= buf.rep.elts();
2506 bufRepLength= (int) buf.rep.length();
2507 for (j= 0; j < bufRepLength; j++)
2508 {
2509 subA1p [j + k] += bufp [j];
2510 subA2p [j + kk] += bufp [j];
2511 }
2512 }
2513 subA1.normalize();
2514 subA2.normalize();
2515}
2516#endif
2517
2518#if (!(HAVE_FLINT && __FLINT_RELEASE >= 20400))
2520reverseSubstReciproFq (const zz_pEX& F, const zz_pEX& G, int d, int k,
2521 const Variable& alpha)
2522{
2523 Variable y= Variable (2);
2524 Variable x= Variable (1);
2525
2526 zz_pEX f= F;
2527 zz_pEX g= G;
2528 int degf= deg(f);
2529 int degg= deg(g);
2530
2531 zz_pEX buf1;
2532 zz_pEX buf2;
2533 zz_pEX buf3;
2534
2535 zz_pE *buf1p;
2536 zz_pE *buf2p;
2537 zz_pE *buf3p;
2538 if (f.rep.length() < (long) d*(k+1)) //zero padding
2539 f.rep.SetLength ((long)d*(k+1));
2540
2541 zz_pE *gp= g.rep.elts();
2542 zz_pE *fp= f.rep.elts();
2544 int i= 0;
2545 int lf= 0;
2546 int lg= d*k;
2547 int degfSubLf= degf;
2548 int deggSubLg= degg-lg;
2549 int repLengthBuf2, repLengthBuf1, ind, tmp;
2550 zz_pE zzpEZero= zz_pE();
2551
2552 while (degf >= lf || lg >= 0)
2553 {
2554 if (degfSubLf >= d)
2555 repLengthBuf1= d;
2556 else if (degfSubLf < 0)
2557 repLengthBuf1= 0;
2558 else
2559 repLengthBuf1= degfSubLf + 1;
2560 buf1.rep.SetLength((long) repLengthBuf1);
2561
2562 buf1p= buf1.rep.elts();
2563 for (ind= 0; ind < repLengthBuf1; ind++)
2564 buf1p [ind]= fp [ind + lf];
2565 buf1.normalize();
2566
2567 repLengthBuf1= buf1.rep.length();
2568
2569 if (deggSubLg >= d - 1)
2570 repLengthBuf2= d - 1;
2571 else if (deggSubLg < 0)
2572 repLengthBuf2= 0;
2573 else
2574 repLengthBuf2= deggSubLg + 1;
2575
2576 buf2.rep.SetLength ((long) repLengthBuf2);
2577 buf2p= buf2.rep.elts();
2578 for (ind= 0; ind < repLengthBuf2; ind++)
2579 buf2p [ind]= gp [ind + lg];
2580 buf2.normalize();
2581
2582 repLengthBuf2= buf2.rep.length();
2583
2584 buf3.rep.SetLength((long) repLengthBuf2 + d);
2585 buf3p= buf3.rep.elts();
2586 buf2p= buf2.rep.elts();
2587 buf1p= buf1.rep.elts();
2588 for (ind= 0; ind < repLengthBuf1; ind++)
2589 buf3p [ind]= buf1p [ind];
2590 for (ind= repLengthBuf1; ind < d; ind++)
2591 buf3p [ind]= zzpEZero;
2592 for (ind= 0; ind < repLengthBuf2; ind++)
2593 buf3p [ind + d]= buf2p [ind];
2594 buf3.normalize();
2595
2596 result += convertNTLzz_pEX2CF (buf3, x, alpha)*power (y, i);
2597 i++;
2598
2599
2600 lf= i*d;
2601 degfSubLf= degf - lf;
2602
2603 lg= d*(k-i);
2604 deggSubLg= degg - lg;
2605
2606 buf1p= buf1.rep.elts();
2607
2608 if (lg >= 0 && deggSubLg > 0)
2609 {
2610 if (repLengthBuf2 > degfSubLf + 1)
2611 degfSubLf= repLengthBuf2 - 1;
2612 tmp= tmin (repLengthBuf1, deggSubLg + 1);
2613 for (ind= 0; ind < tmp; ind++)
2614 gp [ind + lg] -= buf1p [ind];
2615 }
2616
2617 if (lg < 0)
2618 break;
2619
2620 buf2p= buf2.rep.elts();
2621 if (degfSubLf >= 0)
2622 {
2623 for (ind= 0; ind < repLengthBuf2; ind++)
2624 fp [ind + lf] -= buf2p [ind];
2625 }
2626 }
2627
2628 return result;
2629}
2630#endif
2631
2632#ifndef HAVE_FLINT
2634reverseSubstReciproFp (const zz_pX& F, const zz_pX& G, int d, int k)
2635{
2636 Variable y= Variable (2);
2637 Variable x= Variable (1);
2638
2639 zz_pX f= F;
2640 zz_pX g= G;
2641 int degf= deg(f);
2642 int degg= deg(g);
2643
2644 zz_pX buf1;
2645 zz_pX buf2;
2646 zz_pX buf3;
2647
2648 zz_p *buf1p;
2649 zz_p *buf2p;
2650 zz_p *buf3p;
2651
2652 if (f.rep.length() < (long) d*(k+1)) //zero padding
2653 f.rep.SetLength ((long)d*(k+1));
2654
2655 zz_p *gp= g.rep.elts();
2656 zz_p *fp= f.rep.elts();
2658 int i= 0;
2659 int lf= 0;
2660 int lg= d*k;
2661 int degfSubLf= degf;
2662 int deggSubLg= degg-lg;
2663 int repLengthBuf2, repLengthBuf1, ind, tmp;
2664 zz_p zzpZero= zz_p();
2665 while (degf >= lf || lg >= 0)
2666 {
2667 if (degfSubLf >= d)
2668 repLengthBuf1= d;
2669 else if (degfSubLf < 0)
2670 repLengthBuf1= 0;
2671 else
2672 repLengthBuf1= degfSubLf + 1;
2673 buf1.rep.SetLength((long) repLengthBuf1);
2674
2675 buf1p= buf1.rep.elts();
2676 for (ind= 0; ind < repLengthBuf1; ind++)
2677 buf1p [ind]= fp [ind + lf];
2678 buf1.normalize();
2679
2680 repLengthBuf1= buf1.rep.length();
2681
2682 if (deggSubLg >= d - 1)
2683 repLengthBuf2= d - 1;
2684 else if (deggSubLg < 0)
2685 repLengthBuf2= 0;
2686 else
2687 repLengthBuf2= deggSubLg + 1;
2688
2689 buf2.rep.SetLength ((long) repLengthBuf2);
2690 buf2p= buf2.rep.elts();
2691 for (ind= 0; ind < repLengthBuf2; ind++)
2692 buf2p [ind]= gp [ind + lg];
2693
2694 buf2.normalize();
2695
2696 repLengthBuf2= buf2.rep.length();
2697
2698
2699 buf3.rep.SetLength((long) repLengthBuf2 + d);
2700 buf3p= buf3.rep.elts();
2701 buf2p= buf2.rep.elts();
2702 buf1p= buf1.rep.elts();
2703 for (ind= 0; ind < repLengthBuf1; ind++)
2704 buf3p [ind]= buf1p [ind];
2705 for (ind= repLengthBuf1; ind < d; ind++)
2706 buf3p [ind]= zzpZero;
2707 for (ind= 0; ind < repLengthBuf2; ind++)
2708 buf3p [ind + d]= buf2p [ind];
2709 buf3.normalize();
2710
2711 result += convertNTLzzpX2CF (buf3, x)*power (y, i);
2712 i++;
2713
2714
2715 lf= i*d;
2716 degfSubLf= degf - lf;
2717
2718 lg= d*(k-i);
2719 deggSubLg= degg - lg;
2720
2721 buf1p= buf1.rep.elts();
2722
2723 if (lg >= 0 && deggSubLg > 0)
2724 {
2725 if (repLengthBuf2 > degfSubLf + 1)
2726 degfSubLf= repLengthBuf2 - 1;
2727 tmp= tmin (repLengthBuf1, deggSubLg + 1);
2728 for (ind= 0; ind < tmp; ind++)
2729 gp [ind + lg] -= buf1p [ind];
2730 }
2731 if (lg < 0)
2732 break;
2733
2734 buf2p= buf2.rep.elts();
2735 if (degfSubLf >= 0)
2736 {
2737 for (ind= 0; ind < repLengthBuf2; ind++)
2738 fp [ind + lf] -= buf2p [ind];
2739 }
2740 }
2741
2742 return result;
2743}
2744#endif
2745
2746#if (!(HAVE_FLINT && __FLINT_RELEASE >= 20400))
2747CanonicalForm reverseSubstFq (const zz_pEX& F, int d, const Variable& alpha)
2748{
2749 Variable y= Variable (2);
2750 Variable x= Variable (1);
2751
2752 zz_pEX f= F;
2753 zz_pE *fp= f.rep.elts();
2754
2755 zz_pEX buf;
2756 zz_pE *bufp;
2758 int i= 0;
2759 int degf= deg(f);
2760 int k= 0;
2761 int degfSubK, repLength, j;
2762 while (degf >= k)
2763 {
2764 degfSubK= degf - k;
2765 if (degfSubK >= d)
2766 repLength= d;
2767 else
2768 repLength= degfSubK + 1;
2769
2770 buf.rep.SetLength ((long) repLength);
2771 bufp= buf.rep.elts();
2772 for (j= 0; j < repLength; j++)
2773 bufp [j]= fp [j + k];
2774 buf.normalize();
2775
2777 i++;
2778 k= d*i;
2779 }
2780
2781 return result;
2782}
2783#endif
2784
2785#ifndef HAVE_FLINT
2786CanonicalForm reverseSubstFp (const zz_pX& F, int d)
2787{
2788 Variable y= Variable (2);
2789 Variable x= Variable (1);
2790
2791 zz_pX f= F;
2792 zz_p *fp= f.rep.elts();
2793
2794 zz_pX buf;
2795 zz_p *bufp;
2797 int i= 0;
2798 int degf= deg(f);
2799 int k= 0;
2800 int degfSubK, repLength, j;
2801 while (degf >= k)
2802 {
2803 degfSubK= degf - k;
2804 if (degfSubK >= d)
2805 repLength= d;
2806 else
2807 repLength= degfSubK + 1;
2808
2809 buf.rep.SetLength ((long) repLength);
2810 bufp= buf.rep.elts();
2811 for (j= 0; j < repLength; j++)
2812 bufp [j]= fp [j + k];
2813 buf.normalize();
2814
2816 i++;
2817 k= d*i;
2818 }
2819
2820 return result;
2821}
2822
2823// assumes input to be reduced mod M and to be an element of Fp
2825mulMod2NTLFpReci (const CanonicalForm& F, const CanonicalForm& G, const
2827{
2828 int d1= degree (F, 1) + degree (G, 1) + 1;
2829 d1 /= 2;
2830 d1 += 1;
2831
2832 zz_pX F1, F2;
2833 kronSubReciproFp (F1, F2, F, d1);
2834 zz_pX G1, G2;
2835 kronSubReciproFp (G1, G2, G, d1);
2836
2837 int k= d1*degree (M);
2838 MulTrunc (F1, F1, G1, (long) k);
2839
2840 int degtailF= degree (tailcoeff (F), 1);
2841 int degtailG= degree (tailcoeff (G), 1);
2842 int taildegF= taildegree (F);
2843 int taildegG= taildegree (G);
2844 int b= k + degtailF + degtailG - d1*(2+taildegF+taildegG);
2845
2846 reverse (F2, F2);
2847 reverse (G2, G2);
2848 MulTrunc (F2, F2, G2, b + 1);
2849 reverse (F2, F2, b);
2850
2851 int d2= tmax (deg (F2)/d1, deg (F1)/d1);
2852 return reverseSubstReciproFp (F1, F2, d1, d2);
2853}
2854
2855//Kronecker substitution
2857mulMod2NTLFp (const CanonicalForm& F, const CanonicalForm& G, const
2859{
2860 CanonicalForm A= F;
2861 CanonicalForm B= G;
2862
2863 int degAx= degree (A, 1);
2864 int degAy= degree (A, 2);
2865 int degBx= degree (B, 1);
2866 int degBy= degree (B, 2);
2867 int d1= degAx + 1 + degBx;
2868 int d2= tmax (degAy, degBy);
2869
2870 if (d1 > 128 && d2 > 160 && (degAy == degBy) && (2*degAy > degree (M)))
2871 return mulMod2NTLFpReci (A, B, M);
2872
2873 zz_pX NTLA= kronSubFp (A, d1);
2874 zz_pX NTLB= kronSubFp (B, d1);
2875
2876 int k= d1*degree (M);
2877 MulTrunc (NTLA, NTLA, NTLB, (long) k);
2878
2879 A= reverseSubstFp (NTLA, d1);
2880
2881 return A;
2882}
2883#endif
2884
2885#if (!(HAVE_FLINT && __FLINT_RELEASE >= 20400))
2886// assumes input to be reduced mod M and to be an element of Fq not Fp
2888mulMod2NTLFqReci (const CanonicalForm& F, const CanonicalForm& G, const
2889 CanonicalForm& M, const Variable& alpha)
2890{
2891 int d1= degree (F, 1) + degree (G, 1) + 1;
2892 d1 /= 2;
2893 d1 += 1;
2894
2895 zz_pEX F1, F2;
2896 kronSubReciproFq (F1, F2, F, d1, alpha);
2897 zz_pEX G1, G2;
2898 kronSubReciproFq (G1, G2, G, d1, alpha);
2899
2900 int k= d1*degree (M);
2901 MulTrunc (F1, F1, G1, (long) k);
2902
2903 int degtailF= degree (tailcoeff (F), 1);
2904 int degtailG= degree (tailcoeff (G), 1);
2905 int taildegF= taildegree (F);
2906 int taildegG= taildegree (G);
2907 int b= k + degtailF + degtailG - d1*(2+taildegF+taildegG);
2908
2909 reverse (F2, F2);
2910 reverse (G2, G2);
2911 MulTrunc (F2, F2, G2, b + 1);
2912 reverse (F2, F2, b);
2913
2914 int d2= tmax (deg (F2)/d1, deg (F1)/d1);
2915 return reverseSubstReciproFq (F1, F2, d1, d2, alpha);
2916}
2917#endif
2918
2919#ifdef HAVE_FLINT
2921mulMod2FLINTFp (const CanonicalForm& F, const CanonicalForm& G, const
2922 CanonicalForm& M);
2923#endif
2924
2928{
2930 CanonicalForm A= F;
2931 CanonicalForm B= G;
2932
2934 {
2935#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
2936 nmod_poly_t FLINTmipo;
2937 convertFacCF2nmod_poly_t (FLINTmipo, getMipo (alpha));
2938
2939 fq_nmod_ctx_t fq_con;
2940 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
2941
2942 A= mulMod2FLINTFq (A, B, M, alpha, fq_con);
2943 nmod_poly_clear (FLINTmipo);
2945#else
2946 int degAx= degree (A, 1);
2947 int degAy= degree (A, 2);
2948 int degBx= degree (B, 1);
2949 int degBy= degree (B, 2);
2950 int d1= degAx + degBx + 1;
2951 int d2= tmax (degAy, degBy);
2953 {
2955 zz_p::init (getCharacteristic());
2956 }
2957 zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
2958 zz_pE::init (NTLMipo);
2959
2960 int degMipo= degree (getMipo (alpha));
2961 if ((d1 > 128/degMipo) && (d2 > 160/degMipo) && (degAy == degBy) &&
2962 (2*degAy > degree (M)))
2963 return mulMod2NTLFqReci (A, B, M, alpha);
2964
2965 zz_pEX NTLA= kronSubFq (A, d1, alpha);
2966 zz_pEX NTLB= kronSubFq (B, d1, alpha);
2967
2968 int k= d1*degree (M);
2969
2970 MulTrunc (NTLA, NTLA, NTLB, (long) k);
2971
2972 A= reverseSubstFq (NTLA, d1, alpha);
2973#endif
2974 }
2975 else
2976 {
2977#ifdef HAVE_FLINT
2978 A= mulMod2FLINTFp (A, B, M);
2979#else
2980 A= mulMod2NTLFp (A, B, M);
2981#endif
2982 }
2983 return A;
2984}
2985
2987 const CanonicalForm& M)
2988{
2989 if (A.isZero() || B.isZero())
2990 return 0;
2991
2992 ASSERT (M.isUnivariate(), "M must be univariate");
2993
2994 CanonicalForm F= mod (A, M);
2995 CanonicalForm G= mod (B, M);
2996 if (F.inCoeffDomain())
2997 return G*F;
2998 if (G.inCoeffDomain())
2999 return F*G;
3000
3001 Variable y= M.mvar();
3002 int degF= degree (F, y);
3003 int degG= degree (G, y);
3004
3005 if ((degF < 1 && degG < 1) && (F.isUnivariate() && G.isUnivariate()) &&
3006 (F.level() == G.level()))
3007 {
3009 return mod (result, M);
3010 }
3011 else if (degF <= 1 && degG <= 1)
3012 {
3014 return mod (result, M);
3015 }
3016
3017 int sizeF= size (F);
3018 int sizeG= size (G);
3019
3020 int fallBackToNaive= 50;
3021 if (sizeF < fallBackToNaive || sizeG < fallBackToNaive)
3022 {
3023 if (sizeF < sizeG)
3024 return mod (G*F, M);
3025 else
3026 return mod (F*G, M);
3027 }
3028
3029#ifdef HAVE_FLINT
3030 if (getCharacteristic() == 0)
3031 return mulMod2FLINTQa (F, G, M);
3032#endif
3033
3035 (((degF-degG) < 50 && degF > degG) || ((degG-degF) < 50 && degF <= degG)))
3036 return mulMod2NTLFq (F, G, M);
3037
3038 int m= (int) ceil (degree (M)/2.0);
3039 if (degF >= m || degG >= m)
3040 {
3041 CanonicalForm MLo= power (y, m);
3042 CanonicalForm MHi= power (y, degree (M) - m);
3043 CanonicalForm F0= mod (F, MLo);
3044 CanonicalForm F1= div (F, MLo);
3045 CanonicalForm G0= mod (G, MLo);
3046 CanonicalForm G1= div (G, MLo);
3047 CanonicalForm F0G1= mulMod2 (F0, G1, MHi);
3048 CanonicalForm F1G0= mulMod2 (F1, G0, MHi);
3049 CanonicalForm F0G0= mulMod2 (F0, G0, M);
3050 return F0G0 + MLo*(F0G1 + F1G0);
3051 }
3052 else
3053 {
3054 m= (int) ceil (tmax (degF, degG)/2.0);
3055 CanonicalForm yToM= power (y, m);
3056 CanonicalForm F0= mod (F, yToM);
3057 CanonicalForm F1= div (F, yToM);
3058 CanonicalForm G0= mod (G, yToM);
3059 CanonicalForm G1= div (G, yToM);
3060 CanonicalForm H00= mulMod2 (F0, G0, M);
3061 CanonicalForm H11= mulMod2 (F1, G1, M);
3062 CanonicalForm H01= mulMod2 (F0 + F1, G0 + G1, M);
3063 return H11*yToM*yToM + (H01 - H11 - H00)*yToM + H00;
3064 }
3065 DEBOUTLN (cerr, "fatal end in mulMod2");
3066}
3067
3068// end bivariate polys
3069//**********************
3070// multivariate polys
3071
3073{
3074 CanonicalForm A= F;
3075 for (CFListIterator i= M; i.hasItem(); i++)
3076 A= mod (A, i.getItem());
3077 return A;
3078}
3079
3081 const CFList& MOD)
3082{
3083 if (A.isZero() || B.isZero())
3084 return 0;
3085
3086 if (MOD.length() == 1)
3087 return mulMod2 (A, B, MOD.getLast());
3088
3089 CanonicalForm M= MOD.getLast();
3090 CanonicalForm F= mod (A, M);
3091 CanonicalForm G= mod (B, M);
3092 if (F.inCoeffDomain())
3093 return G*F;
3094 if (G.inCoeffDomain())
3095 return F*G;
3096
3097 int sizeF= size (F);
3098 int sizeG= size (G);
3099
3100 if (sizeF / MOD.length() < 100 || sizeG / MOD.length() < 100)
3101 {
3102 if (sizeF < sizeG)
3103 return mod (G*F, MOD);
3104 else
3105 return mod (F*G, MOD);
3106 }
3107
3108 Variable y= M.mvar();
3109 int degF= degree (F, y);
3110 int degG= degree (G, y);
3111
3112 if ((degF <= 1 && F.level() <= M.level()) &&
3113 (degG <= 1 && G.level() <= M.level()))
3114 {
3115 CFList buf= MOD;
3116 buf.removeLast();
3117 if (degF == 1 && degG == 1)
3118 {
3119 CanonicalForm F0= mod (F, y);
3120 CanonicalForm F1= div (F, y);
3121 CanonicalForm G0= mod (G, y);
3122 CanonicalForm G1= div (G, y);
3123 if (degree (M) > 2)
3124 {
3125 CanonicalForm H00= mulMod (F0, G0, buf);
3126 CanonicalForm H11= mulMod (F1, G1, buf);
3127 CanonicalForm H01= mulMod (F0 + F1, G0 + G1, buf);
3128 return H11*y*y + (H01 - H00 - H11)*y + H00;
3129 }
3130 else //here degree (M) == 2
3131 {
3132 buf.append (y);
3133 CanonicalForm F0G1= mulMod (F0, G1, buf);
3134 CanonicalForm F1G0= mulMod (F1, G0, buf);
3135 CanonicalForm F0G0= mulMod (F0, G0, MOD);
3136 CanonicalForm result= F0G0 + y*(F0G1 + F1G0);
3137 return result;
3138 }
3139 }
3140 else if (degF == 1 && degG == 0)
3141 return mulMod (div (F, y), G, buf)*y + mulMod (mod (F, y), G, buf);
3142 else if (degF == 0 && degG == 1)
3143 return mulMod (div (G, y), F, buf)*y + mulMod (mod (G, y), F, buf);
3144 else
3145 return mulMod (F, G, buf);
3146 }
3147 int m= (int) ceil (degree (M)/2.0);
3148 if (degF >= m || degG >= m)
3149 {
3150 CanonicalForm MLo= power (y, m);
3151 CanonicalForm MHi= power (y, degree (M) - m);
3152 CanonicalForm F0= mod (F, MLo);
3153 CanonicalForm F1= div (F, MLo);
3154 CanonicalForm G0= mod (G, MLo);
3155 CanonicalForm G1= div (G, MLo);
3156 CFList buf= MOD;
3157 buf.removeLast();
3158 buf.append (MHi);
3159 CanonicalForm F0G1= mulMod (F0, G1, buf);
3160 CanonicalForm F1G0= mulMod (F1, G0, buf);
3161 CanonicalForm F0G0= mulMod (F0, G0, MOD);
3162 return F0G0 + MLo*(F0G1 + F1G0);
3163 }
3164 else
3165 {
3166 m= (tmax(degF, degG)+1)/2;
3167 CanonicalForm yToM= power (y, m);
3168 CanonicalForm F0= mod (F, yToM);
3169 CanonicalForm F1= div (F, yToM);
3170 CanonicalForm G0= mod (G, yToM);
3171 CanonicalForm G1= div (G, yToM);
3172 CanonicalForm H00= mulMod (F0, G0, MOD);
3173 CanonicalForm H11= mulMod (F1, G1, MOD);
3174 CanonicalForm H01= mulMod (F0 + F1, G0 + G1, MOD);
3175 return H11*yToM*yToM + (H01 - H11 - H00)*yToM + H00;
3176 }
3177 DEBOUTLN (cerr, "fatal end in mulMod");
3178}
3179
3181{
3182 if (L.isEmpty())
3183 return 1;
3184 int l= L.length();
3185 if (l == 1)
3186 return mod (L.getFirst(), M);
3187 else if (l == 2) {
3189 return result;
3190 }
3191 else
3192 {
3193 l /= 2;
3194 CFList tmp1, tmp2;
3195 CFListIterator i= L;
3197 for (int j= 1; j <= l; j++, i++)
3198 tmp1.append (i.getItem());
3199 tmp2= Difference (L, tmp1);
3200 buf1= prodMod (tmp1, M);
3201 buf2= prodMod (tmp2, M);
3203 return result;
3204 }
3205}
3206
3208{
3209 if (L.isEmpty())
3210 return 1;
3211 else if (L.length() == 1)
3212 return L.getFirst();
3213 else if (L.length() == 2)
3214 return mulMod (L.getFirst(), L.getLast(), M);
3215 else
3216 {
3217 int l= L.length()/2;
3218 CFListIterator i= L;
3219 CFList tmp1, tmp2;
3221 for (int j= 1; j <= l; j++, i++)
3222 tmp1.append (i.getItem());
3223 tmp2= Difference (L, tmp1);
3224 buf1= prodMod (tmp1, M);
3225 buf2= prodMod (tmp2, M);
3226 return mulMod (buf1, buf2, M);
3227 }
3228}
3229
3230// end multivariate polys
3231//***************************
3232// division
3233
3235{
3236 if (d == 0)
3237 return F;
3238 CanonicalForm A= F;
3239 Variable y= Variable (2);
3240 Variable x= Variable (1);
3241 if (degree (A, x) > 0)
3242 {
3243 A= swapvar (A, x, y);
3245 CFIterator i= A;
3246 while (d - i.exp() < 0)
3247 i++;
3248
3249 for (; i.hasTerms() && (d - i.exp() >= 0); i++)
3250 result += swapvar (i.coeff(),x,y)*power (x, d - i.exp());
3251 return result;
3252 }
3253 else
3254 return A*power (x, d);
3255}
3256
3258newtonInverse (const CanonicalForm& F, const int n, const CanonicalForm& M)
3259{
3260 int l= ilog2(n);
3261
3262 CanonicalForm g= mod (F, M)[0] [0];
3263
3264 ASSERT (!g.isZero(), "expected a unit");
3265
3267
3268 if (!g.isOne())
3269 g = 1/g;
3270 Variable x= Variable (1);
3272 int exp= 0;
3273 if (n & 1)
3274 {
3275 result= g;
3276 exp= 1;
3277 }
3279
3280 for (int i= 1; i <= l; i++)
3281 {
3282 h= mulMod2 (g, mod (F, power (x, (1 << i))), M);
3283 h= mod (h, power (x, (1 << i)) - 1);
3284 h= div (h, power (x, (1 << (i - 1))));
3285 h= mod (h, M);
3286 g -= power (x, (1 << (i - 1)))*
3287 mod (mulMod2 (g, h, M), power (x, (1 << (i - 1))));
3288
3289 if (n & (1 << i))
3290 {
3291 if (exp)
3292 {
3293 h= mulMod2 (result, mod (F, power (x, exp + (1 << i))), M);
3294 h= mod (h, power (x, exp + (1 << i)) - 1);
3295 h= div (h, power (x, exp));
3296 h= mod (h, M);
3297 result -= power(x, exp)*mod (mulMod2 (g, h, M),
3298 power (x, (1 << i)));
3299 exp += (1 << i);
3300 }
3301 else
3302 {
3303 exp= (1 << i);
3304 result= g;
3305 }
3306 }
3307 }
3308
3309 return result;
3310}
3311
3314 M)
3315{
3316 ASSERT (getCharacteristic() > 0, "positive characteristic expected");
3317
3318 CanonicalForm A= mod (F, M);
3319 CanonicalForm B= mod (G, M);
3320
3321 Variable x= Variable (1);
3322 int degA= degree (A, x);
3323 int degB= degree (B, x);
3324 int m= degA - degB;
3325 if (m < 0)
3326 return 0;
3327
3328 Variable v;
3330 if (degB < 1 || CFFactory::gettype() == GaloisFieldDomain)
3331 {
3333 divrem2 (A, B, Q, R, M);
3334 }
3335 else
3336 {
3337 if (hasFirstAlgVar (A, v) || hasFirstAlgVar (B, v))
3338 {
3339 CanonicalForm R= reverse (A, degA);
3340 CanonicalForm revB= reverse (B, degB);
3341 revB= newtonInverse (revB, m + 1, M);
3342 Q= mulMod2 (R, revB, M);
3343 Q= mod (Q, power (x, m + 1));
3344 Q= reverse (Q, m);
3345 }
3346 else
3347 {
3348 Variable y= Variable (2);
3349#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
3350 nmod_poly_t FLINTmipo;
3351 fq_nmod_ctx_t fq_con;
3352
3353 nmod_poly_init (FLINTmipo, getCharacteristic());
3354 convertFacCF2nmod_poly_t (FLINTmipo, M);
3355
3356 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
3357
3358
3359 fq_nmod_poly_t FLINTA, FLINTB;
3362
3363 fq_nmod_poly_divrem (FLINTA, FLINTB, FLINTA, FLINTB, fq_con);
3364
3365 Q= convertFq_nmod_poly_t2FacCF (FLINTA, x, y, fq_con);
3366
3367 fq_nmod_poly_clear (FLINTA, fq_con);
3368 fq_nmod_poly_clear (FLINTB, fq_con);
3369 nmod_poly_clear (FLINTmipo);
3371#else
3372 bool zz_pEbak= zz_pE::initialized();
3373 zz_pEBak bak;
3374 if (zz_pEbak)
3375 bak.save();
3376 zz_pX mipo= convertFacCF2NTLzzpX (M);
3377 zz_pEX NTLA, NTLB;
3378 NTLA= convertFacCF2NTLzz_pEX (swapvar (A, x, y), mipo);
3379 NTLB= convertFacCF2NTLzz_pEX (swapvar (B, x, y), mipo);
3380 div (NTLA, NTLA, NTLB);
3381 Q= convertNTLzz_pEX2CF (NTLA, x, y);
3382 if (zz_pEbak)
3383 bak.restore();
3384#endif
3385 }
3386 }
3387
3388 return Q;
3389}
3390
3391void
3393 CanonicalForm& R, const CanonicalForm& M)
3394{
3395 CanonicalForm A= mod (F, M);
3396 CanonicalForm B= mod (G, M);
3397 Variable x= Variable (1);
3398 int degA= degree (A, x);
3399 int degB= degree (B, x);
3400 int m= degA - degB;
3401
3402 if (m < 0)
3403 {
3404 R= A;
3405 Q= 0;
3406 return;
3407 }
3408
3409 Variable v;
3410 if (degB <= 1 || CFFactory::gettype() == GaloisFieldDomain)
3411 {
3412 divrem2 (A, B, Q, R, M);
3413 }
3414 else
3415 {
3416 if (hasFirstAlgVar (A, v) || hasFirstAlgVar (B, v))
3417 {
3418 R= reverse (A, degA);
3419
3420 CanonicalForm revB= reverse (B, degB);
3421 revB= newtonInverse (revB, m + 1, M);
3422 Q= mulMod2 (R, revB, M);
3423
3424 Q= mod (Q, power (x, m + 1));
3425 Q= reverse (Q, m);
3426
3427 R= A - mulMod2 (Q, B, M);
3428 }
3429 else
3430 {
3431 Variable y= Variable (2);
3432#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
3433 nmod_poly_t FLINTmipo;
3434 fq_nmod_ctx_t fq_con;
3435
3436 nmod_poly_init (FLINTmipo, getCharacteristic());
3437 convertFacCF2nmod_poly_t (FLINTmipo, M);
3438
3439 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
3440
3441 fq_nmod_poly_t FLINTA, FLINTB;
3444
3445 fq_nmod_poly_divrem (FLINTA, FLINTB, FLINTA, FLINTB, fq_con);
3446
3447 Q= convertFq_nmod_poly_t2FacCF (FLINTA, x, y, fq_con);
3448 R= convertFq_nmod_poly_t2FacCF (FLINTB, x, y, fq_con);
3449
3450 fq_nmod_poly_clear (FLINTA, fq_con);
3451 fq_nmod_poly_clear (FLINTB, fq_con);
3452 nmod_poly_clear (FLINTmipo);
3454#else
3455 zz_pX mipo= convertFacCF2NTLzzpX (M);
3456 zz_pEX NTLA, NTLB;
3457 NTLA= convertFacCF2NTLzz_pEX (swapvar (A, x, y), mipo);
3458 NTLB= convertFacCF2NTLzz_pEX (swapvar (B, x, y), mipo);
3459 zz_pEX NTLQ, NTLR;
3460 DivRem (NTLQ, NTLR, NTLA, NTLB);
3461 Q= convertNTLzz_pEX2CF (NTLQ, x, y);
3462 R= convertNTLzz_pEX2CF (NTLR, x, y);
3463#endif
3464 }
3465 }
3466}
3467
3468static inline
3469CFList split (const CanonicalForm& F, const int m, const Variable& x)
3470{
3471 CanonicalForm A= F;
3472 CanonicalForm buf= 0;
3473 bool swap= false;
3474 if (degree (A, x) <= 0)
3475 return CFList(A);
3476 else if (x.level() != A.level())
3477 {
3478 swap= true;
3479 A= swapvar (A, x, A.mvar());
3480 }
3481
3482 int j= (int) floor ((double) degree (A)/ m);
3483 CFList result;
3484 CFIterator i= A;
3485 for (; j >= 0; j--)
3486 {
3487 while (i.hasTerms() && i.exp() - j*m >= 0)
3488 {
3489 if (swap)
3490 buf += i.coeff()*power (A.mvar(), i.exp() - j*m);
3491 else
3492 buf += i.coeff()*power (x, i.exp() - j*m);
3493 i++;
3494 }
3495 if (swap)
3496 result.append (swapvar (buf, x, F.mvar()));
3497 else
3498 result.append (buf);
3499 buf= 0;
3500 }
3501 return result;
3502}
3503
3504static inline
3505void divrem32 (const CanonicalForm& F, const CanonicalForm& G, CanonicalForm& Q,
3506 CanonicalForm& R, const CFList& M);
3507
3508static inline
3510 CanonicalForm& R, const CFList& M)
3511{
3512 CanonicalForm A= mod (F, M);
3513 CanonicalForm B= mod (G, M);
3514 Variable x= Variable (1);
3515 int degB= degree (B, x);
3516 int degA= degree (A, x);
3517 if (degA < degB)
3518 {
3519 Q= 0;
3520 R= A;
3521 return;
3522 }
3523 if (degB < 1)
3524 {
3525 divrem (A, B, Q, R);
3526 Q= mod (Q, M);
3527 R= mod (R, M);
3528 return;
3529 }
3530 int m= (int) ceil ((double) (degB + 1)/2.0) + 1;
3531 ASSERT (4*m >= degA, "expected degree (F, 1) < 2*degree (G, 1)");
3532 CFList splitA= split (A, m, x);
3533 if (splitA.length() == 3)
3534 splitA.insert (0);
3535 if (splitA.length() == 2)
3536 {
3537 splitA.insert (0);
3538 splitA.insert (0);
3539 }
3540 if (splitA.length() == 1)
3541 {
3542 splitA.insert (0);
3543 splitA.insert (0);
3544 splitA.insert (0);
3545 }
3546
3547 CanonicalForm xToM= power (x, m);
3548
3549 CFListIterator i= splitA;
3550 CanonicalForm H= i.getItem();
3551 i++;
3552 H *= xToM;
3553 H += i.getItem();
3554 i++;
3555 H *= xToM;
3556 H += i.getItem();
3557 i++;
3558
3559 divrem32 (H, B, Q, R, M);
3560
3561 CFList splitR= split (R, m, x);
3562 if (splitR.length() == 1)
3563 splitR.insert (0);
3564
3565 H= splitR.getFirst();
3566 H *= xToM;
3567 H += splitR.getLast();
3568 H *= xToM;
3569 H += i.getItem();
3570
3571 CanonicalForm bufQ;
3572 divrem32 (H, B, bufQ, R, M);
3573
3574 Q *= xToM;
3575 Q += bufQ;
3576 return;
3577}
3578
3579static inline
3581 CanonicalForm& R, const CFList& M)
3582{
3583 CanonicalForm A= mod (F, M);
3584 CanonicalForm B= mod (G, M);
3585 Variable x= Variable (1);
3586 int degB= degree (B, x);
3587 int degA= degree (A, x);
3588 if (degA < degB)
3589 {
3590 Q= 0;
3591 R= A;
3592 return;
3593 }
3594 if (degB < 1)
3595 {
3596 divrem (A, B, Q, R);
3597 Q= mod (Q, M);
3598 R= mod (R, M);
3599 return;
3600 }
3601 int m= (int) ceil ((double) (degB + 1)/ 2.0);
3602 ASSERT (3*m > degA, "expected degree (F, 1) < 3*degree (G, 1)");
3603 CFList splitA= split (A, m, x);
3604 CFList splitB= split (B, m, x);
3605
3606 if (splitA.length() == 2)
3607 {
3608 splitA.insert (0);
3609 }
3610 if (splitA.length() == 1)
3611 {
3612 splitA.insert (0);
3613 splitA.insert (0);
3614 }
3615 CanonicalForm xToM= power (x, m);
3616
3618 CFListIterator i= splitA;
3619 i++;
3620
3621 if (degree (splitA.getFirst(), x) < degree (splitB.getFirst(), x))
3622 {
3623 H= splitA.getFirst()*xToM + i.getItem();
3624 divrem21 (H, splitB.getFirst(), Q, R, M);
3625 }
3626 else
3627 {
3628 R= splitA.getFirst()*xToM + i.getItem() + splitB.getFirst() -
3629 splitB.getFirst()*xToM;
3630 Q= xToM - 1;
3631 }
3632
3633 H= mulMod (Q, splitB.getLast(), M);
3634
3635 R= R*xToM + splitA.getLast() - H;
3636
3637 while (degree (R, x) >= degB)
3638 {
3639 xToM= power (x, degree (R, x) - degB);
3640 Q += LC (R, x)*xToM;
3641 R -= mulMod (LC (R, x), B, M)*xToM;
3642 Q= mod (Q, M);
3643 R= mod (R, M);
3644 }
3645
3646 return;
3647}
3648
3650 CanonicalForm& R, const CanonicalForm& M)
3651{
3652 CanonicalForm A= mod (F, M);
3653 CanonicalForm B= mod (G, M);
3654
3655 if (B.inCoeffDomain())
3656 {
3657 divrem (A, B, Q, R);
3658 return;
3659 }
3660 if (A.inCoeffDomain() && !B.inCoeffDomain())
3661 {
3662 Q= 0;
3663 R= A;
3664 return;
3665 }
3666
3667 if (B.level() < A.level())
3668 {
3669 divrem (A, B, Q, R);
3670 return;
3671 }
3672 if (A.level() > B.level())
3673 {
3674 R= A;
3675 Q= 0;
3676 return;
3677 }
3678 if (B.level() == 1 && B.isUnivariate())
3679 {
3680 divrem (A, B, Q, R);
3681 return;
3682 }
3683
3684 Variable x= Variable (1);
3685 int degB= degree (B, x);
3686 if (degB > degree (A, x))
3687 {
3688 Q= 0;
3689 R= A;
3690 return;
3691 }
3692
3693 CFList splitA= split (A, degB, x);
3694
3695 CanonicalForm xToDegB= power (x, degB);
3696 CanonicalForm H, bufQ;
3697 Q= 0;
3698 CFListIterator i= splitA;
3699 H= i.getItem()*xToDegB;
3700 i++;
3701 H += i.getItem();
3702 CFList buf;
3703 while (i.hasItem())
3704 {
3705 buf= CFList (M);
3706 divrem21 (H, B, bufQ, R, buf);
3707 i++;
3708 if (i.hasItem())
3709 H= R*xToDegB + i.getItem();
3710 Q *= xToDegB;
3711 Q += bufQ;
3712 }
3713 return;
3714}
3715
3717 CanonicalForm& R, const CFList& MOD)
3718{
3719 CanonicalForm A= mod (F, MOD);
3720 CanonicalForm B= mod (G, MOD);
3721 Variable x= Variable (1);
3722 int degB= degree (B, x);
3723 if (degB > degree (A, x))
3724 {
3725 Q= 0;
3726 R= A;
3727 return;
3728 }
3729
3730 if (degB <= 0)
3731 {
3732 divrem (A, B, Q, R);
3733 Q= mod (Q, MOD);
3734 R= mod (R, MOD);
3735 return;
3736 }
3737 CFList splitA= split (A, degB, x);
3738
3739 CanonicalForm xToDegB= power (x, degB);
3740 CanonicalForm H, bufQ;
3741 Q= 0;
3742 CFListIterator i= splitA;
3743 H= i.getItem()*xToDegB;
3744 i++;
3745 H += i.getItem();
3746 while (i.hasItem())
3747 {
3748 divrem21 (H, B, bufQ, R, MOD);
3749 i++;
3750 if (i.hasItem())
3751 H= R*xToDegB + i.getItem();
3752 Q *= xToDegB;
3753 Q += bufQ;
3754 }
3755 return;
3756}
3757
3758bool
3760{
3761 if (B.isZero())
3762 return true;
3763 if (A.isZero())
3764 return false;
3766 return fdivides (A, B);
3767 int p= getCharacteristic();
3768 if (A.inCoeffDomain() || B.inCoeffDomain())
3769 {
3770 if (A.inCoeffDomain())
3771 return true;
3772 else
3773 return false;
3774 }
3775 if (p > 0)
3776 {
3777#if (!defined(HAVE_FLINT) || __FLINT_RELEASE < 20400)
3778 if (fac_NTL_char != p)
3779 {
3780 fac_NTL_char= p;
3781 zz_p::init (p);
3782 }
3783#endif
3786 {
3787#if (HAVE_FLINT && __FLINT_RELEASE >= 20400)
3788 nmod_poly_t FLINTmipo;
3789 fq_nmod_ctx_t fq_con;
3790
3791 nmod_poly_init (FLINTmipo, getCharacteristic());
3792 convertFacCF2nmod_poly_t (FLINTmipo, getMipo (alpha));
3793
3794 fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
3795
3796 fq_nmod_poly_t FLINTA, FLINTB;
3799 int result= fq_nmod_poly_divides (FLINTA, FLINTB, FLINTA, fq_con);
3800 fq_nmod_poly_clear (FLINTA, fq_con);
3801 fq_nmod_poly_clear (FLINTB, fq_con);
3802 nmod_poly_clear (FLINTmipo);
3804 return result;
3805#else
3806 zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
3807 zz_pE::init (NTLMipo);
3808 zz_pEX NTLA= convertFacCF2NTLzz_pEX (A, NTLMipo);
3809 zz_pEX NTLB= convertFacCF2NTLzz_pEX (B, NTLMipo);
3810 return divide (NTLB, NTLA);
3811#endif
3812 }
3813#ifdef HAVE_FLINT
3814 nmod_poly_t FLINTA, FLINTB;
3815 convertFacCF2nmod_poly_t (FLINTA, A);
3816 convertFacCF2nmod_poly_t (FLINTB, B);
3817 nmod_poly_divrem (FLINTB, FLINTA, FLINTB, FLINTA);
3818 bool result= nmod_poly_is_zero (FLINTA);
3819 nmod_poly_clear (FLINTA);
3820 nmod_poly_clear (FLINTB);
3821 return result;
3822#else
3823 zz_pX NTLA= convertFacCF2NTLzzpX (A);
3824 zz_pX NTLB= convertFacCF2NTLzzpX (B);
3825 return divide (NTLB, NTLA);
3826#endif
3827 }
3828#ifdef HAVE_FLINT
3830 bool isRat= isOn (SW_RATIONAL);
3831 if (!isRat)
3832 On (SW_RATIONAL);
3833 if (!hasFirstAlgVar (A, alpha) && !hasFirstAlgVar (B, alpha))
3834 {
3835 fmpq_poly_t FLINTA,FLINTB;
3836 convertFacCF2Fmpq_poly_t (FLINTA, A);
3837 convertFacCF2Fmpq_poly_t (FLINTB, B);
3838 fmpq_poly_rem (FLINTA, FLINTB, FLINTA);
3839 bool result= fmpq_poly_is_zero (FLINTA);
3840 fmpq_poly_clear (FLINTA);
3841 fmpq_poly_clear (FLINTB);
3842 if (!isRat)
3843 Off (SW_RATIONAL);
3844 return result;
3845 }
3846 CanonicalForm Q, R;
3847 newtonDivrem (B, A, Q, R);
3848 if (!isRat)
3849 Off (SW_RATIONAL);
3850 return R.isZero();
3851#else
3852 bool isRat= isOn (SW_RATIONAL);
3853 if (!isRat)
3854 On (SW_RATIONAL);
3855 bool result= fdivides (A, B);
3856 if (!isRat)
3857 Off (SW_RATIONAL);
3858 return result; //maybe NTL?
3859#endif
3860}
3861
3862// end division
3863
3864#else
3866mulNTL (const CanonicalForm& F, const CanonicalForm& G, const modpk& b)
3867{ return F*G; }
3868#endif
CanonicalForm convertFq_poly_t2FacCF(const fq_poly_t p, const Variable &x, const Variable &alpha, const fq_ctx_t ctx)
conversion of a FLINT poly over Fq (for non-word size p) to a CanonicalForm with alg....
void convertFacCF2Fq_t(fq_t result, const CanonicalForm &f, const fq_ctx_t ctx)
conversion of a factory element of F_q (for non-word size p) to a FLINT fq_t
CanonicalForm convertFq_nmod_poly_t2FacCF(const fq_nmod_poly_t p, const Variable &x, const Variable &alpha, const fq_nmod_ctx_t ctx)
conversion of a FLINT poly over Fq to a CanonicalForm with alg. variable alpha and polynomial variabl...
CanonicalForm convertFq_t2FacCF(const fq_t poly, const Variable &alpha)
conversion of a FLINT element of F_q with non-word size p to a CanonicalForm with alg....
CanonicalForm convertFmpq_poly_t2FacCF(const fmpq_poly_t p, const Variable &x)
conversion of a FLINT poly over Q to CanonicalForm
CanonicalForm convertFmpz_mod_poly_t2FacCF(const fmpz_mod_poly_t poly, const Variable &x, const modpk &b)
conversion of a FLINT poly over Z/p (for non word size p) to a CanonicalForm over Z
CanonicalForm convertnmod_poly_t2FacCF(const nmod_poly_t poly, const Variable &x)
conversion of a FLINT poly over Z/p to CanonicalForm
void convertFacCF2Fmpz_mod_poly_t(fmpz_mod_poly_t result, const CanonicalForm &f, const fmpz_t p)
conversion of a factory univariate poly over Z to a FLINT poly over Z/p (for non word size p)
void convertFacCF2Fq_nmod_poly_t(fq_nmod_poly_t result, const CanonicalForm &f, const fq_nmod_ctx_t ctx)
conversion of a factory univariate poly over F_q to a FLINT fq_nmod_poly_t
CanonicalForm convertFmpz_poly_t2FacCF(const fmpz_poly_t poly, const Variable &x)
conversion of a FLINT poly over Z to CanonicalForm
void convertFacCF2Fmpq_poly_t(fmpq_poly_t result, const CanonicalForm &f)
conversion of a factory univariate polynomials over Q to fmpq_poly_t
void convertFacCF2Fmpz_poly_t(fmpz_poly_t result, const CanonicalForm &f)
conversion of a factory univariate polynomial over Z to a fmpz_poly_t
void convertCF2initFmpz(fmpz_t result, const CanonicalForm &f)
conversion of a factory integer to fmpz_t(init.)
void convertFacCF2Fq_poly_t(fq_poly_t result, const CanonicalForm &f, const fq_ctx_t ctx)
conversion of a factory univariate poly over F_q (for non-word size p) to a FLINT fq_poly_t
This file defines functions for conversion to FLINT (www.flintlib.org) and back.
ZZX convertFacCF2NTLZZX(const CanonicalForm &f)
Definition: NTLconvert.cc:691
zz_pEX convertFacCF2NTLzz_pEX(const CanonicalForm &f, const zz_pX &mipo)
Definition: NTLconvert.cc:1064
CanonicalForm convertNTLzz_pEX2CF(const zz_pEX &f, const Variable &x, const Variable &alpha)
Definition: NTLconvert.cc:1092
ZZ_pEX convertFacCF2NTLZZ_pEX(const CanonicalForm &f, const ZZ_pX &mipo)
CanonicalForm in Z_p(a)[X] to NTL ZZ_pEX.
Definition: NTLconvert.cc:1037
CanonicalForm convertNTLzzpX2CF(const zz_pX &poly, const Variable &x)
Definition: NTLconvert.cc:255
CanonicalForm convertNTLZZpX2CF(const ZZ_pX &poly, const Variable &x)
NAME: convertNTLZZpX2CF.
Definition: NTLconvert.cc:250
CanonicalForm convertNTLZZX2CF(const ZZX &polynom, const Variable &x)
Definition: NTLconvert.cc:285
CanonicalForm convertNTLZZ_pEX2CF(const ZZ_pEX &f, const Variable &x, const Variable &alpha)
Definition: NTLconvert.cc:1115
zz_pX convertFacCF2NTLzzpX(const CanonicalForm &f)
Definition: NTLconvert.cc:105
ZZ_pX convertFacCF2NTLZZpX(const CanonicalForm &f)
NAME: convertFacCF2NTLZZpX.
Definition: NTLconvert.cc:64
VAR long fac_NTL_char
Definition: NTLconvert.cc:46
ZZ convertFacCF2NTLZZ(const CanonicalForm &f)
NAME: convertFacCF2NTLZZX.
Definition: NTLconvert.cc:670
Conversion to and from NTL.
#define swap(_i, _j)
bool isOn(int sw)
switches
void On(int sw)
switches
void Off(int sw)
switches
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
Header for factory's main class CanonicalForm.
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition: cf_ops.cc:600
CF_NO_INLINE FACTORY_PUBLIC CanonicalForm div(const CanonicalForm &, const CanonicalForm &)
int ilog2(const CanonicalForm &a)
CanonicalForm tailcoeff(const CanonicalForm &f)
int taildegree(const CanonicalForm &f)
int degree(const CanonicalForm &f)
bool hasFirstAlgVar(const CanonicalForm &f, Variable &a)
check if poly f contains an algebraic variable a
Definition: cf_ops.cc:679
CanonicalForm FACTORY_PUBLIC swapvar(const CanonicalForm &, const Variable &, const Variable &)
swapvar() - swap variables x1 and x2 in f.
Definition: cf_ops.cc:168
CanonicalForm den(const CanonicalForm &f)
List< CanonicalForm > CFList
CanonicalForm LC(const CanonicalForm &f)
int FACTORY_PUBLIC getCharacteristic()
Definition: cf_char.cc:70
int l
Definition: cfEzgcd.cc:100
int m
Definition: cfEzgcd.cc:128
int i
Definition: cfEzgcd.cc:132
int k
Definition: cfEzgcd.cc:99
Variable x
Definition: cfModGcd.cc:4082
int p
Definition: cfModGcd.cc:4078
g
Definition: cfModGcd.cc:4090
CanonicalForm fp
Definition: cfModGcd.cc:4102
CanonicalForm cd(bCommonDen(FF))
Definition: cfModGcd.cc:4089
CanonicalForm gp
Definition: cfModGcd.cc:4102
CanonicalForm b
Definition: cfModGcd.cc:4103
CanonicalForm bCommonDen(const CanonicalForm &f)
CanonicalForm bCommonDen ( const CanonicalForm & f )
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
bool fdivides ( const CanonicalForm & f, const CanonicalForm & g )
declarations of higher level algorithms.
#define ASSERT(expression, message)
Definition: cf_assert.h:99
static const int SW_RATIONAL
set to 1 for computations over Q
Definition: cf_defs.h:31
#define GaloisFieldDomain
Definition: cf_defs.h:18
Iterators for CanonicalForm's.
FILE * f
Definition: checklibs.c:9
static int gettype()
Definition: cf_factory.h:28
class to iterate through CanonicalForm's
Definition: cf_iter.h:44
factory's main class
Definition: canonicalform.h:86
CF_NO_INLINE bool isZero() const
Variable mvar() const
mvar() returns the main variable of CO or Variable() if CO is in a base domain.
bool inCoeffDomain() const
int level() const
level() returns the level of CO.
bool inBaseDomain() const
bool isUnivariate() const
T getFirst() const
Definition: ftmpl_list.cc:279
int length() const
Definition: ftmpl_list.cc:273
void append(const T &)
Definition: ftmpl_list.cc:256
T getLast() const
Definition: ftmpl_list.cc:309
void insert(const T &)
Definition: ftmpl_list.cc:193
int isEmpty() const
Definition: ftmpl_list.cc:267
factory's class for variables
Definition: variable.h:33
int level() const
Definition: variable.h:49
class to do operations mod p^k for int's p and k
Definition: fac_util.h:23
functions to print debug output
#define DEBOUTLN(stream, objects)
Definition: debug.h:49
Variable alpha
Definition: facAbsBiFact.cc:51
return result
Definition: facAbsBiFact.cc:75
const CanonicalForm int const CFList const Variable & y
Definition: facAbsFact.cc:53
CanonicalForm H
Definition: facAbsFact.cc:60
int degg
Definition: facAlgExt.cc:64
CanonicalForm mipo
Definition: facAlgExt.cc:57
CanonicalForm divide(const CanonicalForm &ff, const CanonicalForm &f, const CFList &as)
b *CanonicalForm B
Definition: facBivar.cc:52
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
CFList tmp1
Definition: facFqBivar.cc:72
CanonicalForm buf2
Definition: facFqBivar.cc:73
CFList tmp2
Definition: facFqBivar.cc:72
CanonicalForm buf1
Definition: facFqBivar.cc:73
fq_nmod_ctx_t fq_con
Definition: facHensel.cc:99
int j
Definition: facHensel.cc:110
fq_nmod_ctx_clear(fq_con)
nmod_poly_init(FLINTmipo, getCharacteristic())
fq_nmod_ctx_init_modulus(fq_con, FLINTmipo, "Z")
fq_nmod_poly_init(prod, fq_con)
convertFacCF2nmod_poly_t(FLINTmipo, M)
nmod_poly_clear(FLINTmipo)
fq_nmod_poly_clear(prod, fq_con)
CanonicalForm mod(const CanonicalForm &F, const CFList &M)
reduce F modulo elements in M.
Definition: facMul.cc:3072
CanonicalForm uniReverse(const CanonicalForm &F, int d, const Variable &x)
Definition: facMul.cc:274
void newtonDivrem(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q, CanonicalForm &R)
division with remainder of univariate polynomials over Q and Q(a) using Newton inversion,...
Definition: facMul.cc:346
void kronSubFq(fq_nmod_poly_t result, const CanonicalForm &A, int d, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:1271
CanonicalForm mulNTL(const CanonicalForm &F, const CanonicalForm &G, const modpk &b)
multiplication of univariate polys using FLINT/NTL over F_p, F_q, Z/p^k, Z/p^k[t]/(f),...
Definition: facMul.cc:411
void divrem(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q, CanonicalForm &R, const CFList &MOD)
division with remainder of F by G wrt Variable (1) modulo MOD. Uses an algorithm based on Burnikel,...
Definition: facMul.cc:3716
bool uniFdivides(const CanonicalForm &A, const CanonicalForm &B)
divisibility test for univariate polys
Definition: facMul.cc:3759
CanonicalForm divFLINTQ(const CanonicalForm &F, const CanonicalForm &G)
Definition: facMul.cc:174
void kronSubQa(fmpz_poly_t result, const CanonicalForm &A, int d)
Definition: facMul.cc:46
CanonicalForm reverseSubstFp(const nmod_poly_t F, int d)
Definition: facMul.cc:2054
static CFList split(const CanonicalForm &F, const int m, const Variable &x)
Definition: facMul.cc:3469
static void divrem32(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q, CanonicalForm &R, const CFList &M)
Definition: facMul.cc:3580
CanonicalForm mulMod2FLINTQ(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2272
CanonicalForm reverse(const CanonicalForm &F, int d)
Definition: facMul.cc:3234
CanonicalForm mulMod2(const CanonicalForm &A, const CanonicalForm &B, const CanonicalForm &M)
Karatsuba style modular multiplication for bivariate polynomials.
Definition: facMul.cc:2986
CanonicalForm mulMod2FLINTFqReci(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M, const Variable &alpha, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:2160
CanonicalForm mulFLINTQ(const CanonicalForm &F, const CanonicalForm &G)
Definition: facMul.cc:133
CanonicalForm mulMod2FLINTFpReci(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2090
CanonicalForm mulMod2FLINTFq(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M, const Variable &alpha, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:2202
CanonicalForm reverseSubstReciproFq(const fq_nmod_poly_t F, const fq_nmod_poly_t G, int d, int k, const Variable &alpha, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:1781
CanonicalForm modFLINTQ(const CanonicalForm &F, const CanonicalForm &G)
Definition: facMul.cc:192
void kronSubReciproQ(fmpz_poly_t subA1, fmpz_poly_t subA2, const CanonicalForm &A, int d)
Definition: facMul.cc:1474
void kronSubReciproFq(fq_nmod_poly_t subA1, fq_nmod_poly_t subA2, const CanonicalForm &A, int d, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:1429
CanonicalForm reverseSubstQ(const fmpz_poly_t F, int d)
Definition: facMul.cc:1499
CanonicalForm mulMod2FLINTQReci(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2235
CanonicalForm reverseSubstFq(const fq_nmod_poly_t F, int d, const Variable &alpha, const fq_nmod_ctx_t fq_con)
Definition: facMul.cc:2019
CanonicalForm mulMod(const CanonicalForm &A, const CanonicalForm &B, const CFList &MOD)
Karatsuba style modular multiplication for multivariate polynomials.
Definition: facMul.cc:3080
CanonicalForm mulFLINTQTrunc(const CanonicalForm &F, const CanonicalForm &G, int m)
Definition: facMul.cc:241
CanonicalForm mulFLINTQa(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha)
Definition: facMul.cc:103
CanonicalForm reverseSubstReciproQ(const fmpz_poly_t F, const fmpz_poly_t G, int d, int k)
Definition: facMul.cc:1885
static void divrem21(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q, CanonicalForm &R, const CFList &M)
Definition: facMul.cc:3509
CanonicalForm newtonInverse(const CanonicalForm &F, const int n, const Variable &x)
Definition: facMul.cc:291
void newtonDiv(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q)
Definition: facMul.cc:380
void divrem2(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Q, CanonicalForm &R, const CanonicalForm &M)
division with remainder of F by G wrt Variable (1) modulo M. Uses an algorithm based on Burnikel,...
Definition: facMul.cc:3649
CanonicalForm reverseSubstQa(const fmpz_poly_t F, int d, const Variable &x, const Variable &alpha, const CanonicalForm &den)
Definition: facMul.cc:66
void kronSubReciproFp(nmod_poly_t subA1, nmod_poly_t subA2, const CanonicalForm &A, int d)
Definition: facMul.cc:1388
CanonicalForm divNTL(const CanonicalForm &F, const CanonicalForm &G, const modpk &b)
division of univariate polys using FLINT/NTL over F_p, F_q, Z/p^k, Z/p^k[t]/(f), Z,...
Definition: facMul.cc:936
CanonicalForm mulFLINTQaTrunc(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, int m)
Definition: facMul.cc:210
CanonicalForm modNTL(const CanonicalForm &F, const CanonicalForm &G, const modpk &b)
mod of univariate polys using FLINT/NTL over F_p, F_q, Z/p^k, Z/p^k[t]/(f), Z, Q, Q(a),...
Definition: facMul.cc:731
CanonicalForm prodMod(const CFList &L, const CanonicalForm &M)
product of all elements in L modulo M via divide-and-conquer.
Definition: facMul.cc:3180
CanonicalForm reverseSubstReciproFp(const nmod_poly_t F, const nmod_poly_t G, int d, int k)
Definition: facMul.cc:1662
void kronSubFp(nmod_poly_t result, const CanonicalForm &A, int d)
Definition: facMul.cc:1248
CanonicalForm mulMod2FLINTFp(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2128
CanonicalForm mulMod2NTLFq(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2926
CanonicalForm mulMod2FLINTQa(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M)
Definition: facMul.cc:2332
This file defines functions for fast multiplication and division with remainder.
some useful template functions.
template CanonicalForm tmax(const CanonicalForm &, const CanonicalForm &)
template CanonicalForm tmin(const CanonicalForm &, const CanonicalForm &)
template List< Variable > Difference(const List< Variable > &, const List< Variable > &)
STATIC_VAR TreeM * G
Definition: janet.cc:31
STATIC_VAR Poly * h
Definition: janet.cc:971
void rem(unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
Definition: minpoly.cc:572
gmp_float exp(const gmp_float &a)
Definition: mpr_complex.cc:357
The main handler for Singular numbers which are suitable for Singular polynomials.
int status int void * buf
Definition: si_signals.h:59
#define R
Definition: sirandom.c:27
#define A
Definition: sirandom.c:24
#define M
Definition: sirandom.c:25
#define Q
Definition: sirandom.c:26
int F1(int a1, int &r1)
F1.
void F2(int a2, int &r2)
F2.
bool getReduce(const Variable &alpha)
Definition: variable.cc:232
CanonicalForm getMipo(const Variable &alpha, const Variable &x)
Definition: variable.cc:207