1 /*
2     Copyright (C) 2019 Daniel Schultz
3 
4     This file is part of FLINT.
5 
6     FLINT is free software: you can redistribute it and/or modify it under
7     the terms of the GNU Lesser General Public License (LGPL) as published
8     by the Free Software Foundation; either version 2.1 of the License, or
9     (at your option) any later version.  See <http://www.gnu.org/licenses/>.
10 */
11 
12 #include "fq_nmod_mpoly.h"
13 
14 /*
15     interp_reduce: map from Fp[x] to Fp[x]/poly(x)
16     interp_lift:   map from Fp[x]/poly(x) to Fp[x]
17     interp_crt:    update element of Fp[x] with a new image in Fp[x]/poly(x)
18     interp_mcrt:   same as interp_mcrt, but monomial match, thus easier
19 
20     ..._sm:    poly(x) is x - alpha
21     ..._lg:    poly(x) is modulus of finite field
22 */
23 
24 /****************** bivar - image in F_q *************************************/
25 
26 /*
27     E = A mod alpha(v)
28     A is in Fp[X][][v]
29     E is in (Fp/alpha(v))[X]
30 */
31 
nmod_mpolyn_interp_reduce_lg_poly(fq_nmod_poly_t E,const fq_nmod_ctx_t fqctx,const nmod_mpolyn_t A,const nmod_mpoly_ctx_t ctx)32 void nmod_mpolyn_interp_reduce_lg_poly(
33     fq_nmod_poly_t E,
34     const fq_nmod_ctx_t fqctx,
35     const nmod_mpolyn_t A,
36     const nmod_mpoly_ctx_t ctx)
37 {
38     fq_nmod_t v;
39     slong Ai, Alen, k;
40     nmod_poly_struct * Acoeff;
41     ulong * Aexp;
42     slong N, off, shift;
43 
44     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
45     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
46 
47     fq_nmod_init(v, fqctx);
48 
49     Acoeff = A->coeffs;
50     Aexp = A->exps;
51     Alen = A->length;
52 
53     fq_nmod_poly_zero(E, fqctx);
54     for (Ai = 0; Ai < Alen; Ai++)
55     {
56         k = (Aexp + N*Ai)[off] >> shift;
57         nmod_poly_rem(v, Acoeff + Ai, fqctx->modulus);
58         fq_nmod_poly_set_coeff(E, k, v, fqctx);
59     }
60 
61     fq_nmod_clear(v, fqctx);
62 }
63 
64 
65 /*
66     Convert B to A using the lowest degree representative
67     A is in           Fp [X][][v]
68     B is in (Fp/alpha(v))[X]
69 */
70 
nmod_mpolyn_interp_lift_lg_poly(slong * lastdeg_,nmod_mpolyn_t A,const nmod_mpoly_ctx_t ctx,const fq_nmod_poly_t B,const fq_nmod_ctx_t fqctx)71 void nmod_mpolyn_interp_lift_lg_poly(
72     slong * lastdeg_,
73     nmod_mpolyn_t A,
74     const nmod_mpoly_ctx_t ctx,
75     const fq_nmod_poly_t B,
76     const fq_nmod_ctx_t fqctx)
77 {
78     slong Bexp;
79     slong Blen = fq_nmod_poly_length(B, fqctx);
80     fq_nmod_struct * Bcoeff = B->coeffs;
81     nmod_poly_struct * Acoeff;
82     ulong * Aexp;
83     slong Ai;
84     slong lastdeg = -WORD(1);
85     slong N, off, shift;
86 
87     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
88     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
89 
90     nmod_mpolyn_fit_length(A, Blen, ctx);
91     Acoeff = A->coeffs;
92     Aexp = A->exps;
93 
94     Ai = 0;
95     for (Bexp = Blen - 1; Bexp >= 0; Bexp--)
96     {
97         if (!fq_nmod_is_zero(Bcoeff + Bexp, fqctx))
98         {
99             FLINT_ASSERT(Ai < A->alloc);
100             nmod_poly_set(Acoeff + Ai, Bcoeff + Bexp);
101             lastdeg = FLINT_MAX(lastdeg, nmod_poly_degree(Bcoeff + Bexp));
102             mpoly_monomial_zero(Aexp + N*Ai, N);
103             (Aexp + N*Ai)[off] = Bexp << shift;
104             Ai++;
105         }
106     }
107     A->length = Ai;
108 
109     *lastdeg_ = lastdeg;
110 }
111 
112 
113 /*
114     F = F + modulus*((A - F(alpha))/(modulus(alpha)))
115     no assumptions about matching monomials
116     F is in Fp[X][v]
117     A is in (Fp/alpha(v))[X]    alpha(v) is fqctx->modulus(v)
118 */
119 
nmod_mpolyn_interp_crt_lg_poly(slong * lastdeg_,nmod_mpolyn_t F,nmod_mpolyn_t T,nmod_poly_t modulus,const nmod_mpoly_ctx_t ctx,fq_nmod_poly_t A,const fq_nmod_ctx_t fqctx)120 int nmod_mpolyn_interp_crt_lg_poly(
121     slong * lastdeg_,
122     nmod_mpolyn_t F,
123     nmod_mpolyn_t T,
124     nmod_poly_t modulus,
125     const nmod_mpoly_ctx_t ctx,
126     fq_nmod_poly_t A,
127     const fq_nmod_ctx_t fqctx)
128 {
129     int changed = 0;
130     slong lastdeg = -WORD(1);
131     fq_nmod_t u, v;
132     nmod_poly_t w;
133     slong Fi, Ti, Aexp;
134     fq_nmod_struct * Acoeff = A->coeffs;
135     slong Flen = F->length;
136     nmod_poly_struct * Fcoeffs = F->coeffs;
137     ulong * Fexps = F->exps;
138     nmod_poly_struct * Tcoeffs;
139     ulong * Texps;
140     nmod_poly_t tp;
141     fq_nmod_t inv_m_eval;
142     slong N, off, shift;
143 
144     FLINT_ASSERT(T->bits == F->bits);
145 
146     N = mpoly_words_per_exp_sp(F->bits, ctx->minfo);
147     mpoly_gen_offset_shift_sp(&off, &shift, 0, F->bits, ctx->minfo);
148 
149     fq_nmod_init(inv_m_eval, fqctx);
150     nmod_poly_rem(inv_m_eval, modulus, fqctx->modulus);
151     fq_nmod_inv(inv_m_eval, inv_m_eval, fqctx);
152 
153     Fi = 0;
154 
155     fq_nmod_init(u, fqctx);
156     fq_nmod_init(v, fqctx);
157     nmod_poly_init(w, fqctx->modulus->mod.n);
158 
159     Aexp = fq_nmod_poly_degree(A, fqctx);
160 
161     nmod_poly_init(tp, ctx->ffinfo->mod.n);
162 
163     nmod_mpolyn_fit_length(T, Flen + Aexp + 1, ctx);
164     Tcoeffs = T->coeffs;
165     Texps = T->exps;
166     Ti = 0;
167 
168     while (Fi < Flen || Aexp >= 0)
169     {
170         FLINT_ASSERT(Ti < T->alloc);
171 
172         if (Fi < Flen)
173         {
174             FLINT_ASSERT(!nmod_poly_is_zero(Fcoeffs + Fi));
175             FLINT_ASSERT(nmod_poly_degree(Fcoeffs + Fi) < nmod_poly_degree(modulus));
176         }
177 
178         if (Aexp >= 0)
179         {
180             FLINT_ASSERT(!fq_nmod_is_zero(Acoeff + Aexp, fqctx));
181         }
182 
183         mpoly_monomial_zero(Texps + N*Ti, N);
184 
185         if (Fi < Flen && Aexp >= 0 && ((Fexps + N*Fi)[off]>>shift) == Aexp)
186         {
187             /* F term ok, A term ok */
188             nmod_poly_rem(u, Fcoeffs + Fi, fqctx->modulus);
189             fq_nmod_sub(v, Acoeff + Aexp, u, fqctx);
190             if (!fq_nmod_is_zero(v, fqctx))
191             {
192                 changed = 1;
193                 fq_nmod_mul(u, v, inv_m_eval, fqctx);
194                 nmod_poly_mul(w, modulus, u);
195                 nmod_poly_add(Tcoeffs + Ti, Fcoeffs + Fi, w);
196             }
197             else
198             {
199                 nmod_poly_set(Tcoeffs + Ti, Fcoeffs + Fi);
200             }
201 
202             (Texps + N*Ti)[off] = (Fexps + N*Fi)[off];
203             Fi++;
204             do {
205                 Aexp--;
206             } while (Aexp >= 0 && fq_nmod_is_zero(Acoeff + Aexp, fqctx));
207         }
208         else if (Fi < Flen && (Aexp < 0 || ((Fexps + N*Fi)[off]>>shift) > Aexp))
209         {
210             /* F term ok, A term missing */
211             nmod_poly_rem(v, Fcoeffs + Fi, fqctx->modulus);
212             if (!fq_nmod_is_zero(v, fqctx))
213             {
214                 changed = 1;
215                 fq_nmod_mul(u, v, inv_m_eval, fqctx);
216                 nmod_poly_mul(w, u, modulus);
217                 nmod_poly_sub(Tcoeffs + 0, Fcoeffs + 0, w);
218             }
219             else
220             {
221                 nmod_poly_set(Tcoeffs + Ti, Fcoeffs + Fi);
222             }
223 
224             (Texps + N*Ti)[off] = (Fexps + N*Fi)[off];
225             Fi++;
226         }
227         else if (Aexp >= 0 && (Fi >= Flen || ((Fexps + N*Fi)[off]>>shift) < Aexp))
228         {
229             /* F term missing, A term ok */
230             changed = 1;
231             fq_nmod_mul(u, Acoeff + Aexp, inv_m_eval, fqctx);
232             nmod_poly_mul(Tcoeffs + Ti, modulus, u);
233             (Texps + N*Ti)[off] = Aexp << shift;
234             do {
235                 Aexp--;
236             } while (Aexp >= 0 && fq_nmod_is_zero(Acoeff + Aexp, fqctx));
237         }
238         else
239         {
240             FLINT_ASSERT(0);
241         }
242 
243         FLINT_ASSERT(!nmod_poly_is_zero(Tcoeffs + Ti));
244         lastdeg = FLINT_MAX(lastdeg, nmod_poly_degree(Tcoeffs + Ti));
245         Ti++;
246     }
247     T->length = Ti;
248 
249     if (changed)
250     {
251         nmod_mpolyn_swap(T, F);
252     }
253 
254     fq_nmod_clear(u, fqctx);
255     fq_nmod_clear(v, fqctx);
256     nmod_poly_clear(w);
257 
258     fq_nmod_clear(inv_m_eval, fqctx);
259 
260     *lastdeg_ = lastdeg;
261     return changed;
262 }
263 
264 
265 /*********************** multivar - image in F_q *****************************/
266 
267 /*
268     E = A mod alpha(x_var)
269     A is in                Fp[x_0, ..., x_(var-2), x_(var-1)][x_var]
270     E is in (Fp/alpha(x_var))[x_0, ..., x_(var-2)][x_(var-1)]
271     alpha is ectx->modulus
272 */
273 
nmod_mpolyn_interp_reduce_lg_mpolyn(fq_nmod_mpolyn_t E,fq_nmod_mpoly_ctx_t ectx,nmod_mpolyn_t A,slong var,const nmod_mpoly_ctx_t ctx)274 void nmod_mpolyn_interp_reduce_lg_mpolyn(
275     fq_nmod_mpolyn_t E,
276     fq_nmod_mpoly_ctx_t ectx,
277     nmod_mpolyn_t A,
278     slong var,
279     const nmod_mpoly_ctx_t ctx)
280 {
281     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
282     slong offset, shift, k;
283     ulong mask;
284     fq_nmod_t v;
285     nmod_poly_struct * Acoeff = A->coeffs;
286     ulong * Aexp = A->exps;
287     slong Alen = A->length;
288     slong Ai;
289     fq_nmod_poly_struct * Ecoeff;
290     ulong * Eexp;
291     slong Ei;
292 
293     fq_nmod_init(v, ectx->fqctx);
294 
295     FLINT_ASSERT(var > 0);
296     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
297     mask = (-UWORD(1)) >> (FLINT_BITS - A->bits);
298 
299     Ecoeff = E->coeffs;
300     Eexp = E->exps;
301     Ei = 0;
302     for (Ai = 0; Ai < Alen; Ai++)
303     {
304         nmod_poly_rem(v, Acoeff + Ai, ectx->fqctx->modulus);
305         k = ((Aexp + N*Ai)[offset] >> shift) & mask;
306         if (fq_nmod_is_zero(v, ectx->fqctx))
307         {
308             continue;
309         }
310 
311         if (Ei > 0 && mpoly_monomial_equal_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)))
312         {
313             /* append to previous */
314             fq_nmod_poly_set_coeff(Ecoeff + Ei - 1, k, v, ectx->fqctx);
315         }
316         else
317         {
318             FLINT_ASSERT(Ei == 0 || mpoly_monomial_gt_nomask_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)));
319 
320             /* create new */
321             if (Ei >= E->alloc)
322             {
323                 fq_nmod_mpolyn_fit_length(E, Ei + 1, ectx);
324                 Ecoeff = E->coeffs;
325                 Eexp = E->exps;
326             }
327             mpoly_monomial_set_extra(Eexp + N*Ei, Aexp + N*Ai, N, offset, -(k << shift));
328             fq_nmod_poly_zero(Ecoeff + Ei, ectx->fqctx);
329             fq_nmod_poly_set_coeff(Ecoeff + Ei, k, v, ectx->fqctx);
330             Ei++;
331         }
332     }
333     E->length = Ei;
334 
335     fq_nmod_clear(v, ectx->fqctx);
336 }
337 
338 
339 /*
340     A = B using lowest degree representative
341     A is in                      Fp [x_0, ..., x_(var-1), x_(var-1)][x_var]
342     B is in (Fp[x_var]/alpha(x_var))[x_0, ..., x_(var-2)][x_(var-1)]
343 */
nmod_mpolyn_interp_lift_lg_mpolyn(slong * lastdeg_,nmod_mpolyn_t A,slong var,const nmod_mpoly_ctx_t ctx,fq_nmod_mpolyn_t B,const fq_nmod_mpoly_ctx_t ectx)344 void nmod_mpolyn_interp_lift_lg_mpolyn(
345     slong * lastdeg_,
346     nmod_mpolyn_t A,
347     slong var,
348     const nmod_mpoly_ctx_t ctx,
349     fq_nmod_mpolyn_t B,
350     const fq_nmod_mpoly_ctx_t ectx)
351 {
352     slong N = mpoly_words_per_exp_sp(B->bits, ctx->minfo);
353     slong offset, shift;
354     slong vi;
355     fq_nmod_poly_struct * Bcoeff = B->coeffs;
356     ulong * Bexp = B->exps;
357     slong Blen = B->length;
358     slong Bi;
359     nmod_poly_struct * Acoeff;
360     ulong * Aexp;
361     slong Ai;
362     slong lastdeg = -WORD(1);
363 
364     FLINT_ASSERT(A->bits == B->bits);
365 
366     nmod_mpolyn_fit_length(A, Blen, ctx);
367     Acoeff = A->coeffs;
368     Aexp = A->exps;
369 
370     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
371 
372     Ai = 0;
373     for (Bi = 0; Bi < Blen; Bi++)
374     {
375         if (Ai + (Bcoeff + Bi)->length >= A->alloc)
376         {
377             nmod_mpolyn_fit_length(A, Ai + (Bcoeff + Bi)->length, ctx);
378             Acoeff = A->coeffs;
379             Aexp = A->exps;
380         }
381         for (vi = (Bcoeff + Bi)->length - 1; vi >= 0; vi--)
382         {
383             if (!fq_nmod_is_zero((Bcoeff + Bi)->coeffs + vi, ectx->fqctx))
384             {
385                 mpoly_monomial_set_extra(Aexp + N*Ai, Bexp + N*Bi, N, offset, vi << shift);
386                 nmod_poly_set(Acoeff + Ai, (Bcoeff + Bi)->coeffs + vi);
387                 lastdeg = FLINT_MAX(lastdeg, nmod_poly_degree(Acoeff + Ai));
388                 Ai++;
389             }
390         }
391     }
392     A->length = Ai;
393 
394     *lastdeg_ = lastdeg;
395 }
396 
397 
nmod_mpolyn_interp_crt_lg_mpolyn(slong * lastdeg_,nmod_mpolyn_t F,nmod_mpolyn_t T,nmod_poly_t modulus,slong var,const nmod_mpoly_ctx_t ctx,fq_nmod_mpolyn_t A,const fq_nmod_mpoly_ctx_t ectx)398 int nmod_mpolyn_interp_crt_lg_mpolyn(
399     slong * lastdeg_,
400     nmod_mpolyn_t F,
401     nmod_mpolyn_t T,
402     nmod_poly_t modulus,
403     slong var,
404     const nmod_mpoly_ctx_t ctx,
405     fq_nmod_mpolyn_t A,
406     const fq_nmod_mpoly_ctx_t ectx)
407 {
408     int changed = 0;
409     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
410     slong offset, shift;
411     slong lastdeg = -WORD(1);
412     slong vi;
413     fq_nmod_t u, v;
414     nmod_poly_t w;
415 
416     nmod_poly_struct * Tcoeff;
417     ulong * Texp;
418     slong Ti;
419 
420     fq_nmod_poly_struct * Acoeff = A->coeffs;
421     slong Alen = A->length;
422     ulong * Aexp = A->exps;
423     slong Ai;
424 
425     nmod_poly_struct * Fcoeff = F->coeffs;
426     slong Flen = F->length;
427     ulong * Fexp = F->exps;
428     slong Fi;
429     fq_nmod_t inv_m_eval;
430 
431     fq_nmod_init(inv_m_eval, ectx->fqctx);
432     nmod_poly_rem(inv_m_eval, modulus, ectx->fqctx->modulus);
433     fq_nmod_inv(inv_m_eval, inv_m_eval, ectx->fqctx);
434 
435     fq_nmod_init(u, ectx->fqctx);
436     fq_nmod_init(v, ectx->fqctx);
437     nmod_poly_init(w, ctx->ffinfo->mod.n);
438 
439     FLINT_ASSERT(var > 0);
440     FLINT_ASSERT(T->bits == A->bits);
441     FLINT_ASSERT(F->bits == A->bits);
442 
443     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
444 
445     Flen = F->length;
446 
447     nmod_mpolyn_fit_length(T, FLINT_MAX(Flen, Alen), ctx);
448     Tcoeff = T->coeffs;
449     Texp = T->exps;
450     Ti = 0;
451 
452     Fi = Ai = vi = 0;
453     if (Ai < Alen)
454     {
455         vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
456     }
457     while (Fi < Flen || Ai < Alen)
458     {
459         if (Ti >= T->alloc)
460         {
461             nmod_mpolyn_fit_length(T, Ti + FLINT_MAX(Flen - Fi, Alen - Ai), ctx);
462             Tcoeff = T->coeffs;
463             Texp = T->exps;
464         }
465 
466         if (Fi < Flen && Ai < Alen && mpoly_monomial_equal_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift))
467         {
468             /* F term ok, A term ok */
469             nmod_poly_rem(u, Fcoeff + Fi, ectx->fqctx->modulus);
470             fq_nmod_sub(v, (Acoeff + Ai)->coeffs + vi, u, ectx->fqctx);
471             if (!fq_nmod_is_zero(v, ectx->fqctx))
472             {
473                 changed = 1;
474                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
475                 nmod_poly_mul(w, modulus, u);
476                 nmod_poly_add(Tcoeff + Ti, Fcoeff + Fi, w);
477             }
478             else
479             {
480                 nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi);
481             }
482             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
483 
484             Fi++;
485             do {
486                 vi--;
487             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ectx->fqctx));
488             if (vi < 0)
489             {
490                 Ai++;
491                 if (Ai < Alen)
492                 {
493                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
494                 }
495             }
496         }
497         else if (Fi < Flen && (Ai >= Alen || mpoly_monomial_gt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)))
498         {
499             /* F term ok, A term missing */
500             nmod_poly_rem(v, Fcoeff + Fi, ectx->fqctx->modulus);
501             if (!fq_nmod_is_zero(v, ectx->fqctx))
502             {
503                 changed = 1;
504                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
505                 nmod_poly_mul(w, u, modulus);
506                 nmod_poly_sub(Tcoeff + Ti, Fcoeff + Fi, w);
507             }
508             else
509             {
510                 nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi);
511             }
512             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
513 
514             Fi++;
515         }
516         else
517         {
518             FLINT_ASSERT(Ai < Alen && (Fi >= Flen || mpoly_monomial_lt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)));
519 
520             /* F term missing, A term ok */
521             changed = 1;
522             fq_nmod_mul(u, (Acoeff + Ai)->coeffs + vi, inv_m_eval, ectx->fqctx);
523             nmod_poly_mul(Tcoeff + Ti, modulus, u);
524             mpoly_monomial_set_extra(Texp + N*Ti, Aexp + N*Ai, N, offset, vi << shift);
525 
526             do {
527                 vi--;
528             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ectx->fqctx));
529             if (vi < 0)
530             {
531                 Ai++;
532                 if (Ai < Alen)
533                 {
534                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
535                 }
536             }
537         }
538 
539         FLINT_ASSERT(!nmod_poly_is_zero(Tcoeff + Ti));
540         lastdeg = FLINT_MAX(lastdeg, nmod_poly_degree(Tcoeff + Ti));
541         Ti++;
542     }
543     T->length = Ti;
544 
545     if (changed)
546     {
547         nmod_mpolyn_swap(T, F);
548     }
549 
550     fq_nmod_clear(inv_m_eval, ectx->fqctx);
551     fq_nmod_clear(u, ectx->fqctx);
552     fq_nmod_clear(v, ectx->fqctx);
553     nmod_poly_clear(w);
554 
555     FLINT_ASSERT(nmod_mpolyn_is_canonical(F, ctx));
556 
557     *lastdeg_ = lastdeg;
558 
559     return changed;
560 }
561 
562 
563 /******************** multivar - image in Fq *********************************/
564 
565 /*
566     A = B mod alpha(x_var)
567     B is in               Fp [x_0, ..., x_(var-2), x_(var-1)][x_var]
568     A is in (Fp/alpha(x_var))[x_0, ..., x_(var-2), x_(var-1)]
569     alpha is ectx->modulus
570 */
nmod_mpolyn_interp_reduce_lg_mpoly(fq_nmod_mpoly_t A,nmod_mpolyn_t B,const fq_nmod_mpoly_ctx_t ffctx,const nmod_mpoly_ctx_t ctx)571 void nmod_mpolyn_interp_reduce_lg_mpoly(
572     fq_nmod_mpoly_t A,
573     nmod_mpolyn_t B,
574     const fq_nmod_mpoly_ctx_t ffctx,
575     const nmod_mpoly_ctx_t ctx)
576 {
577     slong i, k, N;
578 
579     FLINT_ASSERT(A->bits == B->bits);
580     FLINT_ASSERT(B->bits <= FLINT_BITS);
581     FLINT_ASSERT(ctx->minfo->ord == ORD_LEX);
582     FLINT_ASSERT(ffctx->minfo->ord == ORD_LEX);
583     FLINT_ASSERT(ffctx->minfo->nvars == ctx->minfo->nvars);
584 
585     N = mpoly_words_per_exp(B->bits, ctx->minfo);
586 
587     k = 0;
588     fq_nmod_mpoly_fit_length(A, B->length, ffctx);
589     for (i = 0; i < B->length; i++)
590     {
591         mpoly_monomial_set(A->exps + N*k, B->exps + N*i, N);
592         nmod_poly_rem(A->coeffs + k, B->coeffs + i, ffctx->fqctx->modulus);
593         k += !fq_nmod_is_zero(A->coeffs + k, ffctx->fqctx);
594     }
595     A->length = k;
596 }
597 
598 /*
599     A = B mod alpha(x_var)
600     B is in               Fp [X][x_0, ..., x_(var-2), x_(var-1)][x_var]
601     A is in (Fp/alpha(x_var))[X][x_0, ..., x_(var-2), x_(var-1)]
602     alpha is ectx->modulus
603 */
nmod_mpolyun_interp_reduce_lg_mpolyu(fq_nmod_mpolyu_t A,nmod_mpolyun_t B,const fq_nmod_mpoly_ctx_t ffctx,const nmod_mpoly_ctx_t ctx)604 void nmod_mpolyun_interp_reduce_lg_mpolyu(
605     fq_nmod_mpolyu_t A,
606     nmod_mpolyun_t B,
607     const fq_nmod_mpoly_ctx_t ffctx,
608     const nmod_mpoly_ctx_t ctx)
609 {
610     slong i, k, Blen;
611     fq_nmod_mpoly_struct * Acoeff;
612     nmod_mpolyn_struct * Bcoeff;
613     ulong * Aexp, * Bexp;
614 
615     Blen = B->length;
616     fq_nmod_mpolyu_fit_length(A, Blen, ffctx);
617     Acoeff = A->coeffs;
618     Bcoeff = B->coeffs;
619     Aexp = A->exps;
620     Bexp = B->exps;
621 
622     k = 0;
623     for (i = 0; i < Blen; i++)
624     {
625         nmod_mpolyn_interp_reduce_lg_mpoly(Acoeff + k, Bcoeff + i, ffctx, ctx);
626         Aexp[k] = Bexp[i];
627         k += !fq_nmod_mpoly_is_zero(Acoeff + k, ffctx);
628     }
629     A->length = k;
630 }
631 
632 
633 /*
634     A = Ap using lowest degree
635     A  is in               Fp [x_0, ..., x_(var-2), x_(var-1)][x_var]
636     Ap is in (Fp/alpha(x_var))[x_0, ..., x_(var-2), x_(var-1)]
637     alpha is ectx->modulus
638 */
nmod_mpolyn_interp_lift_lg_mpoly(nmod_mpolyn_t A,const nmod_mpoly_ctx_t ctx,fq_nmod_mpoly_t Ap,const fq_nmod_mpoly_ctx_t ctxp)639 void nmod_mpolyn_interp_lift_lg_mpoly(
640     nmod_mpolyn_t A,
641     const nmod_mpoly_ctx_t ctx,
642     fq_nmod_mpoly_t Ap,
643     const fq_nmod_mpoly_ctx_t ctxp)
644 {
645     slong i, N;
646 
647     FLINT_ASSERT(Ap->bits == A->bits);
648     nmod_mpolyn_fit_length(A, Ap->length, ctx);
649     N = mpoly_words_per_exp(A->bits, ctx->minfo);
650     for (i = 0; i < Ap->length; i++)
651     {
652         mpoly_monomial_set(A->exps + N*i, Ap->exps + N*i, N);
653         nmod_poly_set(A->coeffs + i, Ap->coeffs + i);
654     }
655     A->length = Ap->length;
656 }
657 
658 /*
659     A = Ap using lowest degree
660     A  is in               Fp [X][x_0, ..., x_(var-2), x_(var-1)][x_var]
661     Ap is in (Fp/alpha(x_var))[X][x_0, ..., x_(var-2), x_(var-1)]
662     alpha is ectx->modulus
663 */
nmod_mpolyun_interp_lift_lg_mpolyu(nmod_mpolyun_t A,const nmod_mpoly_ctx_t ctx,fq_nmod_mpolyu_t Ap,const fq_nmod_mpoly_ctx_t ctxp)664 void nmod_mpolyun_interp_lift_lg_mpolyu(
665     nmod_mpolyun_t A,
666     const nmod_mpoly_ctx_t ctx,
667     fq_nmod_mpolyu_t Ap,
668     const fq_nmod_mpoly_ctx_t ctxp)
669 {
670     slong i;
671 
672     FLINT_ASSERT(Ap->bits == A->bits);
673     nmod_mpolyun_fit_length(A, Ap->length, ctx);
674     for (i = 0; i < Ap->length; i++)
675     {
676         A->exps[i] = Ap->exps[i];
677         nmod_mpolyn_interp_lift_lg_mpoly(A->coeffs + i, ctx, Ap->coeffs + i, ctxp);
678     }
679     A->length = Ap->length;
680 }
681 
682 /*
683     F = F + modulus*((A - F(alpha))/(modulus(alpha)))
684     no assumptions about matching monomials
685 */
nmod_mpolyn_interp_crt_lg_mpoly(slong * lastdeg,nmod_mpolyn_t F,nmod_mpolyn_t T,nmod_poly_t m,const nmod_mpoly_ctx_t ctx,fq_nmod_mpoly_t A,fq_nmod_t inv_m_eval,const fq_nmod_mpoly_ctx_t ffctx)686 int nmod_mpolyn_interp_crt_lg_mpoly(
687     slong * lastdeg,
688     nmod_mpolyn_t F,
689     nmod_mpolyn_t T,
690     nmod_poly_t m,
691     const nmod_mpoly_ctx_t ctx,
692     fq_nmod_mpoly_t A,
693     fq_nmod_t inv_m_eval,
694     const fq_nmod_mpoly_ctx_t ffctx)
695 {
696     int changed = 0;
697     slong i, j, k;
698     slong N;
699     fq_nmod_t u, v;
700     nmod_poly_t w;
701     flint_bitcnt_t bits = A->bits;
702     slong Flen = F->length, Alen = A->length;
703     ulong * Fexp = F->exps, * Aexp = A->exps;
704     ulong * Texp;
705     fq_nmod_struct * Acoeff = A->coeffs;
706     nmod_poly_struct * Fcoeff = F->coeffs;
707     nmod_poly_struct * Tcoeff;
708 
709     FLINT_ASSERT(F->bits == bits);
710     FLINT_ASSERT(ctx->minfo->ord == ORD_LEX);
711 
712     fq_nmod_init(u, ffctx->fqctx);
713     fq_nmod_init(v, ffctx->fqctx);
714     nmod_poly_init(w, ffctx->fqctx->modulus->mod.n);
715 
716     nmod_mpolyn_fit_length(T, Flen + Alen, ctx);
717     Texp = T->exps;
718     Tcoeff = T->coeffs;
719 
720     N = mpoly_words_per_exp(bits, ctx->minfo);
721 
722     i = j = k = 0;
723     while (i < Flen || j < Alen)
724     {
725         if (i < Flen && (j >= Alen
726                         || mpoly_monomial_gt_nomask(Fexp + N*i, Aexp + N*j, N)))
727         {
728             FLINT_ASSERT(!nmod_poly_is_zero(Fcoeff + i));
729             FLINT_ASSERT(nmod_poly_degree(Fcoeff + i) < nmod_poly_degree(m));
730 
731             /* F term ok, A term missing */
732             nmod_poly_rem(v, Fcoeff + i, ffctx->fqctx->modulus);
733             if (!fq_nmod_is_zero(v, ffctx->fqctx))
734             {
735                 changed = 1;
736                 fq_nmod_mul(u, v, inv_m_eval, ffctx->fqctx);
737                 nmod_poly_mul(w, u, m);
738                 nmod_poly_sub(Tcoeff + k, Fcoeff + i, w);
739             } else {
740                 nmod_poly_set(Tcoeff + k, Fcoeff + i);
741             }
742             *lastdeg = FLINT_MAX(*lastdeg, nmod_poly_degree(Tcoeff + k));
743 
744             mpoly_monomial_set(Texp + N*k, Fexp + N*i, N);
745             FLINT_ASSERT(!nmod_poly_is_zero(Tcoeff + k));
746             k++;
747             i++;
748         }
749         else if (j < Alen && (i >= Flen
750                         || mpoly_monomial_lt_nomask(Fexp + N*i, Aexp + N*j, N)))
751         {
752             /* F term missing, A term ok */
753             if (!fq_nmod_is_zero(Acoeff + j, ffctx->fqctx))
754             {
755                 changed = 1;
756                 fq_nmod_mul(u, Acoeff + j, inv_m_eval, ffctx->fqctx);
757                 nmod_poly_mul(Tcoeff + k, m, u);
758                 *lastdeg = FLINT_MAX(*lastdeg, nmod_poly_degree(Tcoeff + k));
759                 mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
760                 k++;
761             }
762             j++;
763         }
764         else if (i < Flen && j < Alen
765                              && mpoly_monomial_equal(Fexp + N*i, Aexp + N*j, N))
766         {
767             FLINT_ASSERT(!nmod_poly_is_zero(Fcoeff + i));
768             FLINT_ASSERT(nmod_poly_degree(Fcoeff + i) < nmod_poly_degree(m));
769 
770             /* F term ok, A term ok */
771             nmod_poly_rem(u, Fcoeff + i, ffctx->fqctx->modulus);
772             fq_nmod_sub(v, Acoeff + j, u, ffctx->fqctx);
773             if (!fq_nmod_is_zero(v, ffctx->fqctx))
774             {
775                 changed = 1;
776                 fq_nmod_mul(u, v, inv_m_eval, ffctx->fqctx);
777                 nmod_poly_mul(w, m, u);
778                 nmod_poly_add(Tcoeff + k, Fcoeff + i, w);
779             } else {
780                 nmod_poly_set(Tcoeff + k, Fcoeff + i);
781             }
782             *lastdeg = FLINT_MAX(*lastdeg, nmod_poly_degree(Tcoeff + k));
783             mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
784             FLINT_ASSERT(!nmod_poly_is_zero(Tcoeff + k));
785             k++;
786             i++;
787             j++;
788         }
789         else
790         {
791             FLINT_ASSERT(0);
792         }
793     }
794 
795     nmod_mpolyn_set_length(T, k, ctx);
796 
797     if (changed)
798     {
799         nmod_mpolyn_swap(T, F);
800     }
801 
802     fq_nmod_clear(u, ffctx->fqctx);
803     fq_nmod_clear(v, ffctx->fqctx);
804     nmod_poly_clear(w);
805 
806     return changed;
807 }
808 
809 /*
810     F = F + modulus*((A - F(alpha))/(modulus(alpha)))
811     no assumptions about matching monomials
812 */
nmod_mpolyun_interp_crt_lg_mpolyu(slong * lastdeg,nmod_mpolyun_t F,nmod_mpolyun_t T,nmod_poly_t m,const nmod_mpoly_ctx_t ctx,fq_nmod_mpolyu_t A,const fq_nmod_mpoly_ctx_t ffctx)813 int nmod_mpolyun_interp_crt_lg_mpolyu(
814     slong * lastdeg,
815     nmod_mpolyun_t F,
816     nmod_mpolyun_t T,
817     nmod_poly_t m,
818     const nmod_mpoly_ctx_t ctx,
819     fq_nmod_mpolyu_t A,
820     const fq_nmod_mpoly_ctx_t ffctx)
821 {
822     int changed = 0;
823     slong i, j, k;
824     ulong * Texp;
825     ulong * Fexp;
826     ulong * Aexp;
827     slong Flen;
828     slong Alen;
829     nmod_mpolyn_t S;
830     nmod_mpolyn_struct * Tcoeff;
831     nmod_mpolyn_struct * Fcoeff;
832     fq_nmod_mpoly_struct  * Acoeff;
833     fq_nmod_mpoly_t zero;
834     fq_nmod_t inv_m_eval;
835 
836     *lastdeg = -WORD(1);
837 
838     FLINT_ASSERT(F->bits == T->bits);
839     FLINT_ASSERT(T->bits == A->bits);
840 
841     nmod_mpolyn_init(S, F->bits, ctx);
842 
843     Flen = F->length;
844     Alen = A->length;
845     nmod_mpolyun_fit_length(T, Flen + Alen, ctx);
846 
847     Tcoeff = T->coeffs;
848     Fcoeff = F->coeffs;
849     Acoeff = A->coeffs;
850     Texp = T->exps;
851     Fexp = F->exps;
852     Aexp = A->exps;
853 
854     fq_nmod_mpoly_init(zero, ffctx);
855     fq_nmod_mpoly_fit_bits(zero, A->bits, ffctx);
856     zero->bits = A->bits;
857 
858     fq_nmod_init(inv_m_eval, ffctx->fqctx);
859     nmod_poly_rem(inv_m_eval, m, ffctx->fqctx->modulus);
860     fq_nmod_inv(inv_m_eval, inv_m_eval, ffctx->fqctx);
861 
862     i = j = k = 0;
863     while (i < Flen || j < Alen)
864     {
865         if (i < Flen && (j >= Alen || Fexp[i] > Aexp[j]))
866         {
867             /* F term ok, A term missing */
868             nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
869             changed |= nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
870                                            S, m, ctx, zero, inv_m_eval, ffctx);
871             Texp[k] = Fexp[i];
872             k++;
873             i++;
874         }
875         else if (j < Alen && (i >= Flen || Aexp[j] > Fexp[i]))
876         {
877             /* F term missing, A term ok */
878             nmod_mpolyn_zero(Tcoeff + k, ctx);
879             changed |= nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
880                                      S, m, ctx, Acoeff + j, inv_m_eval, ffctx);
881             Texp[k] = Aexp[j];
882             k++;
883             j++;
884         }
885         else if (i < Flen && j < Alen && (Fexp[i] == Aexp[j]))
886         {
887             /* F term ok, A term ok */
888             nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
889             changed |= nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
890                                      S, m, ctx, Acoeff + j, inv_m_eval, ffctx);
891             Texp[k] = Aexp[j];
892             FLINT_ASSERT(!nmod_mpolyn_is_zero(Tcoeff + k, ctx));
893             k++;
894             i++;
895             j++;
896         }
897         else
898         {
899             FLINT_ASSERT(0);
900         }
901     }
902     T->length = k;
903 
904     if (changed)
905     {
906         nmod_mpolyun_swap(T, F);
907     }
908 
909     fq_nmod_clear(inv_m_eval, ffctx->fqctx);
910 
911     nmod_mpolyn_clear(S, ctx);
912     fq_nmod_mpoly_clear(zero, ffctx);
913     return changed;
914 }
915 
916 
917 /*
918     Update H so that it does not change mod m, and is now A mod p
919     It is asserted that the monomials in H and A match
920 */
nmod_mpolyn_CRT_fq_nmod_mpoly(slong * lastdeg,nmod_mpolyn_t H,const nmod_mpoly_ctx_t ctx,nmod_poly_t m,fq_nmod_t inv_m_eval,fq_nmod_mpoly_t A,const fq_nmod_mpoly_ctx_t ctxp)921 int nmod_mpolyn_CRT_fq_nmod_mpoly(
922     slong * lastdeg,
923     nmod_mpolyn_t H,
924     const nmod_mpoly_ctx_t ctx,
925     nmod_poly_t m,
926     fq_nmod_t inv_m_eval,
927     fq_nmod_mpoly_t A,
928     const fq_nmod_mpoly_ctx_t ctxp)
929 {
930     slong i;
931 #if WANT_ASSERT
932     slong N;
933 #endif
934     int changed = 0;
935     fq_nmod_t u, v;
936     nmod_poly_t w;
937 
938     fq_nmod_init(u, ctxp->fqctx);
939     fq_nmod_init(v, ctxp->fqctx);
940     nmod_poly_init(w, ctxp->fqctx->modulus->mod.n);
941 
942     FLINT_ASSERT(H->length == A->length);
943     FLINT_ASSERT(H->bits == A->bits);
944 #if WANT_ASSERT
945     N = mpoly_words_per_exp(A->bits, ctx->minfo);
946 #endif
947     for (i = 0; i < A->length; i++)
948     {
949         FLINT_ASSERT(mpoly_monomial_equal(H->exps + N*i, A->exps + N*i, N));
950         nmod_poly_rem(u, H->coeffs + i, ctxp->fqctx->modulus);
951         fq_nmod_sub(v, A->coeffs + i, u, ctxp->fqctx);
952         if (!fq_nmod_is_zero(v, ctxp->fqctx))
953         {
954             changed = 1;
955             fq_nmod_mul(u, v, inv_m_eval, ctxp->fqctx);
956             nmod_poly_mul(w, u, m);
957             nmod_poly_add(H->coeffs + i, H->coeffs + i, w);
958         }
959 
960         *lastdeg = FLINT_MAX(*lastdeg, nmod_poly_degree(H->coeffs + i));
961 
962         FLINT_ASSERT(nmod_poly_degree(H->coeffs + i) < nmod_poly_degree(m)
963                                                      + nmod_poly_degree(ctxp->fqctx->modulus));
964     }
965 
966     fq_nmod_clear(u, ctxp->fqctx);
967     fq_nmod_clear(v, ctxp->fqctx);
968     nmod_poly_clear(w);
969 
970     return changed;
971 }
972 
nmod_mpolyun_interp_mcrt_lg_mpolyu(slong * lastdeg,nmod_mpolyun_t H,const nmod_mpoly_ctx_t ctx,nmod_poly_t m,fq_nmod_mpolyu_t A,const fq_nmod_mpoly_ctx_t ctxp)973 int nmod_mpolyun_interp_mcrt_lg_mpolyu(
974     slong * lastdeg,
975     nmod_mpolyun_t H,
976     const nmod_mpoly_ctx_t ctx,
977     nmod_poly_t m,
978     fq_nmod_mpolyu_t A,
979     const fq_nmod_mpoly_ctx_t ctxp)
980 {
981     slong i;
982     int changed = 0;
983     fq_nmod_t inv_m_eval;
984 
985     *lastdeg = -WORD(1);
986 
987     fq_nmod_init(inv_m_eval, ctxp->fqctx);
988     nmod_poly_rem(inv_m_eval, m, ctxp->fqctx->modulus);
989     fq_nmod_inv(inv_m_eval, inv_m_eval, ctxp->fqctx);
990 
991     FLINT_ASSERT(H->bits == A->bits);
992     FLINT_ASSERT(H->length == A->length);
993     for (i = 0; i < A->length; i++)
994     {
995         FLINT_ASSERT(H->exps[i] == A->exps[i]);
996         changed |= nmod_mpolyn_CRT_fq_nmod_mpoly(lastdeg, H->coeffs + i, ctx, m,
997                                               inv_m_eval, A->coeffs + i, ctxp);
998     }
999     H->length = A->length;
1000     fq_nmod_clear(inv_m_eval, ctxp->fqctx);
1001     return changed;
1002 }
1003 
1004 /*****************************************************************************/
1005 
1006 /*
1007     E = A(v = alpha)
1008     A is in Fq[X][v]
1009     E is in Fq[X]
1010 */
1011 
fq_nmod_mpolyn_interp_reduce_sm_poly(fq_nmod_poly_t E,const fq_nmod_mpolyn_t A,const fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)1012 void fq_nmod_mpolyn_interp_reduce_sm_poly(
1013     fq_nmod_poly_t E,
1014     const fq_nmod_mpolyn_t A,
1015     const fq_nmod_t alpha,
1016     const fq_nmod_mpoly_ctx_t ctx)
1017 {
1018     fq_nmod_t v;
1019     slong Ai, Alen, k;
1020     fq_nmod_poly_struct * Acoeff;
1021     ulong * Aexp;
1022     slong N, off, shift;
1023 
1024     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1025     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
1026 
1027     fq_nmod_init(v, ctx->fqctx);
1028 
1029     Acoeff = A->coeffs;
1030     Aexp = A->exps;
1031     Alen = A->length;
1032     Ai = 0;
1033     fq_nmod_poly_zero(E, ctx->fqctx);
1034     for (Ai = 0; Ai < Alen; Ai++)
1035     {
1036         fq_nmod_poly_evaluate_fq_nmod(v, Acoeff + Ai, alpha, ctx->fqctx);
1037         k = (Aexp + N*Ai)[off] >> shift;
1038         fq_nmod_poly_set_coeff(E, k, v, ctx->fqctx);
1039     }
1040 
1041     fq_nmod_clear(v, ctx->fqctx);
1042 }
1043 
1044 /*
1045     A = B
1046     A is in Fq[X][v]  (no v appears)
1047     B is in Fq[X]
1048 */
1049 
fq_nmod_mpolyn_interp_lift_sm_poly(fq_nmod_mpolyn_t A,const fq_nmod_poly_t B,const fq_nmod_mpoly_ctx_t ctx)1050 void fq_nmod_mpolyn_interp_lift_sm_poly(
1051     fq_nmod_mpolyn_t A,
1052     const fq_nmod_poly_t B,
1053     const fq_nmod_mpoly_ctx_t ctx)
1054 {
1055     slong Bi;
1056     slong Blen = fq_nmod_poly_length(B, ctx->fqctx);
1057     fq_nmod_struct * Bcoeff = B->coeffs;
1058     fq_nmod_poly_struct * Acoeff;
1059     ulong * Aexp;
1060     slong Ai;
1061     slong N, off, shift;
1062 
1063     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1064     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
1065 
1066     fq_nmod_mpolyn_fit_length(A, Blen, ctx);
1067     Acoeff = A->coeffs;
1068     Aexp = A->exps;
1069 
1070     Ai = 0;
1071     for (Bi = Blen - 1; Bi >= 0; Bi--)
1072     {
1073         if (!fq_nmod_is_zero(Bcoeff + Bi, ctx->fqctx))
1074         {
1075             FLINT_ASSERT(Ai < A->alloc);
1076 
1077             fq_nmod_poly_set_fq_nmod(Acoeff + Ai, Bcoeff + Bi, ctx->fqctx);
1078             mpoly_monomial_zero(Aexp + N*Ai, N);
1079             (Aexp + N*Ai)[off] = Bi << shift;
1080             Ai++;
1081         }
1082     }
1083     A->length = Ai;
1084 }
1085 
1086 /*
1087     F = F + modulus*(A - F(v = alpha))
1088     no assumptions about matching monomials
1089     F is in Fq[X][v]
1090     A is in Fq[X]
1091     it is expected that modulus(alpha) == 1
1092 */
1093 
fq_nmod_mpolyn_interp_crt_sm_poly(slong * lastdeg_,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,const fq_nmod_poly_t A,const fq_nmod_poly_t modulus,const fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)1094 int fq_nmod_mpolyn_interp_crt_sm_poly(
1095     slong * lastdeg_,
1096     fq_nmod_mpolyn_t F,
1097     fq_nmod_mpolyn_t T,
1098     const fq_nmod_poly_t A,
1099     const fq_nmod_poly_t modulus,
1100     const fq_nmod_t alpha,
1101     const fq_nmod_mpoly_ctx_t ctx)
1102 {
1103     int changed = 0;
1104     slong lastdeg = -WORD(1);
1105     fq_nmod_t u, v;
1106     slong Fi, Ti, Ai;
1107     fq_nmod_struct * Acoeff = A->coeffs;
1108     slong Flen = F->length;
1109     fq_nmod_poly_struct * Fcoeffs = F->coeffs;
1110     ulong * Fexps = F->exps;
1111     fq_nmod_poly_struct * Tcoeffs;
1112     ulong * Texps;
1113     fq_nmod_poly_t tp;
1114     slong N, off, shift;
1115 
1116     N = mpoly_words_per_exp_sp(F->bits, ctx->minfo);
1117     mpoly_gen_offset_shift_sp(&off, &shift, 0, F->bits, ctx->minfo);
1118 
1119     Fi = 0;
1120     Ai = fq_nmod_poly_degree(A, ctx->fqctx);
1121 
1122     fq_nmod_init(u, ctx->fqctx);
1123     fq_nmod_init(v, ctx->fqctx);
1124     fq_nmod_poly_init(tp, ctx->fqctx);
1125 
1126     fq_nmod_mpolyn_fit_length(T, Flen + Ai + 1, ctx);
1127     Tcoeffs = T->coeffs;
1128     Texps = T->exps;
1129     Ti = 0;
1130 
1131     while (Fi < Flen || Ai >= 0)
1132     {
1133         FLINT_ASSERT(Ti < T->alloc);
1134 
1135         if (Fi < Flen)
1136         {
1137             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeffs + Fi, ctx->fqctx));
1138             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeffs + Fi, ctx->fqctx) < fq_nmod_poly_degree(modulus, ctx->fqctx));
1139         }
1140 
1141         if (Ai >= 0)
1142         {
1143             FLINT_ASSERT(!fq_nmod_is_zero(Acoeff + Ai, ctx->fqctx));
1144         }
1145 
1146         mpoly_monomial_zero(Texps + N*Ti, N);
1147 
1148         if (Fi < Flen && Ai >= 0 && ((Fexps + N*Fi)[off]>>shift) == Ai)
1149         {
1150             /* F term ok, A term ok */
1151             fq_nmod_poly_evaluate_fq_nmod(u, Fcoeffs + Fi, alpha, ctx->fqctx);
1152             fq_nmod_sub(v, Acoeff + Ai, u, ctx->fqctx);
1153             if (!fq_nmod_is_zero(v, ctx->fqctx))
1154             {
1155                 changed = 1;
1156                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
1157                 fq_nmod_poly_add(Tcoeffs + Ti, Fcoeffs + Fi, tp, ctx->fqctx);
1158             }
1159             else
1160             {
1161                 fq_nmod_poly_set(Tcoeffs + Ti, Fcoeffs + Fi, ctx->fqctx);
1162             }
1163 
1164             (Texps + N*Ti)[off] = Ai << shift;
1165             Fi++;
1166             do {
1167                 Ai--;
1168             } while (Ai >= 0 && fq_nmod_is_zero(Acoeff + Ai, ctx->fqctx));
1169         }
1170         else if (Fi < Flen && (Ai < 0 || ((Fexps + N*Fi)[off]>>shift) > Ai))
1171         {
1172             /* F term ok, A term missing */
1173             fq_nmod_poly_evaluate_fq_nmod(v, Fcoeffs + Fi, alpha, ctx->fqctx);
1174             if (!fq_nmod_is_zero(v, ctx->fqctx))
1175             {
1176                 changed = 1;
1177                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
1178                 fq_nmod_poly_sub(Tcoeffs + Ti, Fcoeffs + Fi, tp, ctx->fqctx);
1179             }
1180             else
1181             {
1182                 fq_nmod_poly_set(Tcoeffs + Ti, Fcoeffs + Fi, ctx->fqctx);
1183             }
1184 
1185             (Texps + N*Ti)[off] = (Fexps + N*Fi)[off];
1186             Fi++;
1187         }
1188         else if (Ai >= 0 && (Fi >= Flen || ((Fexps + N*Fi)[off]>>shift) < Ai))
1189         {
1190             /* F term missing, A term ok */
1191             changed = 1;
1192             fq_nmod_poly_scalar_mul_fq_nmod(Tcoeffs + Ti, modulus, Acoeff + Ai, ctx->fqctx);
1193 
1194             (Texps + N*Ti)[off] = Ai << shift;
1195             do {
1196                 Ai--;
1197             } while (Ai >= 0 && fq_nmod_is_zero(Acoeff + Ai, ctx->fqctx));
1198         }
1199         else
1200         {
1201             FLINT_ASSERT(0);
1202         }
1203 
1204         FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeffs + Ti, ctx->fqctx));
1205         lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Tcoeffs + Ti, ctx->fqctx));
1206 
1207         Ti++;
1208     }
1209     T->length = Ti;
1210 
1211     fq_nmod_clear(u, ctx->fqctx);
1212     fq_nmod_clear(v, ctx->fqctx);
1213     fq_nmod_poly_clear(tp, ctx->fqctx);
1214 
1215     if (changed)
1216     {
1217         fq_nmod_mpolyn_swap(T, F);
1218     }
1219 
1220     *lastdeg_ = lastdeg;
1221     return changed;
1222 }
1223 
1224 
1225 /****************************************************************************/
1226 
1227 /*
1228     E = A(x_var = alpha)
1229     A is in Fq[x_0, ..., x_(var-2), x_(var-1)][x_var]
1230     E is in Fq[x_0, ..., x_(var-2)][x_(var-1)]
1231 */
fq_nmod_mpolyn_interp_reduce_sm_mpolyn(fq_nmod_mpolyn_t E,fq_nmod_mpolyn_t A,slong var,fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)1232 void fq_nmod_mpolyn_interp_reduce_sm_mpolyn(
1233     fq_nmod_mpolyn_t E,
1234     fq_nmod_mpolyn_t A,
1235     slong var,
1236     fq_nmod_t alpha,
1237     const fq_nmod_mpoly_ctx_t ctx)
1238 {
1239     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1240     slong offset, shift, k;
1241     ulong mask;
1242     fq_nmod_t v;
1243     fq_nmod_poly_struct * Acoeff = A->coeffs;
1244     ulong * Aexp = A->exps;
1245     slong Alen = A->length;
1246     slong Ai;
1247     fq_nmod_poly_struct * Ecoeff;
1248     ulong * Eexp;
1249     slong Ei;
1250 
1251     fq_nmod_init(v, ctx->fqctx);
1252 
1253     FLINT_ASSERT(var > 0);
1254     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1255     mask = (-UWORD(1)) >> (FLINT_BITS - A->bits);
1256 
1257     Ecoeff = E->coeffs;
1258     Eexp = E->exps;
1259     Ei = 0;
1260     for (Ai = 0; Ai < Alen; Ai++)
1261     {
1262         fq_nmod_poly_evaluate_fq_nmod(v, Acoeff + Ai, alpha, ctx->fqctx);
1263         k = ((Aexp + N*Ai)[offset] >> shift) & mask;
1264         if (fq_nmod_is_zero(v, ctx->fqctx))
1265         {
1266             continue;
1267         }
1268 
1269         if (Ei > 0 && mpoly_monomial_equal_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)))
1270         {
1271             /* append to previous */
1272             fq_nmod_poly_set_coeff(Ecoeff + Ei - 1, k, v, ctx->fqctx);
1273         }
1274         else
1275         {
1276             FLINT_ASSERT(Ei == 0 || mpoly_monomial_gt_nomask_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)));
1277 
1278             /* create new */
1279             if (Ei >= E->alloc)
1280             {
1281                 fq_nmod_mpolyn_fit_length(E, Ei + 1, ctx);
1282                 Ecoeff = E->coeffs;
1283                 Eexp = E->exps;
1284             }
1285             mpoly_monomial_set_extra(Eexp + N*Ei, Aexp + N*Ai, N, offset, -(k << shift));
1286             fq_nmod_poly_zero(Ecoeff + Ei, ctx->fqctx);
1287             fq_nmod_poly_set_coeff(Ecoeff + Ei, k, v, ctx->fqctx);
1288             Ei++;
1289         }
1290     }
1291     E->length = Ei;
1292 
1293     fq_nmod_clear(v, ctx->fqctx);
1294 }
1295 
1296 
1297 /*
1298     A = B
1299     A is in Fq[x_0, ..., x_(var-1), x_(var-1)][x_var]
1300     B is in Fq[x_0, ..., x_(var-2)][x_(var-1)]
1301 */
fq_nmod_mpolyn_interp_lift_sm_mpolyn(fq_nmod_mpolyn_t A,fq_nmod_mpolyn_t B,slong var,const fq_nmod_mpoly_ctx_t ctx)1302 void fq_nmod_mpolyn_interp_lift_sm_mpolyn(
1303     fq_nmod_mpolyn_t A,
1304     fq_nmod_mpolyn_t B,
1305     slong var,
1306     const fq_nmod_mpoly_ctx_t ctx)
1307 {
1308     slong N = mpoly_words_per_exp_sp(B->bits, ctx->minfo);
1309     slong offset, shift;
1310     slong vi;
1311 
1312     fq_nmod_poly_struct * Bcoeff = B->coeffs;
1313     ulong * Bexp = B->exps;
1314     slong Blen = B->length;
1315     slong Bi;
1316 
1317     fq_nmod_poly_struct * Acoeff;
1318     ulong * Aexp;
1319     slong Ai;
1320 
1321     fq_nmod_mpolyn_fit_length(A, Blen, ctx);
1322     Acoeff = A->coeffs;
1323     Aexp = A->exps;
1324 
1325     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1326 
1327     Ai = 0;
1328     for (Bi = 0; Bi < Blen; Bi++)
1329     {
1330         if (Ai + (Bcoeff + Bi)->length >= A->alloc)
1331         {
1332             fq_nmod_mpolyn_fit_length(A, Ai + (Bcoeff + Bi)->length, ctx);
1333             Acoeff = A->coeffs;
1334             Aexp = A->exps;
1335         }
1336         for (vi = (Bcoeff + Bi)->length - 1; vi >= 0; vi--)
1337         {
1338             if (!fq_nmod_is_zero((Bcoeff + Bi)->coeffs + vi, ctx->fqctx))
1339             {
1340                 mpoly_monomial_set_extra(Aexp + N*Ai, Bexp + N*Bi, N, offset, vi << shift);
1341                 fq_nmod_poly_zero(Acoeff + Ai, ctx->fqctx);
1342                 fq_nmod_poly_set_coeff(Acoeff + Ai, 0, (Bcoeff + Bi)->coeffs + vi, ctx->fqctx);
1343                 Ai++;
1344             }
1345         }
1346     }
1347     A->length = Ai;
1348 }
1349 
1350 
1351 /*
1352     T = F + modulus*(A - F(x_var = alpha))
1353     no assumptions about matching monomials
1354     F is in Fq[x_0, ..., x_(var-1), x_(var-1)][x_var]
1355     A is in Fq[x_0, ..., x_(var-2)][x_(var-1)]
1356     in order to fxn correctly, modulus(alpha) should be 1
1357 */
fq_nmod_mpolyn_interp_crt_sm_mpolyn(slong * lastdeg_,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,fq_nmod_mpolyn_t A,slong var,fq_nmod_poly_t modulus,const fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)1358 int fq_nmod_mpolyn_interp_crt_sm_mpolyn(
1359     slong * lastdeg_,
1360     fq_nmod_mpolyn_t F,
1361     fq_nmod_mpolyn_t T,
1362     fq_nmod_mpolyn_t A,
1363     slong var,
1364     fq_nmod_poly_t modulus,
1365     const fq_nmod_t alpha,
1366     const fq_nmod_mpoly_ctx_t ctx)
1367 {
1368     int changed = 0;
1369     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1370     slong lastdeg = -WORD(1);
1371     slong offset, shift;
1372     slong vi;
1373     fq_nmod_t v;
1374     fq_nmod_poly_t tp;
1375 
1376     fq_nmod_poly_struct * Tcoeff;
1377     ulong * Texp;
1378     slong Ti;
1379 
1380     fq_nmod_poly_struct * Acoeff = A->coeffs;
1381     slong Alen = A->length;
1382     ulong * Aexp = A->exps;
1383     slong Ai;
1384 
1385     fq_nmod_poly_struct * Fcoeff = F->coeffs;
1386     slong Flen = F->length;
1387     ulong * Fexp = F->exps;
1388     slong Fi;
1389 
1390     fq_nmod_poly_init(tp, ctx->fqctx);
1391     fq_nmod_init(v, ctx->fqctx);
1392 
1393     FLINT_ASSERT(var > 0);
1394     FLINT_ASSERT(T->bits == A->bits);
1395     FLINT_ASSERT(F->bits == A->bits);
1396 
1397     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1398 
1399     Flen = F->length;
1400 
1401     fq_nmod_mpolyn_fit_length(T, FLINT_MAX(Flen, Alen), ctx);
1402     Tcoeff = T->coeffs;
1403     Texp = T->exps;
1404     Ti = 0;
1405 
1406     Fi = Ai = vi = 0;
1407     if (Ai < Alen)
1408     {
1409         vi = fq_nmod_poly_degree(A->coeffs + Ai, ctx->fqctx);
1410     }
1411     while (Fi < Flen || Ai < Alen)
1412     {
1413         if (Ti >= T->alloc)
1414         {
1415             fq_nmod_mpolyn_fit_length(T, Ti + FLINT_MAX(Flen - Fi, Alen - Ai), ctx);
1416             Tcoeff = T->coeffs;
1417             Texp = T->exps;
1418         }
1419 
1420         if (Ai < Alen)
1421         {
1422             FLINT_ASSERT(!fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ctx->fqctx));
1423         }
1424 
1425         if (Fi < Flen && Ai < Alen && mpoly_monomial_equal_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift))
1426         {
1427             /* F term ok, A term ok */
1428             fq_nmod_poly_evaluate_fq_nmod(v, Fcoeff + Fi, alpha, ctx->fqctx);
1429             fq_nmod_sub(v, (Acoeff + Ai)->coeffs + vi, v, ctx->fqctx);
1430             if (!fq_nmod_is_zero(v, ctx->fqctx))
1431             {
1432                 changed = 1;
1433                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
1434                 fq_nmod_poly_add(Tcoeff + Ti, Fcoeff + Fi, tp, ctx->fqctx);
1435             }
1436             else
1437             {
1438                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
1439             }
1440             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
1441 
1442             Fi++;
1443             do {
1444                 vi--;
1445             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ctx->fqctx));
1446             if (vi < 0)
1447             {
1448                 Ai++;
1449                 if (Ai < Alen)
1450                 {
1451                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ctx->fqctx);
1452                 }
1453             }
1454         }
1455         else if (Fi < Flen && (Ai >= Alen || mpoly_monomial_gt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)))
1456         {
1457             /* F term ok, A term missing */
1458             fq_nmod_poly_evaluate_fq_nmod(v, Fcoeff + Fi, alpha, ctx->fqctx);
1459             if (!fq_nmod_is_zero(v, ctx->fqctx))
1460             {
1461                 changed = 1;
1462                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
1463                 fq_nmod_poly_sub(Tcoeff + Ti, Fcoeff + Fi, tp, ctx->fqctx);
1464             }
1465             else
1466             {
1467                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
1468             }
1469             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
1470 
1471             Fi++;
1472         }
1473         else
1474         {
1475             FLINT_ASSERT(Ai < Alen && (Fi >= Flen || mpoly_monomial_lt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)));
1476 
1477             /* F term missing, A term ok */
1478             changed = 1;
1479             fq_nmod_poly_scalar_mul_fq_nmod(Tcoeff + Ti, modulus, (Acoeff + Ai)->coeffs + vi, ctx->fqctx);
1480             mpoly_monomial_set_extra(Texp + N*Ti, Aexp + N*Ai, N, offset, vi << shift);
1481 
1482             do {
1483                 vi--;
1484             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ctx->fqctx));
1485             if (vi < 0)
1486             {
1487                 Ai++;
1488                 if (Ai < Alen)
1489                 {
1490                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ctx->fqctx);
1491                 }
1492             }
1493         }
1494 
1495         FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + Ti, ctx->fqctx));
1496         lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Tcoeff + Ti, ctx->fqctx));
1497         Ti++;
1498     }
1499     T->length = Ti;
1500 
1501     if (changed)
1502     {
1503         fq_nmod_mpolyn_swap(T, F);
1504     }
1505 
1506     fq_nmod_poly_clear(tp, ctx->fqctx);
1507     fq_nmod_clear(v, ctx->fqctx);
1508 
1509     FLINT_ASSERT(fq_nmod_mpolyn_is_canonical(F, ctx));
1510 
1511     *lastdeg_ = lastdeg;
1512     return changed;
1513 }
1514 
1515 /*****************************************************************************/
1516 
1517 /*
1518     A is in              Fq [X][v]
1519     E is in (Fq[v]/alpha(v))[X]
1520 */
fq_nmod_mpolyn_interp_reduce_lg_poly(fq_nmod_poly_t E,const fq_nmod_mpoly_ctx_t ectx,fq_nmod_mpolyn_t A,const fq_nmod_mpoly_ctx_t ctx,const bad_fq_nmod_embed_t emb)1521 void fq_nmod_mpolyn_interp_reduce_lg_poly(
1522     fq_nmod_poly_t E,
1523     const fq_nmod_mpoly_ctx_t ectx,
1524     fq_nmod_mpolyn_t A,
1525     const fq_nmod_mpoly_ctx_t ctx,
1526     const bad_fq_nmod_embed_t emb)
1527 {
1528     fq_nmod_t v;
1529     slong Ai, Alen, k;
1530     fq_nmod_poly_struct * Acoeff;
1531     ulong * Aexp;
1532     slong N, off, shift;
1533 
1534     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1535     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
1536 
1537     fq_nmod_init(v, ectx->fqctx);
1538 
1539     Acoeff = A->coeffs;
1540     Aexp = A->exps;
1541     Alen = A->length;
1542 
1543     fq_nmod_poly_zero(E, ectx->fqctx);
1544     for (Ai = 0; Ai < Alen; Ai++)
1545     {
1546         k = (Aexp + N*Ai)[off] >> shift;
1547         bad_fq_nmod_embed_sm_to_lg(v, Acoeff + Ai, emb);
1548         fq_nmod_poly_set_coeff(E, k, v, ectx->fqctx);
1549     }
1550 
1551     fq_nmod_clear(v, ectx->fqctx);
1552 }
1553 
1554 
1555 /*
1556     A = B
1557     A is in              Fq [X][][v]
1558     B is in (Fq[v]/alpha(v))[X]
1559 */
fq_nmod_mpolyn_interp_lift_lg_poly(slong * lastdeg_,fq_nmod_mpolyn_t A,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_poly_t B,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)1560 void fq_nmod_mpolyn_interp_lift_lg_poly(
1561     slong * lastdeg_,
1562     fq_nmod_mpolyn_t A,
1563     const fq_nmod_mpoly_ctx_t ctx,
1564     fq_nmod_poly_t B,
1565     const fq_nmod_mpoly_ctx_t ectx,
1566     const bad_fq_nmod_embed_t emb)
1567 {
1568     slong Bexp;
1569     slong Blen = fq_nmod_poly_length(B, ectx->fqctx);
1570     fq_nmod_struct * Bcoeff = B->coeffs;
1571     fq_nmod_poly_struct * Acoeff;
1572     ulong * Aexp;
1573     slong Ai;
1574     slong lastdeg = -WORD(1);
1575     slong N, off, shift;
1576 
1577     N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1578     mpoly_gen_offset_shift_sp(&off, &shift, 0, A->bits, ctx->minfo);
1579 
1580     fq_nmod_mpolyn_fit_length(A, Blen, ctx);
1581     Acoeff = A->coeffs;
1582     Aexp = A->exps;
1583 
1584     Ai = 0;
1585     for (Bexp = Blen - 1; Bexp >= 0; Bexp--)
1586     {
1587         if (!fq_nmod_is_zero(Bcoeff + Bexp, ectx->fqctx))
1588         {
1589             FLINT_ASSERT(Ai < A->alloc);
1590             bad_fq_nmod_embed_lg_to_sm(Acoeff + Ai, Bcoeff + Bexp, emb);
1591             FLINT_ASSERT(!fq_nmod_poly_is_zero(Acoeff + Ai, ctx->fqctx));
1592             lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Acoeff + Ai, ctx->fqctx));
1593             mpoly_monomial_zero(Aexp + N*Ai, N);
1594             (Aexp + N*Ai)[off] = Bexp << shift;
1595             Ai++;
1596         }
1597     }
1598     A->length = Ai;
1599 
1600     *lastdeg_ = lastdeg;
1601 }
1602 
1603 
1604 /*
1605     F = F + modulus*((A - F(alpha))/(modulus(alpha)))
1606     F is in              Fq [X][][v]
1607     A is in (Fq[v]/alpha(v))[X]
1608     alpha(v) is ectx->fqctx->modulus(v)
1609 */
fq_nmod_mpolyn_interp_crt_lg_poly(slong * lastdeg_,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,fq_nmod_poly_t modulus,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_poly_t A,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)1610 int fq_nmod_mpolyn_interp_crt_lg_poly(
1611     slong * lastdeg_,
1612     fq_nmod_mpolyn_t F,
1613     fq_nmod_mpolyn_t T,
1614     fq_nmod_poly_t modulus,
1615     const fq_nmod_mpoly_ctx_t ctx,
1616     fq_nmod_poly_t A,
1617     const fq_nmod_mpoly_ctx_t ectx,
1618     const bad_fq_nmod_embed_t emb)
1619 {
1620     int changed = 0;
1621     slong lastdeg = -WORD(1);
1622     fq_nmod_t u, v;
1623     fq_nmod_poly_t u_sm, w;
1624     slong Fi, Ti, Aexp;
1625     fq_nmod_struct * Acoeff = A->coeffs;
1626     slong Flen = F->length;
1627     fq_nmod_poly_struct * Fcoeff = F->coeffs;
1628     ulong * Fexp = F->exps;
1629     fq_nmod_poly_struct * Tcoeff;
1630     ulong * Texp;
1631     fq_nmod_poly_t tp;
1632     fq_nmod_t inv_m_eval;
1633     slong N, off, shift;
1634 
1635     FLINT_ASSERT(T->bits == F->bits);
1636 
1637     N = mpoly_words_per_exp_sp(F->bits, ctx->minfo);
1638     mpoly_gen_offset_shift_sp(&off, &shift, 0, F->bits, ctx->minfo);
1639 
1640     fq_nmod_init(inv_m_eval, ectx->fqctx);
1641     bad_fq_nmod_embed_sm_to_lg(inv_m_eval, modulus, emb);
1642     fq_nmod_inv(inv_m_eval, inv_m_eval, ectx->fqctx);
1643 
1644     fq_nmod_init(u, ectx->fqctx);
1645     fq_nmod_init(v, ectx->fqctx);
1646     fq_nmod_poly_init(u_sm, ctx->fqctx);
1647     fq_nmod_poly_init(w, ctx->fqctx);
1648 
1649     Fi = 0;
1650 
1651     Aexp = fq_nmod_poly_degree(A, ectx->fqctx);
1652 
1653     fq_nmod_poly_init(tp, ctx->fqctx);
1654 
1655     fq_nmod_mpolyn_fit_length(T, Flen + Aexp + 1, ctx);
1656     Tcoeff = T->coeffs;
1657     Texp = T->exps;
1658     Ti = 0;
1659 
1660     while (Fi < Flen || Aexp >= 0)
1661     {
1662         FLINT_ASSERT(Ti < T->alloc);
1663 
1664         if (Fi < Flen)
1665         {
1666             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeff + Fi, ctx->fqctx));
1667             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeff + Fi, ctx->fqctx) < fq_nmod_poly_degree(modulus, ctx->fqctx));
1668         }
1669 
1670         if (Aexp >= 0)
1671         {
1672             FLINT_ASSERT(!fq_nmod_is_zero(Acoeff + Aexp, ectx->fqctx));
1673         }
1674 
1675         mpoly_monomial_zero(Texp + N*Ti, N);
1676 
1677         if (Fi < Flen && Aexp >= 0 && ((Fexp + N*Fi)[off]>>shift) == Aexp)
1678         {
1679             /* F term ok, A term ok */
1680             bad_fq_nmod_embed_sm_to_lg(u, Fcoeff + Fi, emb);
1681             fq_nmod_sub(v, Acoeff + Aexp, u, ectx->fqctx);
1682             if (!fq_nmod_is_zero(v, ectx->fqctx))
1683             {
1684                 changed = 1;
1685                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
1686                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
1687                 fq_nmod_poly_mul(w, modulus, u_sm, ctx->fqctx);
1688                 fq_nmod_poly_add(Tcoeff + Ti, Fcoeff + Fi, w, ctx->fqctx);
1689             }
1690             else
1691             {
1692                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
1693             }
1694 
1695             (Texp + N*Ti)[off] = (Fexp + N*Fi)[off];
1696             Fi++;
1697             do {
1698                 Aexp--;
1699             } while (Aexp >= 0 && fq_nmod_is_zero(Acoeff + Aexp, ectx->fqctx));
1700         }
1701         else if (Fi < Flen && (Aexp < 0 || ((Fexp + N*Fi)[off]>>shift) > Aexp))
1702         {
1703             /* F term ok, A term missing */
1704             bad_fq_nmod_embed_sm_to_lg(v, Fcoeff + Fi, emb);
1705             if (!fq_nmod_is_zero(v, ectx->fqctx))
1706             {
1707                 changed = 1;
1708                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
1709                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
1710                 fq_nmod_poly_mul(w, u_sm, modulus, ctx->fqctx);
1711                 fq_nmod_poly_sub(Tcoeff + Ti, Fcoeff + Fi, w, ctx->fqctx);
1712             }
1713             else
1714             {
1715                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
1716             }
1717 
1718             (Texp + N*Ti)[off] = (Fexp + N*Fi)[off];
1719             Fi++;
1720         }
1721         else if (Aexp >= 0 && (Fi >= Flen || ((Fexp + N*Fi)[off]>>shift) < Aexp))
1722         {
1723             /* F term missing, A term ok */
1724             changed = 1;
1725             fq_nmod_mul(u, Acoeff + Aexp, inv_m_eval, ectx->fqctx);
1726             bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
1727             fq_nmod_poly_mul(Tcoeff + Ti, modulus, u_sm, ctx->fqctx);
1728 
1729             (Texp + N*Ti)[off] = Aexp << shift;
1730             do {
1731                 Aexp--;
1732             } while (Aexp >= 0 && fq_nmod_is_zero(Acoeff + Aexp, ectx->fqctx));
1733         }
1734         else
1735         {
1736             FLINT_ASSERT(0);
1737         }
1738 
1739         FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + Ti, ctx->fqctx));
1740         lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Tcoeff + Ti, ctx->fqctx));
1741         Ti++;
1742     }
1743     T->length = Ti;
1744 
1745     if (changed)
1746     {
1747         fq_nmod_mpolyn_swap(T, F);
1748     }
1749 
1750     fq_nmod_clear(u, ectx->fqctx);
1751     fq_nmod_clear(v, ectx->fqctx);
1752     fq_nmod_poly_clear(u_sm, ctx->fqctx);
1753     fq_nmod_poly_clear(w, ctx->fqctx);
1754 
1755     fq_nmod_clear(inv_m_eval, ectx->fqctx);
1756 
1757     *lastdeg_ = lastdeg;
1758     return changed;
1759 }
1760 
1761 
1762 
1763 /*****************************************************************************/
1764 
1765 /*
1766     A is in                      Fq [x_0,...,x_(var-2), x_(var-1)][x_var]
1767     E is in (Fq[x_var]/alpha(x_var))[x_0,...,x_(var-2)][x_(var-1)]
1768 */
1769 
fq_nmod_mpolyn_interp_reduce_lg_mpolyn(fq_nmod_mpolyn_t E,const fq_nmod_mpoly_ctx_t ectx,fq_nmod_mpolyn_t A,slong var,const fq_nmod_mpoly_ctx_t ctx,const bad_fq_nmod_embed_t emb)1770 void fq_nmod_mpolyn_interp_reduce_lg_mpolyn(
1771     fq_nmod_mpolyn_t E,
1772     const fq_nmod_mpoly_ctx_t ectx,
1773     fq_nmod_mpolyn_t A,
1774     slong var,
1775     const fq_nmod_mpoly_ctx_t ctx,
1776     const bad_fq_nmod_embed_t emb)
1777 {
1778     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1779     slong offset, shift, k;
1780     ulong mask;
1781     fq_nmod_t v;
1782     fq_nmod_poly_struct * Acoeff = A->coeffs;
1783     ulong * Aexp = A->exps;
1784     slong Alen = A->length;
1785     slong Ai;
1786     fq_nmod_poly_struct * Ecoeff;
1787     ulong * Eexp;
1788     slong Ei;
1789 
1790     fq_nmod_init(v, ectx->fqctx);
1791 
1792     FLINT_ASSERT(var > 0);
1793     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1794     mask = (-UWORD(1)) >> (FLINT_BITS - A->bits);
1795 
1796     Ecoeff = E->coeffs;
1797     Eexp = E->exps;
1798     Ei = 0;
1799     for (Ai = 0; Ai < Alen; Ai++)
1800     {
1801         bad_fq_nmod_embed_sm_to_lg(v, Acoeff + Ai, emb);
1802         k = ((Aexp + N*Ai)[offset] >> shift) & mask;
1803         if (fq_nmod_is_zero(v, ectx->fqctx))
1804         {
1805             continue;
1806         }
1807 
1808         if (Ei > 0 && mpoly_monomial_equal_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)))
1809         {
1810             /* append to previous */
1811             fq_nmod_poly_set_coeff(Ecoeff + Ei - 1, k, v, ectx->fqctx);
1812         }
1813         else
1814         {
1815             FLINT_ASSERT(Ei == 0 || mpoly_monomial_gt_nomask_extra(Eexp + N*(Ei - 1), Aexp + N*Ai, N, offset, -(k << shift)));
1816 
1817             /* create new */
1818             if (Ei >= E->alloc)
1819             {
1820                 fq_nmod_mpolyn_fit_length(E, Ei + 1, ectx);
1821                 Ecoeff = E->coeffs;
1822                 Eexp = E->exps;
1823             }
1824             mpoly_monomial_set_extra(Eexp + N*Ei, Aexp + N*Ai, N, offset, -(k << shift));
1825             fq_nmod_poly_zero(Ecoeff + Ei, ectx->fqctx);
1826             fq_nmod_poly_set_coeff(Ecoeff + Ei, k, v, ectx->fqctx);
1827             Ei++;
1828         }
1829     }
1830     E->length = Ei;
1831 
1832     fq_nmod_clear(v, ectx->fqctx);
1833 }
1834 
1835 
1836 /*
1837     A = B using lowest degree representative
1838     A is in                      Fq [x_0,...,x_(var-2), x_(var-1)][x_var]
1839     B is in (Fq[x_var]/alpha(x_var))[x_0,...,x_(var-2)][x_(var-1)]
1840     alpha is emb->h
1841 */
fq_nmod_mpolyn_interp_lift_lg_mpolyn(slong * lastdeg_,fq_nmod_mpolyn_t A,slong var,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpolyn_t B,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)1842 void fq_nmod_mpolyn_interp_lift_lg_mpolyn(
1843     slong * lastdeg_,
1844     fq_nmod_mpolyn_t A,
1845     slong var,
1846     const fq_nmod_mpoly_ctx_t ctx,
1847     fq_nmod_mpolyn_t B,
1848     const fq_nmod_mpoly_ctx_t ectx,
1849     const bad_fq_nmod_embed_t emb)
1850 {
1851     slong N = mpoly_words_per_exp_sp(B->bits, ctx->minfo);
1852     slong offset, shift;
1853     slong vi;
1854     fq_nmod_poly_struct * Bcoeff = B->coeffs;
1855     ulong * Bexp = B->exps;
1856     slong Blen = B->length;
1857     slong Bi;
1858     fq_nmod_poly_struct * Acoeff;
1859     ulong * Aexp;
1860     slong Ai;
1861     slong lastdeg = -WORD(1);
1862 
1863     FLINT_ASSERT(A->bits == B->bits);
1864 
1865     fq_nmod_mpolyn_fit_length(A, Blen, ctx);
1866     Acoeff = A->coeffs;
1867     Aexp = A->exps;
1868 
1869     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1870 
1871     Ai = 0;
1872     for (Bi = 0; Bi < Blen; Bi++)
1873     {
1874         if (Ai + (Bcoeff + Bi)->length >= A->alloc)
1875         {
1876             fq_nmod_mpolyn_fit_length(A, Ai + (Bcoeff + Bi)->length, ctx);
1877             Acoeff = A->coeffs;
1878             Aexp = A->exps;
1879         }
1880         for (vi = (Bcoeff + Bi)->length - 1; vi >= 0; vi--)
1881         {
1882             if (!fq_nmod_is_zero((Bcoeff + Bi)->coeffs + vi, ectx->fqctx))
1883             {
1884                 mpoly_monomial_set_extra(Aexp + N*Ai, Bexp + N*Bi, N, offset, vi << shift);
1885                 bad_fq_nmod_embed_lg_to_sm(Acoeff + Ai, (Bcoeff + Bi)->coeffs + vi, emb);
1886                 FLINT_ASSERT(!fq_nmod_poly_is_zero(Acoeff + Ai, ctx->fqctx));
1887                 lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Acoeff + Ai, ctx->fqctx));
1888                 Ai++;
1889             }
1890         }
1891     }
1892     A->length = Ai;
1893 
1894     *lastdeg_ = lastdeg;
1895 }
1896 
1897 
1898 /*
1899     F = F + modulus*((A - F(alpha))/modulus(alpha))
1900     F is in                      Fq [x_0,...,x_(var-2), x_(var-1)][x_var]
1901     A is in (Fq[x_var]/alpha(x_var))[x_0,...,x_(var-2)][x_(var-1)]
1902 */
fq_nmod_mpolyn_interp_crt_lg_mpolyn(slong * lastdeg_,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,fq_nmod_poly_t modulus,slong var,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpolyn_t A,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)1903 int fq_nmod_mpolyn_interp_crt_lg_mpolyn(
1904     slong * lastdeg_,
1905     fq_nmod_mpolyn_t F,
1906     fq_nmod_mpolyn_t T,
1907     fq_nmod_poly_t modulus,
1908     slong var,
1909     const fq_nmod_mpoly_ctx_t ctx,
1910     fq_nmod_mpolyn_t A,
1911     const fq_nmod_mpoly_ctx_t ectx,
1912     const bad_fq_nmod_embed_t emb)
1913 {
1914     int changed = 0;
1915     slong N = mpoly_words_per_exp_sp(A->bits, ctx->minfo);
1916     slong offset, shift;
1917     slong lastdeg = -WORD(1);
1918     slong vi;
1919     fq_nmod_t u, v;
1920     fq_nmod_poly_t u_sm, w;
1921     fq_nmod_t inv_m_eval;
1922 
1923     fq_nmod_poly_struct * Tcoeff;
1924     ulong * Texp;
1925     slong Ti;
1926 
1927     fq_nmod_poly_struct * Acoeff = A->coeffs;
1928     slong Alen = A->length;
1929     ulong * Aexp = A->exps;
1930     slong Ai;
1931 
1932     fq_nmod_poly_struct * Fcoeff = F->coeffs;
1933     slong Flen = F->length;
1934     ulong * Fexp = F->exps;
1935     slong Fi;
1936 
1937     fq_nmod_init(inv_m_eval, ectx->fqctx);
1938     bad_fq_nmod_embed_sm_to_lg(inv_m_eval, modulus, emb);
1939     fq_nmod_inv(inv_m_eval, inv_m_eval, ectx->fqctx);
1940 
1941     fq_nmod_init(u, ectx->fqctx);
1942     fq_nmod_init(v, ectx->fqctx);
1943     fq_nmod_poly_init(u_sm, ctx->fqctx);
1944     fq_nmod_poly_init(w, ctx->fqctx);
1945 
1946     FLINT_ASSERT(var > 0);
1947     FLINT_ASSERT(T->bits == A->bits);
1948     FLINT_ASSERT(F->bits == A->bits);
1949 
1950     mpoly_gen_offset_shift_sp(&offset, &shift, var - 1, A->bits, ctx->minfo);
1951 
1952     Flen = F->length;
1953 
1954     fq_nmod_mpolyn_fit_length(T, FLINT_MAX(Flen, Alen), ctx);
1955     Tcoeff = T->coeffs;
1956     Texp = T->exps;
1957     Ti = 0;
1958 
1959     Fi = Ai = vi = 0;
1960     if (Ai < Alen)
1961     {
1962         vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
1963     }
1964     while (Fi < Flen || Ai < Alen)
1965     {
1966         if (Ti >= T->alloc)
1967         {
1968             fq_nmod_mpolyn_fit_length(T, Ti + FLINT_MAX(Flen - Fi, Alen - Ai), ctx);
1969             Tcoeff = T->coeffs;
1970             Texp = T->exps;
1971         }
1972 
1973         if (Fi < Flen && Ai < Alen && mpoly_monomial_equal_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift))
1974         {
1975             /* F term ok, A term ok */
1976             bad_fq_nmod_embed_sm_to_lg(u, Fcoeff + Fi, emb);
1977             fq_nmod_sub(v, (Acoeff + Ai)->coeffs + vi, u, ectx->fqctx);
1978             if (!fq_nmod_is_zero(v, ectx->fqctx))
1979             {
1980                 changed = 1;
1981                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
1982                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
1983                 fq_nmod_poly_mul(w, modulus, u_sm, ctx->fqctx);
1984                 fq_nmod_poly_add(Tcoeff + Ti, Fcoeff + Fi, w, ctx->fqctx);
1985             }
1986             else
1987             {
1988                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
1989             }
1990             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
1991 
1992             Fi++;
1993             do {
1994                 vi--;
1995             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ectx->fqctx));
1996             if (vi < 0)
1997             {
1998                 Ai++;
1999                 if (Ai < Alen)
2000                 {
2001                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
2002                 }
2003             }
2004         }
2005         else if (Fi < Flen && (Ai >= Alen || mpoly_monomial_gt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)))
2006         {
2007             /* F term ok, A term missing */
2008             bad_fq_nmod_embed_sm_to_lg(v, Fcoeff + Fi, emb);
2009             if (!fq_nmod_is_zero(v, ectx->fqctx))
2010             {
2011                 changed = 1;
2012                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
2013                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2014                 fq_nmod_poly_mul(w, u_sm, modulus, ctx->fqctx);
2015                 fq_nmod_poly_sub(Tcoeff + Ti, Fcoeff + Fi, w, ctx->fqctx);
2016             }
2017             else
2018             {
2019                 fq_nmod_poly_set(Tcoeff + Ti, Fcoeff + Fi, ctx->fqctx);
2020             }
2021             mpoly_monomial_set(Texp + N*Ti, Fexp + N*Fi, N);
2022 
2023             Fi++;
2024         }
2025         else
2026         {
2027             FLINT_ASSERT(Ai < Alen && (Fi >= Flen || mpoly_monomial_lt_nomask_extra(Fexp + N*Fi, Aexp + N*Ai, N, offset, vi << shift)));
2028 
2029             /* F term missing, A term ok */
2030             changed = 1;
2031             fq_nmod_mul(u, (Acoeff + Ai)->coeffs + vi, inv_m_eval, ectx->fqctx);
2032             bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2033             fq_nmod_poly_mul(Tcoeff + Ti, modulus, u_sm, ctx->fqctx);
2034             mpoly_monomial_set_extra(Texp + N*Ti, Aexp + N*Ai, N, offset, vi << shift);
2035 
2036             do {
2037                 vi--;
2038             } while (vi >= 0 && fq_nmod_is_zero((Acoeff + Ai)->coeffs + vi, ectx->fqctx));
2039             if (vi < 0)
2040             {
2041                 Ai++;
2042                 if (Ai < Alen)
2043                 {
2044                     vi = fq_nmod_poly_degree(A->coeffs + Ai, ectx->fqctx);
2045                 }
2046             }
2047         }
2048 
2049         FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + Ti, ctx->fqctx));
2050         lastdeg = FLINT_MAX(lastdeg, fq_nmod_poly_degree(Tcoeff + Ti, ctx->fqctx));
2051         Ti++;
2052     }
2053     T->length = Ti;
2054 
2055     if (changed)
2056     {
2057         fq_nmod_mpolyn_swap(T, F);
2058     }
2059 
2060     fq_nmod_clear(inv_m_eval, ectx->fqctx);
2061     fq_nmod_clear(u, ectx->fqctx);
2062     fq_nmod_clear(v, ectx->fqctx);
2063     fq_nmod_poly_clear(u_sm, ctx->fqctx);
2064     fq_nmod_poly_clear(w, ctx->fqctx);
2065 
2066     FLINT_ASSERT(fq_nmod_mpolyn_is_canonical(F, ctx));
2067 
2068     *lastdeg_ = lastdeg;
2069     return changed;
2070 }
2071 
2072 
2073 
2074 /*****************************************************************************/
2075 
2076 /* evaluate A at lastvar = alpha */
fq_nmod_mpolyn_interp_reduce_sm_mpoly(fq_nmod_mpoly_t B,fq_nmod_mpolyn_t A,fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)2077 void fq_nmod_mpolyn_interp_reduce_sm_mpoly(
2078     fq_nmod_mpoly_t B,
2079     fq_nmod_mpolyn_t A,
2080     fq_nmod_t alpha,
2081     const fq_nmod_mpoly_ctx_t ctx)
2082 {
2083     slong i, k, N;
2084     FLINT_ASSERT(B->bits == A->bits);
2085 
2086     fq_nmod_mpoly_fit_length(B, A->length, ctx);
2087     N = mpoly_words_per_exp(A->bits, ctx->minfo);
2088     k = 0;
2089     for (i = 0; i < A->length; i++)
2090     {
2091         mpoly_monomial_set(B->exps + N*k, A->exps + N*i, N);
2092         fq_nmod_poly_evaluate_fq_nmod(B->coeffs + k, A->coeffs + i, alpha,
2093                                                                    ctx->fqctx);
2094         k += !fq_nmod_is_zero(B->coeffs + k, ctx->fqctx);
2095     }
2096     B->length = k;
2097 }
2098 
fq_nmod_mpolyun_interp_reduce_sm_mpolyu(fq_nmod_mpolyu_t B,fq_nmod_mpolyun_t A,fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)2099 void fq_nmod_mpolyun_interp_reduce_sm_mpolyu(
2100     fq_nmod_mpolyu_t B,
2101     fq_nmod_mpolyun_t A,
2102     fq_nmod_t alpha,
2103     const fq_nmod_mpoly_ctx_t ctx)
2104 {
2105     slong i, k;
2106 
2107     FLINT_ASSERT(B->bits == A->bits);
2108 
2109     fq_nmod_mpolyu_fit_length(B, A->length, ctx);
2110     k = 0;
2111     for (i = 0; i < A->length; i++)
2112     {
2113         B->exps[k] = A->exps[i];
2114         fq_nmod_mpolyn_interp_reduce_sm_mpoly(B->coeffs + k, A->coeffs + i, alpha, ctx);
2115         k += !fq_nmod_mpoly_is_zero(B->coeffs + k, ctx);
2116     }
2117     B->length = k;
2118 }
2119 
2120 
fq_nmod_mpolyn_interp_lift_sm_mpoly(fq_nmod_mpolyn_t A,const fq_nmod_mpoly_t B,const fq_nmod_mpoly_ctx_t ctx)2121 void fq_nmod_mpolyn_interp_lift_sm_mpoly(
2122     fq_nmod_mpolyn_t A,
2123     const fq_nmod_mpoly_t B,
2124     const fq_nmod_mpoly_ctx_t ctx)
2125 {
2126     slong i;
2127     slong N;
2128     fq_nmod_poly_struct * Acoeff;
2129     fq_nmod_struct * Bcoeff;
2130     ulong * Aexp, * Bexp;
2131     slong Blen;
2132 
2133     fq_nmod_mpolyn_fit_bits(A, B->bits, ctx);
2134     A->bits = B->bits;
2135 
2136     Blen = B->length;
2137     fq_nmod_mpolyn_fit_length(A, Blen, ctx);
2138     Acoeff = A->coeffs;
2139     Bcoeff = B->coeffs;
2140     Aexp = A->exps;
2141     Bexp = B->exps;
2142 
2143     N = mpoly_words_per_exp(B->bits, ctx->minfo);
2144 
2145     for (i = 0; i < Blen; i++)
2146     {
2147         fq_nmod_poly_zero(Acoeff + i, ctx->fqctx);
2148         fq_nmod_poly_set_coeff(Acoeff + i, 0, Bcoeff + i, ctx->fqctx);
2149         mpoly_monomial_set(Aexp + N*i, Bexp + N*i, N);
2150     }
2151     A->length = Blen;
2152 }
2153 
fq_nmod_mpolyun_interp_lift_sm_mpolyu(fq_nmod_mpolyun_t A,const fq_nmod_mpolyu_t B,const fq_nmod_mpoly_ctx_t ctx)2154 void fq_nmod_mpolyun_interp_lift_sm_mpolyu(
2155     fq_nmod_mpolyun_t A,
2156     const fq_nmod_mpolyu_t B,
2157     const fq_nmod_mpoly_ctx_t ctx)
2158 {
2159     slong i;
2160 
2161     FLINT_ASSERT(A->bits == B->bits);
2162     fq_nmod_mpolyun_fit_length(A, B->length, ctx);
2163     for (i = 0; i < B->length; i++)
2164     {
2165         A->exps[i] = B->exps[i];
2166         fq_nmod_mpolyn_interp_lift_sm_mpoly(A->coeffs + i, B->coeffs + i, ctx);
2167         FLINT_ASSERT((A->coeffs + i)->bits == B->bits);
2168     }
2169     A->length = B->length;
2170 }
2171 
2172 
2173 /*
2174     F = F + modulus*(A - F(alpha))
2175 */
fq_nmod_mpolyn_interp_crt_sm_mpoly(slong * lastdeg,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,fq_nmod_mpoly_t A,fq_nmod_poly_t modulus,fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)2176 int fq_nmod_mpolyn_interp_crt_sm_mpoly(
2177     slong * lastdeg,
2178     fq_nmod_mpolyn_t F,
2179     fq_nmod_mpolyn_t T,
2180     fq_nmod_mpoly_t A,
2181     fq_nmod_poly_t modulus,
2182     fq_nmod_t alpha,
2183     const fq_nmod_mpoly_ctx_t ctx)
2184 {
2185     int changed = 0;
2186     slong i, j, k;
2187     slong N;
2188     fq_nmod_t v;
2189     flint_bitcnt_t bits = A->bits;
2190     slong Flen = F->length, Alen = A->length;
2191     ulong * Fexp = F->exps, * Aexp = A->exps;
2192     ulong * Texp;
2193     fq_nmod_struct * Acoeff = A->coeffs;
2194     fq_nmod_poly_struct * Fcoeff = F->coeffs;
2195     fq_nmod_poly_struct * Tcoeff;
2196     fq_nmod_poly_t tp;
2197 
2198     FLINT_ASSERT(F->bits == bits);
2199     FLINT_ASSERT(ctx->minfo->ord == ORD_LEX);
2200 
2201     fq_nmod_init(v, ctx->fqctx);
2202     fq_nmod_poly_init(tp, ctx->fqctx);
2203 
2204     fq_nmod_mpolyn_fit_length(T, Flen + Alen, ctx);
2205     Texp = T->exps;
2206     Tcoeff = T->coeffs;
2207 
2208     N = mpoly_words_per_exp(bits, ctx->minfo);
2209 
2210     i = j = k = 0;
2211     while (i < Flen || j < Alen)
2212     {
2213         if (i < Flen && (j >= Alen
2214                         || mpoly_monomial_gt_nomask(Fexp + N*i, Aexp + N*j, N)))
2215         {
2216             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeff + i, ctx->fqctx));
2217             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeff + i, ctx->fqctx)
2218                                    < fq_nmod_poly_degree(modulus, ctx->fqctx));
2219 
2220             /* F term ok, A term missing */
2221             fq_nmod_poly_evaluate_fq_nmod(v, Fcoeff + i, alpha, ctx->fqctx);
2222             if (!fq_nmod_is_zero(v, ctx->fqctx))
2223             {
2224                 changed = 1;
2225                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
2226                 fq_nmod_poly_sub(Tcoeff + k, Fcoeff + i, tp, ctx->fqctx);
2227             } else {
2228                 fq_nmod_poly_set(Tcoeff + k, Fcoeff + i, ctx->fqctx);
2229             }
2230             lastdeg[0] = FLINT_MAX(lastdeg[0],
2231                                   fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2232 
2233             mpoly_monomial_set(Texp + N*k, Fexp + N*i, N);
2234             FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + k, ctx->fqctx));
2235             k++;
2236             i++;
2237         }
2238         else if (j < Alen && (i >= Flen
2239                         || mpoly_monomial_lt_nomask(Fexp + N*i, Aexp + N*j, N)))
2240         {
2241             /* F term missing, A term ok */
2242             if (!fq_nmod_is_zero(Acoeff + j, ctx->fqctx))
2243             {
2244                 changed = 1;
2245                 fq_nmod_poly_zero(Tcoeff + k, ctx->fqctx);
2246                 fq_nmod_poly_scalar_mul_fq_nmod(Tcoeff + k, modulus, Acoeff + j,
2247                                                                    ctx->fqctx);
2248                 lastdeg[0] = FLINT_MAX(lastdeg[0],
2249                                   fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2250                 mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
2251                 k++;
2252             }
2253             j++;
2254         }
2255         else if (i < Flen && j < Alen
2256                              && mpoly_monomial_equal(Fexp + N*i, Aexp + N*j, N))
2257         {
2258             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeff + i, ctx->fqctx));
2259             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeff + i, ctx->fqctx) < fq_nmod_poly_degree(modulus, ctx->fqctx));
2260 
2261             /* F term ok, A term ok */
2262             fq_nmod_poly_evaluate_fq_nmod(v, Fcoeff + i, alpha, ctx->fqctx);
2263             fq_nmod_sub(v, Acoeff + j, v, ctx->fqctx);
2264             if (!fq_nmod_is_zero(v, ctx->fqctx))
2265             {
2266                 changed = 1;
2267                 fq_nmod_poly_scalar_mul_fq_nmod(tp, modulus, v, ctx->fqctx);
2268                 fq_nmod_poly_add(Tcoeff + k, Fcoeff + i, tp, ctx->fqctx);
2269             } else {
2270                 fq_nmod_poly_set(Tcoeff + k, Fcoeff + i, ctx->fqctx);
2271             }
2272             lastdeg[0] = FLINT_MAX(lastdeg[0],
2273                                   fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2274             mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
2275             FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + k, ctx->fqctx));
2276             k++;
2277             i++;
2278             j++;
2279         } else
2280         {
2281             FLINT_ASSERT(0);
2282         }
2283     }
2284     T->length = k;
2285 
2286     if (changed)
2287     {
2288         fq_nmod_mpolyn_swap(T, F);
2289     }
2290 
2291     fq_nmod_poly_clear(tp, ctx->fqctx);
2292     fq_nmod_clear(v, ctx->fqctx);
2293 
2294     return changed;
2295 }
2296 
2297 
fq_nmod_mpolyun_interp_crt_sm_mpolyu(slong * lastdeg,fq_nmod_mpolyun_t F,fq_nmod_mpolyun_t T,fq_nmod_mpolyu_t A,fq_nmod_poly_t modulus,fq_nmod_t alpha,const fq_nmod_mpoly_ctx_t ctx)2298 int fq_nmod_mpolyun_interp_crt_sm_mpolyu(
2299     slong * lastdeg,
2300     fq_nmod_mpolyun_t F,
2301     fq_nmod_mpolyun_t T,
2302     fq_nmod_mpolyu_t A,
2303     fq_nmod_poly_t modulus,
2304     fq_nmod_t alpha,
2305     const fq_nmod_mpoly_ctx_t ctx)
2306 {
2307     int changed = 0;
2308     slong i, j, k;
2309     ulong * Texp;
2310     ulong * Fexp;
2311     ulong * Aexp;
2312     slong Flen;
2313     slong Alen;
2314     fq_nmod_mpolyn_t S;
2315     fq_nmod_mpolyn_struct * Tcoeff;
2316     fq_nmod_mpolyn_struct * Fcoeff;
2317     fq_nmod_mpoly_struct  * Acoeff;
2318     fq_nmod_mpoly_t zero;
2319 
2320     lastdeg[0] = -WORD(1);
2321 
2322     FLINT_ASSERT(F->bits == T->bits);
2323     FLINT_ASSERT(T->bits == A->bits);
2324 
2325     fq_nmod_mpolyn_init(S, F->bits, ctx);
2326 
2327     Flen = F->length;
2328     Alen = A->length;
2329     fq_nmod_mpolyun_fit_length(T, Flen + Alen, ctx);
2330 
2331     Tcoeff = T->coeffs;
2332     Fcoeff = F->coeffs;
2333     Acoeff = A->coeffs;
2334     Texp = T->exps;
2335     Fexp = F->exps;
2336     Aexp = A->exps;
2337 
2338     fq_nmod_mpoly_init(zero, ctx);
2339     fq_nmod_mpoly_fit_bits(zero, A->bits, ctx);
2340     zero->bits = A->bits;
2341 
2342     i = j = k = 0;
2343     while (i < Flen || j < Alen)
2344     {
2345         if (i < Flen && (j >= Alen || Fexp[i] > Aexp[j]))
2346         {
2347             /* F term ok, A term missing */
2348             fq_nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
2349             changed |= fq_nmod_mpolyn_interp_crt_sm_mpoly(lastdeg, Tcoeff + k, S, zero, modulus, alpha, ctx);
2350             Texp[k] = Fexp[i];
2351             k++;
2352             i++;
2353         }
2354         else if (j < Alen && (i >= Flen || Aexp[j] > Fexp[i]))
2355         {
2356             /* F term missing, A term ok */
2357             fq_nmod_mpolyn_zero(Tcoeff + k, ctx);
2358             changed |= fq_nmod_mpolyn_interp_crt_sm_mpoly(lastdeg, Tcoeff + k, S, Acoeff + j, modulus, alpha, ctx);
2359             Texp[k] = Aexp[j];
2360             k++;
2361             j++;
2362         }
2363         else if (i < Flen && j < Alen && (Fexp[i] == Aexp[j]))
2364         {
2365             /* F term ok, A term ok */
2366             fq_nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
2367             changed |= fq_nmod_mpolyn_interp_crt_sm_mpoly(lastdeg, Tcoeff + k, S, Acoeff + j, modulus, alpha, ctx);
2368             Texp[k] = Aexp[j];
2369             FLINT_ASSERT(!fq_nmod_mpolyn_is_zero(Tcoeff + k, ctx));
2370             k++;
2371             i++;
2372             j++;
2373         } else
2374         {
2375             FLINT_ASSERT(0);
2376         }
2377     }
2378     T->length = k;
2379 
2380     if (changed)
2381     {
2382         fq_nmod_mpolyun_swap(T, F);
2383     }
2384 
2385     fq_nmod_mpolyn_clear(S, ctx);
2386     fq_nmod_mpoly_clear(zero, ctx);
2387     return changed;
2388 }
2389 
2390 
2391 /*****************************************************************************/
2392 
2393 /* reduce B via the map F_q[x] -> F_q^n */
fq_nmod_mpolyn_interp_reduce_lg_mpoly(fq_nmod_mpoly_t A,fq_nmod_mpolyn_t B,const fq_nmod_mpoly_ctx_t ectx,const fq_nmod_mpoly_ctx_t ctx,const bad_fq_nmod_embed_t emb)2394 void fq_nmod_mpolyn_interp_reduce_lg_mpoly(
2395     fq_nmod_mpoly_t A,
2396     fq_nmod_mpolyn_t B,
2397     const fq_nmod_mpoly_ctx_t ectx,
2398     const fq_nmod_mpoly_ctx_t ctx,
2399     const bad_fq_nmod_embed_t emb)
2400 {
2401     slong i, k, N;
2402 
2403     FLINT_ASSERT(A->bits == B->bits);
2404     FLINT_ASSERT(B->bits <= FLINT_BITS);
2405     FLINT_ASSERT(ctx->minfo->ord == ORD_LEX);
2406     FLINT_ASSERT(ectx->minfo->ord == ORD_LEX);
2407     FLINT_ASSERT(ectx->minfo->nvars == ctx->minfo->nvars);
2408 
2409     N = mpoly_words_per_exp_sp(B->bits, ctx->minfo);
2410 
2411     k = 0;
2412     fq_nmod_mpoly_fit_length(A, k + 1, ectx);
2413     for (i = 0; i < B->length; i++)
2414     {
2415         fq_nmod_mpoly_fit_length(A, k + 1, ectx);
2416         mpoly_monomial_set(A->exps + N*k, B->exps + N*i, N);
2417         bad_fq_nmod_embed_sm_to_lg(A->coeffs + k, B->coeffs + i, emb);
2418         k += !fq_nmod_is_zero(A->coeffs + k, ectx->fqctx);
2419     }
2420 
2421     _fq_nmod_mpoly_set_length(A, k, ectx);
2422 }
2423 
2424 /* reduce B via the map from F_q[x] -> F_q^n */
fq_nmod_mpolyun_interp_reduce_lg_mpolyu(fq_nmod_mpolyu_t A,fq_nmod_mpolyun_t B,const fq_nmod_mpoly_ctx_t ectx,const fq_nmod_mpoly_ctx_t ctx,const bad_fq_nmod_embed_t emb)2425 void fq_nmod_mpolyun_interp_reduce_lg_mpolyu(
2426     fq_nmod_mpolyu_t A,
2427     fq_nmod_mpolyun_t B,
2428     const fq_nmod_mpoly_ctx_t ectx,
2429     const fq_nmod_mpoly_ctx_t ctx,
2430     const bad_fq_nmod_embed_t emb)
2431 {
2432     slong i, k, Blen;
2433     fq_nmod_mpoly_struct * Acoeff;
2434     fq_nmod_mpolyn_struct * Bcoeff;
2435     ulong * Aexp, * Bexp;
2436 
2437     Blen = B->length;
2438     fq_nmod_mpolyu_fit_length(A, Blen, ectx);
2439     Acoeff = A->coeffs;
2440     Bcoeff = B->coeffs;
2441     Aexp = A->exps;
2442     Bexp = B->exps;
2443 
2444     k = 0;
2445     for (i = 0; i < Blen; i++)
2446     {
2447         fq_nmod_mpolyn_interp_reduce_lg_mpoly(Acoeff + k, Bcoeff + i, ectx, ctx, emb);
2448         Aexp[k] = Bexp[i];
2449         k += !fq_nmod_mpoly_is_zero(Acoeff + k, ectx);
2450     }
2451     A->length = k;
2452 }
2453 
2454 /* Convert B to A using the map  F_q^n -> F_q[x] */
fq_nmod_mpolyn_interp_lift_lg_mpoly(fq_nmod_mpolyn_t A,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpoly_t B,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)2455 void fq_nmod_mpolyn_interp_lift_lg_mpoly(
2456     fq_nmod_mpolyn_t A,
2457     const fq_nmod_mpoly_ctx_t ctx,
2458     fq_nmod_mpoly_t B,
2459     const fq_nmod_mpoly_ctx_t ectx,
2460     const bad_fq_nmod_embed_t emb)
2461 {
2462     slong i, N;
2463 
2464     FLINT_ASSERT(B->bits == A->bits);
2465     fq_nmod_mpolyn_fit_length(A, B->length, ctx);
2466     N = mpoly_words_per_exp(A->bits, ctx->minfo);
2467     for (i = 0; i < B->length; i++)
2468     {
2469         mpoly_monomial_set(A->exps + N*i, B->exps + N*i, N);
2470         bad_fq_nmod_embed_lg_to_sm(A->coeffs + i, B->coeffs + i, emb);
2471         FLINT_ASSERT(!fq_nmod_poly_is_zero(A->coeffs + i, ctx->fqctx));
2472     }
2473     A->length = B->length;
2474 }
2475 
fq_nmod_mpolyun_interp_lift_lg_mpolyu(fq_nmod_mpolyun_t A,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpolyu_t B,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)2476 void fq_nmod_mpolyun_interp_lift_lg_mpolyu(
2477     fq_nmod_mpolyun_t A,
2478     const fq_nmod_mpoly_ctx_t ctx,
2479     fq_nmod_mpolyu_t B,
2480     const fq_nmod_mpoly_ctx_t ectx,
2481     const bad_fq_nmod_embed_t emb)
2482 {
2483     slong i;
2484 
2485     FLINT_ASSERT(B->bits == A->bits);
2486     fq_nmod_mpolyun_fit_length(A, B->length, ctx);
2487     for (i = 0; i < B->length; i++)
2488     {
2489         A->exps[i] = B->exps[i];
2490         fq_nmod_mpolyn_interp_lift_lg_mpoly(A->coeffs + i, ctx, B->coeffs + i, ectx, emb);
2491         FLINT_ASSERT(!fq_nmod_mpolyn_is_zero(A->coeffs + i, ctx));
2492     }
2493     A->length = B->length;
2494 }
2495 
2496 
2497 /*
2498     update F so that it doesn't change mod m and is A mod emb->h
2499 
2500     F = F + m*((A - F(alpha))/(m(alpha)))
2501 
2502     no assumptions about matching monomials
2503 */
fq_nmod_mpolyn_interp_crt_lg_mpoly(slong * lastdeg,fq_nmod_mpolyn_t F,fq_nmod_mpolyn_t T,fq_nmod_poly_t m,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpoly_t A,fq_nmod_t inv_m_eval,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)2504 int fq_nmod_mpolyn_interp_crt_lg_mpoly(
2505     slong * lastdeg,
2506     fq_nmod_mpolyn_t F,
2507     fq_nmod_mpolyn_t T,
2508     fq_nmod_poly_t m,
2509     const fq_nmod_mpoly_ctx_t ctx,
2510     fq_nmod_mpoly_t A,
2511     fq_nmod_t inv_m_eval,
2512     const fq_nmod_mpoly_ctx_t ectx,
2513     const bad_fq_nmod_embed_t emb)
2514 {
2515     int changed = 0;
2516     slong i, j, k;
2517     slong N;
2518     fq_nmod_t u, v;
2519     fq_nmod_poly_t u_sm, w;
2520     flint_bitcnt_t bits = A->bits;
2521     slong Flen = F->length, Alen = A->length;
2522     ulong * Fexp = F->exps, * Aexp = A->exps;
2523     ulong * Texp;
2524     fq_nmod_struct * Acoeff = A->coeffs;
2525     fq_nmod_poly_struct * Fcoeff = F->coeffs;
2526     fq_nmod_poly_struct * Tcoeff;
2527 
2528     FLINT_ASSERT(bits <= FLINT_BITS);
2529     FLINT_ASSERT(F->bits == bits);
2530     FLINT_ASSERT(ctx->minfo->ord == ORD_LEX);
2531 
2532     fq_nmod_init(u, ectx->fqctx);
2533     fq_nmod_init(v, ectx->fqctx);
2534     fq_nmod_poly_init(u_sm, ctx->fqctx);
2535     fq_nmod_poly_init(w, ctx->fqctx);
2536 
2537     fq_nmod_mpolyn_fit_length(T, Flen + Alen, ctx);
2538     Texp = T->exps;
2539     Tcoeff = T->coeffs;
2540 
2541     N = mpoly_words_per_exp(bits, ctx->minfo);
2542 
2543     i = j = k = 0;
2544     while (i < Flen || j < Alen)
2545     {
2546         if (i < Flen && (j >= Alen
2547                         || mpoly_monomial_gt_nomask(Fexp + N*i, Aexp + N*j, N)))
2548         {
2549             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeff + i, ctx->fqctx));
2550             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeff + i, ctx->fqctx)
2551                                   < fq_nmod_poly_degree(m, ctx->fqctx));
2552 
2553             /* F term ok, A term missing */
2554             bad_fq_nmod_embed_sm_to_lg(v, Fcoeff + i, emb);
2555             if (!fq_nmod_is_zero(v, ectx->fqctx))
2556             {
2557                 changed = 1;
2558                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
2559                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2560                 fq_nmod_poly_mul(w, u_sm, m, ctx->fqctx);
2561                 fq_nmod_poly_sub(Tcoeff + k, Fcoeff + i, w, ctx->fqctx);
2562             }
2563             else
2564             {
2565                 fq_nmod_poly_set(Tcoeff + k, Fcoeff + i, ctx->fqctx);
2566             }
2567             lastdeg[0] = FLINT_MAX(lastdeg[0], fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2568 
2569             mpoly_monomial_set(Texp + N*k, Fexp + N*i, N);
2570             FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + k, ctx->fqctx));
2571             k++;
2572             i++;
2573         }
2574         else if (j < Alen && (i >= Flen
2575                         || mpoly_monomial_lt_nomask(Fexp + N*i, Aexp + N*j, N)))
2576         {
2577             /* F term missing, A term ok */
2578             if (!fq_nmod_is_zero(Acoeff + j, ectx->fqctx))
2579             {
2580                 changed = 1;
2581                 fq_nmod_mul(u, Acoeff + j, inv_m_eval, ectx->fqctx);
2582                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2583                 fq_nmod_poly_mul(Tcoeff + k, m, u_sm, ctx->fqctx);
2584                 lastdeg[0] = FLINT_MAX(lastdeg[0], fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2585                 mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
2586                 k++;
2587             }
2588             j++;
2589         }
2590         else if (i < Flen && j < Alen
2591                              && mpoly_monomial_equal(Fexp + N*i, Aexp + N*j, N))
2592         {
2593             FLINT_ASSERT(!fq_nmod_poly_is_zero(Fcoeff + i, ctx->fqctx));
2594             FLINT_ASSERT(fq_nmod_poly_degree(Fcoeff + i, ctx->fqctx)
2595                                         < fq_nmod_poly_degree(m, ctx->fqctx));
2596 
2597             /* F term ok, A term ok */
2598             bad_fq_nmod_embed_sm_to_lg(u, Fcoeff + i, emb);
2599             fq_nmod_sub(v, Acoeff + j, u, ectx->fqctx);
2600             if (!fq_nmod_is_zero(v, ectx->fqctx))
2601             {
2602                 changed = 1;
2603                 fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
2604                 bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2605                 fq_nmod_poly_mul(w, m, u_sm, ctx->fqctx);
2606                 fq_nmod_poly_add(Tcoeff + k, Fcoeff + i, w, ctx->fqctx);
2607             }
2608             else
2609             {
2610                 fq_nmod_poly_set(Tcoeff + k, Fcoeff + i, ctx->fqctx);
2611             }
2612             lastdeg[0] = FLINT_MAX(lastdeg[0], fq_nmod_poly_degree(Tcoeff + k, ctx->fqctx));
2613             mpoly_monomial_set(Texp + N*k, Aexp + N*j, N);
2614             FLINT_ASSERT(!fq_nmod_poly_is_zero(Tcoeff + k, ctx->fqctx));
2615             k++;
2616             i++;
2617             j++;
2618         }
2619         else
2620         {
2621             FLINT_ASSERT(0);
2622         }
2623     }
2624 
2625     fq_nmod_mpolyn_set_length(T, k, ctx);
2626 
2627     if (changed)
2628     {
2629         fq_nmod_mpolyn_swap(T, F);
2630     }
2631 
2632     fq_nmod_clear(u, ectx->fqctx);
2633     fq_nmod_clear(v, ectx->fqctx);
2634     fq_nmod_poly_clear(u_sm, ctx->fqctx);
2635     fq_nmod_poly_clear(w, ctx->fqctx);
2636 
2637     return changed;
2638 }
2639 
2640 
fq_nmod_mpolyun_interp_crt_lg_mpolyu(slong * lastdeg,fq_nmod_mpolyun_t F,fq_nmod_mpolyun_t T,fq_nmod_poly_t m,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_mpolyu_t A,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)2641 int fq_nmod_mpolyun_interp_crt_lg_mpolyu(
2642     slong * lastdeg,
2643     fq_nmod_mpolyun_t F,
2644     fq_nmod_mpolyun_t T,
2645     fq_nmod_poly_t m,
2646     const fq_nmod_mpoly_ctx_t ctx,
2647     fq_nmod_mpolyu_t A,
2648     const fq_nmod_mpoly_ctx_t ectx,
2649     const bad_fq_nmod_embed_t emb)
2650 {
2651     int changed = 0;
2652     slong i, j, k;
2653     ulong * Texp;
2654     ulong * Fexp;
2655     ulong * Aexp;
2656     slong Flen;
2657     slong Alen;
2658     fq_nmod_mpolyn_t S;
2659     fq_nmod_mpolyn_struct * Tcoeff;
2660     fq_nmod_mpolyn_struct * Fcoeff;
2661     fq_nmod_mpoly_struct  * Acoeff;
2662     fq_nmod_mpoly_t zero;
2663     fq_nmod_t inv_m_eval;
2664 
2665     lastdeg[0] = -WORD(1);
2666 
2667     FLINT_ASSERT(F->bits == T->bits);
2668     FLINT_ASSERT(T->bits == A->bits);
2669 
2670     fq_nmod_mpolyn_init(S, F->bits, ctx);
2671 
2672     Flen = F->length;
2673     Alen = A->length;
2674     fq_nmod_mpolyun_fit_length(T, Flen + Alen, ctx);
2675 
2676     Tcoeff = T->coeffs;
2677     Fcoeff = F->coeffs;
2678     Acoeff = A->coeffs;
2679     Texp = T->exps;
2680     Fexp = F->exps;
2681     Aexp = A->exps;
2682 
2683     fq_nmod_mpoly_init(zero, ectx);
2684     fq_nmod_mpoly_fit_bits(zero, A->bits, ectx);
2685     zero->bits = A->bits;
2686 
2687     fq_nmod_init(inv_m_eval, ectx->fqctx);
2688     bad_fq_nmod_embed_sm_to_lg(inv_m_eval, m, emb);
2689     fq_nmod_inv(inv_m_eval, inv_m_eval, ectx->fqctx);
2690 
2691     i = j = k = 0;
2692     while (i < Flen || j < Alen)
2693     {
2694         if (i < Flen && (j >= Alen || Fexp[i] > Aexp[j]))
2695         {
2696             /* F term ok, A term missing */
2697             fq_nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
2698             changed |= fq_nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
2699                                        S, m, ctx, zero, inv_m_eval, ectx, emb);
2700             Texp[k] = Fexp[i];
2701             FLINT_ASSERT(!fq_nmod_mpolyn_is_zero(Tcoeff + k, ctx));
2702             k++;
2703             i++;
2704         }
2705         else if (j < Alen && (i >= Flen || Aexp[j] > Fexp[i]))
2706         {
2707             /* F term missing, A term ok */
2708             fq_nmod_mpolyn_zero(Tcoeff + k, ctx);
2709             changed |= fq_nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
2710                                  S, m, ctx, Acoeff + j, inv_m_eval, ectx, emb);
2711             Texp[k] = Aexp[j];
2712             FLINT_ASSERT(!fq_nmod_mpolyn_is_zero(Tcoeff + k, ctx));
2713             k++;
2714             j++;
2715         }
2716         else if (i < Flen && j < Alen && (Fexp[i] == Aexp[j]))
2717         {
2718             /* F term ok, A term ok */
2719             fq_nmod_mpolyn_set(Tcoeff + k, Fcoeff + i, ctx);
2720             changed |= fq_nmod_mpolyn_interp_crt_lg_mpoly(lastdeg, Tcoeff + k,
2721                                  S, m, ctx, Acoeff + j, inv_m_eval, ectx, emb);
2722             Texp[k] = Aexp[j];
2723             FLINT_ASSERT(!fq_nmod_mpolyn_is_zero(Tcoeff + k, ctx));
2724             k++;
2725             i++;
2726             j++;
2727         }
2728         else
2729         {
2730             FLINT_ASSERT(0);
2731         }
2732     }
2733 
2734     T->length = k;
2735 
2736     if (changed)
2737     {
2738         fq_nmod_mpolyun_swap(T, F);
2739     }
2740 
2741     fq_nmod_clear(inv_m_eval, ectx->fqctx);
2742 
2743     fq_nmod_mpolyn_clear(S, ctx);
2744     fq_nmod_mpoly_clear(zero, ectx);
2745     return changed;
2746 }
2747 
2748 /*
2749     Update H so that it does not change mod m, and is now A mod p
2750     It is asserted that the monomials in H and A match
2751 */
fq_nmod_mpolyn_interp_mcrt_lg_mpoly(slong * lastdeg,fq_nmod_mpolyn_t H,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_poly_t m,fq_nmod_t inv_m_eval,fq_nmod_mpoly_t A,const fq_nmod_mpoly_ctx_t ectx,const bad_fq_nmod_embed_t emb)2752 int fq_nmod_mpolyn_interp_mcrt_lg_mpoly(
2753     slong * lastdeg,
2754     fq_nmod_mpolyn_t H,
2755     const fq_nmod_mpoly_ctx_t ctx,
2756     fq_nmod_poly_t m,
2757     fq_nmod_t inv_m_eval,
2758     fq_nmod_mpoly_t A,
2759     const fq_nmod_mpoly_ctx_t ectx,
2760     const bad_fq_nmod_embed_t emb)
2761 {
2762     slong i;
2763 #if WANT_ASSERT
2764     slong N;
2765 #endif
2766     int changed = 0;
2767     fq_nmod_t u, v;
2768     fq_nmod_poly_t w, u_sm;
2769 
2770     fq_nmod_init(u, ectx->fqctx);
2771     fq_nmod_init(v, ectx->fqctx);
2772     fq_nmod_poly_init(w, ctx->fqctx);
2773     fq_nmod_poly_init(u_sm, ctx->fqctx);
2774 
2775     FLINT_ASSERT(H->length == A->length);
2776     FLINT_ASSERT(H->bits == A->bits);
2777 #if WANT_ASSERT
2778     N = mpoly_words_per_exp(A->bits, ctx->minfo);
2779 #endif
2780     for (i = 0; i < A->length; i++)
2781     {
2782         FLINT_ASSERT(mpoly_monomial_equal(H->exps + N*i, A->exps + N*i, N));
2783         bad_fq_nmod_embed_sm_to_lg(u, H->coeffs + i, emb);
2784         fq_nmod_sub(v, A->coeffs + i, u, ectx->fqctx);
2785         if (!fq_nmod_is_zero(v, ectx->fqctx))
2786         {
2787             changed = 1;
2788             fq_nmod_mul(u, v, inv_m_eval, ectx->fqctx);
2789             bad_fq_nmod_embed_lg_to_sm(u_sm, u, emb);
2790             fq_nmod_poly_mul(w, u_sm, m, ctx->fqctx);
2791             fq_nmod_poly_add(H->coeffs + i, H->coeffs + i, w, ctx->fqctx);
2792         }
2793 
2794         *lastdeg = FLINT_MAX(*lastdeg,
2795                                fq_nmod_poly_degree(H->coeffs + i, ctx->fqctx));
2796 
2797         FLINT_ASSERT(fq_nmod_poly_degree(H->coeffs + i, ctx->fqctx)
2798                          <  fq_nmod_poly_degree(m, ctx->fqctx)
2799                           + fq_nmod_poly_degree(emb->h, ctx->fqctx));
2800     }
2801 
2802     fq_nmod_clear(u, ectx->fqctx);
2803     fq_nmod_clear(v, ectx->fqctx);
2804     fq_nmod_poly_clear(w, ctx->fqctx);
2805     fq_nmod_poly_clear(u_sm, ctx->fqctx);
2806 
2807     return changed;
2808 }
2809 
fq_nmod_mpolyun_interp_mcrt_lg_mpolyu(slong * lastdeg,fq_nmod_mpolyun_t H,const fq_nmod_mpoly_ctx_t ctx,fq_nmod_poly_t m,fq_nmod_mpolyu_t A,const fq_nmod_mpoly_ctx_t ectx,bad_fq_nmod_embed_t emb)2810 int fq_nmod_mpolyun_interp_mcrt_lg_mpolyu(
2811     slong * lastdeg,
2812     fq_nmod_mpolyun_t H,
2813     const fq_nmod_mpoly_ctx_t ctx,
2814     fq_nmod_poly_t m,
2815     fq_nmod_mpolyu_t A,
2816     const fq_nmod_mpoly_ctx_t ectx,
2817     bad_fq_nmod_embed_t emb)
2818 {
2819     slong i;
2820     int changed = 0;
2821     fq_nmod_t inv_m_eval;
2822 
2823     *lastdeg = -WORD(1);
2824 
2825     fq_nmod_init(inv_m_eval, ectx->fqctx);
2826     bad_fq_nmod_embed_sm_to_lg(inv_m_eval, m, emb);
2827     fq_nmod_inv(inv_m_eval, inv_m_eval, ectx->fqctx);
2828 
2829     FLINT_ASSERT(H->bits == A->bits);
2830     FLINT_ASSERT(H->length == A->length);
2831     for (i = 0; i < A->length; i++)
2832     {
2833         FLINT_ASSERT(H->exps[i] == A->exps[i]);
2834         changed |= fq_nmod_mpolyn_interp_mcrt_lg_mpoly(lastdeg, H->coeffs + i,
2835                                  ctx, m, inv_m_eval, A->coeffs + i, ectx, emb);
2836     }
2837     H->length = A->length;
2838     fq_nmod_clear(inv_m_eval, ectx->fqctx);
2839     return changed;
2840 }
2841 
2842 
2843