1 /*
2 Copyright (C) 2016 William Hart
3 Copyright (C) 2017 Daniel Schultz
4
5 This file is part of FLINT.
6
7 FLINT is free software: you can redistribute it and/or modify it under
8 the terms of the GNU Lesser General Public License (LGPL) as published
9 by the Free Software Foundation; either version 2.1 of the License, or
10 (at your option) any later version. See <http://www.gnu.org/licenses/>.
11 */
12
13 #include "fmpz_mpoly.h"
14
_fmpz_mpoly_pow_fps1(fmpz ** poly1,ulong ** exp1,slong * alloc,const fmpz * poly2,const ulong * exp2,slong len2,ulong k,ulong maskhi)15 slong _fmpz_mpoly_pow_fps1(fmpz ** poly1, ulong ** exp1, slong * alloc,
16 const fmpz * poly2, const ulong * exp2, slong len2, ulong k,
17 ulong maskhi)
18 {
19 const slong topbit = (WORD(1) << (FLINT_BITS - 1));
20 const slong mask = ~topbit;
21 slong i, rnext, g_alloc, gnext;
22 slong next_loc;
23 slong next_free, Q_len = 0, heap_len = 2; /* heap zero index unused */
24 mpoly_heap1_s * heap;
25 mpoly_heap_t * chain;
26 mpoly_heap_t ** Q, ** reuse;
27 mpoly_heap_t * x;
28 fmpz * p1 = *poly1, * gc = NULL;
29 ulong * e1 = *exp1, * ge, * fik;
30 ulong exp, finalexp, temp2;
31 slong * largest;
32 fmpz_t t1, C, S, temp1;
33 int first;
34 TMP_INIT;
35
36 TMP_START;
37
38 next_loc = len2 + 4; /* something bigger than heap can ever be */
39 heap = (mpoly_heap1_s *) TMP_ALLOC((len2 + 1)*sizeof(mpoly_heap1_s));
40 /* 2x as we pull from heap and insert more before processing pulled ones */
41 chain = (mpoly_heap_t *) TMP_ALLOC(2*len2*sizeof(mpoly_heap_t));
42 reuse = (mpoly_heap_t **) TMP_ALLOC(2*len2*sizeof(mpoly_heap_t *));
43 Q = (mpoly_heap_t **) TMP_ALLOC(len2*sizeof(mpoly_heap_t *));
44 /* we add 1 to all entries of largest to free up the value 0 */
45 largest = (slong *) TMP_ALLOC(len2*sizeof(slong));
46
47 fmpz_init(t1);
48 fmpz_init(C);
49 fmpz_init(S);
50 fmpz_init(temp1);
51
52 for (i = 0; i < 2*len2; i++)
53 reuse[i] = chain + i;
54
55 g_alloc = k*(len2 - 1) + 1;
56 ge = (ulong *) flint_malloc(g_alloc*sizeof(ulong));
57 gnext = 0;
58 rnext = 0;
59
60 ge[0] = exp2[0]*(k - 1);
61
62 e1[0] = exp2[0]*k;
63
64 gc = (fmpz *) flint_calloc(g_alloc, sizeof(fmpz));
65 fmpz_pow_ui(gc + 0, poly2 + 0, k - 1);
66 fmpz_mul(p1 + 0, gc + 0, poly2 + 0);
67
68 next_free = 0;
69
70 x = reuse[next_free++];
71 x->i = 1;
72 x->j = 0;
73 x->next = NULL;
74
75 HEAP_ASSIGN(heap[1], exp2[1] + ge[0], x);
76
77 for (i = 0; i < len2; i++)
78 largest[i] = topbit;
79 largest[1] = 1;
80
81 fik = (ulong *) TMP_ALLOC(len2*sizeof(ulong));
82
83 for (i = 0; i < len2; i++)
84 fik[i] = exp2[i]*(k - 1);
85
86 finalexp = exp2[0];
87
88 while (heap_len > 1)
89 {
90 exp = heap[1].exp;
91
92 rnext++;
93 gnext++;
94
95 if (rnext >= *alloc)
96 {
97 p1 = (fmpz *) flint_realloc(p1, 2*sizeof(fmpz)*(*alloc));
98 e1 = (ulong *) flint_realloc(e1, 2*sizeof(ulong)*(*alloc));
99 flint_mpn_zero(p1 + *alloc, *alloc);
100 (*alloc) *= 2;
101 }
102
103 if (gnext >= g_alloc)
104 {
105 ge = (ulong *) flint_realloc(ge, 2*sizeof(ulong)*g_alloc);
106 gc = (fmpz *) flint_realloc(gc, 2*sizeof(fmpz)*g_alloc);
107 flint_mpn_zero(gc + g_alloc, g_alloc);
108 g_alloc *= 2;
109 }
110
111 first = 1;
112
113 fmpz_zero(C);
114 fmpz_zero(S);
115
116 while (heap_len > 1 && heap[1].exp == exp)
117 {
118 x = _mpoly_heap_pop1(heap, &heap_len, maskhi);
119
120 largest[x->i] |= topbit;
121
122 fmpz_mul(t1, poly2 + x->i, gc + x->j);
123 fmpz_add(S, S, t1);
124
125 if ((exp^maskhi) >= (finalexp^maskhi))
126 {
127 temp2 = fik[x->i] - ge[x->j];
128
129 if ((slong) temp2 < 0)
130 fmpz_submul_ui(C, t1, -temp2);
131 else
132 fmpz_addmul_ui(C, t1, temp2);
133 }
134
135 if (first)
136 {
137 ge[gnext] = exp - exp2[0];
138 first = 0;
139 }
140
141 Q[Q_len++] = x;
142
143 while ((x = x->next) != NULL)
144 {
145 largest[x->i] |= topbit;
146
147 fmpz_mul(t1, poly2 + x->i, gc + x->j);
148 fmpz_add(S, S, t1);
149
150 if ((exp^maskhi) >= (finalexp^maskhi))
151 {
152 temp2 = fik[x->i] - ge[x->j];
153
154 if (0 > (slong) temp2)
155 fmpz_submul_ui(C, t1, -temp2);
156 else
157 fmpz_addmul_ui(C, t1, temp2);
158 }
159
160 Q[Q_len++] = x;
161 }
162 }
163
164 while (Q_len > 0)
165 {
166 slong i, j;
167
168 x = Q[--Q_len];
169 i = x->i;
170 j = x->j;
171
172 if (i < len2 - 1 && largest[i + 1] == (j | topbit))
173 {
174 x->i++;
175 x->next = NULL;
176
177 _mpoly_heap_insert1(heap, exp2[i + 1] + ge[j], x,
178 &next_loc, &heap_len, maskhi);
179 largest[i + 1] = j + 1;
180 } else
181 reuse[--next_free] = x;
182
183 if (j < gnext - 1 && (largest[i] & mask) < j + 2)
184 {
185 x = reuse[next_free++];
186
187 x->i = i;
188 x->j = j + 1;
189 x->next = NULL;
190
191 _mpoly_heap_insert1(heap, exp2[i] + ge[j + 1], x,
192 &next_loc, &heap_len, maskhi);
193 largest[i] = j + 2;
194 }
195 }
196
197 if (!fmpz_is_zero(C))
198 {
199 slong t2 = exp - k*exp2[0];
200
201 if (t2 < 0)
202 {
203 fmpz_divexact_ui(temp1, C, -t2);
204 fmpz_neg(temp1, temp1);
205 } else
206 fmpz_divexact_ui(temp1, C, t2);
207
208 fmpz_add(S, S, temp1);
209 fmpz_divexact(gc + gnext, temp1, poly2 + 0);
210
211 if ((largest[1] & topbit) != 0)
212 {
213 x = reuse[next_free++];
214
215 x->i = 1;
216 x->j = gnext;
217 x->next = NULL;
218
219 _mpoly_heap_insert1(heap, exp2[1] + ge[gnext], x,
220 &next_loc, &heap_len, maskhi);
221
222 largest[1] = gnext + 1;
223 }
224 }
225
226 if (!fmpz_is_zero(S))
227 {
228 fmpz_set(p1 + rnext, S);
229 e1[rnext] = ge[gnext] + exp2[0];
230 } else
231 rnext--;
232
233 if (fmpz_is_zero(C))
234 gnext--;
235 }
236
237 rnext++;
238
239 (*poly1) = p1;
240 (*exp1) = e1;
241
242 fmpz_clear(t1);
243 fmpz_clear(C);
244 fmpz_clear(S);
245 fmpz_clear(temp1);
246
247 flint_free(ge);
248 for (i = 0; i < g_alloc; i++)
249 fmpz_clear(gc + i);
250 flint_free(gc);
251
252 TMP_END;
253
254 return rnext;
255 }
256
fmpz_set_mpn(fmpz_t f,ulong * c_in,slong n)257 void fmpz_set_mpn(fmpz_t f, ulong * c_in, slong n)
258 {
259 slong i;
260 ulong * c;
261
262 TMP_INIT;
263
264 TMP_START;
265
266 c = (ulong *) TMP_ALLOC(n*sizeof(ulong));
267
268 for (i = 0; i < n; i++)
269 c[i] = c_in[i];
270
271 while (n > 0 && c[n - 1] == 0)
272 n--;
273
274 if (n <= 1)
275 fmpz_set_ui(f, c[0]);
276 else
277 {
278 __mpz_struct * mpz = _fmpz_promote(f);
279
280 mpz_realloc2(mpz, n*FLINT_BITS);
281
282 mpn_copyi(mpz->_mp_d, c, n);
283 mpz->_mp_size = n;
284 }
285
286 TMP_END;
287 }
288
289
fmpz_set_mpn_signed(fmpz_t f,ulong * c_in,slong n)290 void fmpz_set_mpn_signed(fmpz_t f, ulong * c_in, slong n)
291 {
292 slong i;
293 ulong * c;
294 int neg;
295
296 TMP_INIT;
297
298 TMP_START;
299
300 c = (ulong *) TMP_ALLOC(n*sizeof(ulong));
301
302 for (i = 0; i < n; i++)
303 c[i] = c_in[i];
304
305 neg = 0 > (slong) c[n - 1];
306
307 if (neg)
308 mpn_neg_n(c, c, n);
309
310 while (n > 0 && c[n - 1] == 0)
311 n--;
312
313 if (n <= 1)
314 fmpz_set_ui(f, c[0]);
315 else
316 {
317 __mpz_struct * mpz = _fmpz_promote(f);
318
319 mpz_realloc2(mpz, n*FLINT_BITS);
320
321 mpn_copyi(mpz->_mp_d, c, n);
322 mpz->_mp_size = n;
323 }
324
325 if (neg)
326 fmpz_neg(f, f);
327
328 TMP_END;
329 }
330
331
_fmpz_mpoly_pow_fps(fmpz ** poly1,ulong ** exp1,slong * alloc,const fmpz * poly2,const ulong * exp2,slong len2,ulong k,flint_bitcnt_t bits,slong N,const ulong * cmpmask)332 slong _fmpz_mpoly_pow_fps(fmpz ** poly1, ulong ** exp1, slong * alloc,
333 const fmpz * poly2, const ulong * exp2, slong len2, ulong k,
334 flint_bitcnt_t bits, slong N, const ulong * cmpmask)
335 {
336 const slong topbit = (WORD(1) << (FLINT_BITS - 1));
337 const slong mask = ~topbit;
338 slong i, rnext, g_alloc, gnext, exp_next;
339 slong next_loc;
340 slong next_free, Q_len = 0, heap_len = 2; /* heap zero index unused */
341 mpoly_heap_s * heap;
342 mpoly_heap_t * chain;
343 mpoly_heap_t ** Q, ** reuse;
344 mpoly_heap_t * x;
345 fmpz * p1 = *poly1, * gc = NULL;
346 ulong * e1 = *exp1, * ge, * fik, * exp, * exps, * exp_copy;
347 ulong ** exp_list;
348 ulong * finalexp, * temp2;
349 slong * largest;
350 fmpz_t t1, t2, C, S, temp1;
351 int first;
352 TMP_INIT;
353
354 if (N == 1)
355 return _fmpz_mpoly_pow_fps1(poly1, exp1, alloc, poly2, exp2, len2, k,
356 cmpmask[0]);
357
358 TMP_START;
359
360 next_loc = len2 + 4; /* something bigger than heap can ever be */
361 heap = (mpoly_heap_s *) TMP_ALLOC((len2 + 1)*sizeof(mpoly_heap_s));
362 /* 2x as we pull from heap and insert more before processing pulled ones */
363 chain = (mpoly_heap_t *) TMP_ALLOC(2*len2*sizeof(mpoly_heap_t));
364 reuse = (mpoly_heap_t **) TMP_ALLOC(2*len2*sizeof(mpoly_heap_t *));
365 Q = (mpoly_heap_t **) TMP_ALLOC(len2*sizeof(mpoly_heap_t *));
366 /* we add 1 to all entries of largest to free up the value 0 */
367 largest = (slong *) TMP_ALLOC(len2*sizeof(slong));
368 exps = (ulong *) TMP_ALLOC((len2 + 1)*N*sizeof(ulong));
369 exp_list = (ulong **) TMP_ALLOC((len2 + 1)*sizeof(ulong *));
370 finalexp = (ulong *) TMP_ALLOC(N*sizeof(ulong));
371 temp2 = (ulong *) TMP_ALLOC(N*sizeof(ulong));
372 exp_copy = (ulong *) TMP_ALLOC(N*sizeof(ulong));
373
374 fmpz_init(t1);
375 fmpz_init(t2);
376 fmpz_init(C);
377 fmpz_init(S);
378 fmpz_init(temp1);
379
380 for (i = 0; i < len2 + 1; i++)
381 exp_list[i] = exps + i*N;
382
383 exp_next = 0;
384
385 for (i = 0; i < 2*len2; i++)
386 reuse[i] = chain + i;
387
388 g_alloc = k*(len2 - 1) + 1;
389 ge = (ulong *) flint_malloc(g_alloc*sizeof(ulong)*N);
390 gnext = 0;
391 rnext = 0;
392
393 if (bits <= FLINT_BITS)
394 {
395 mpoly_monomial_mul_ui(ge + 0, exp2 + 0, N, k - 1);
396 mpoly_monomial_mul_ui(e1 + 0, exp2 + 0, N, k);
397 } else
398 {
399 mpoly_monomial_mul_ui_mp(ge + 0, exp2 + 0, N, k - 1);
400 mpoly_monomial_mul_ui_mp(e1 + 0, exp2 + 0, N, k);
401 }
402
403 gc = (fmpz *) flint_calloc(g_alloc, sizeof(fmpz));
404 fmpz_pow_ui(gc + 0, poly2 + 0, k - 1);
405 fmpz_mul(p1 + 0, gc + 0, poly2 + 0);
406
407 next_free = 0;
408
409 x = reuse[next_free++];
410 x->i = 1;
411 x->j = 0;
412 x->next = NULL;
413
414 heap[1].next = x;
415 heap[1].exp = exp_list[exp_next++];
416
417 if (bits <= FLINT_BITS)
418 mpoly_monomial_add(heap[1].exp, exp2 + N, ge + 0, N);
419 else
420 mpoly_monomial_add_mp(heap[1].exp, exp2 + N, ge + 0, N);
421
422 for (i = 0; i < len2; i++)
423 largest[i] = topbit;
424 largest[1] = 1;
425
426 fik = (ulong *) TMP_ALLOC(N*len2*sizeof(ulong));
427
428 for (i = 0; i < len2; i++)
429 {
430 if (bits <= FLINT_BITS)
431 mpoly_monomial_mul_ui(fik + i*N, exp2 + i*N, N, k - 1);
432 else
433 mpoly_monomial_mul_ui_mp(fik + i*N, exp2 + i*N, N, k - 1);
434 }
435
436 mpoly_monomial_set(finalexp, exp2, N);
437
438 while (heap_len > 1)
439 {
440 exp = heap[1].exp;
441 mpoly_monomial_set(exp_copy, exp, N);
442
443 rnext++;
444 gnext++;
445
446 if (rnext >= *alloc)
447 {
448 p1 = (fmpz *) flint_realloc(p1, 2*sizeof(fmpz)*(*alloc));
449 e1 = (ulong *) flint_realloc(e1, 2*N*sizeof(ulong)*(*alloc));
450 flint_mpn_zero(p1 + *alloc, *alloc);
451 (*alloc) *= 2;
452 }
453
454 if (gnext >= g_alloc)
455 {
456 ge = (ulong *) flint_realloc(ge, 2*N*sizeof(ulong)*g_alloc);
457 gc = (fmpz *) flint_realloc(gc, 2*sizeof(fmpz)*g_alloc);
458 flint_mpn_zero(gc + g_alloc, g_alloc);
459 g_alloc *= 2;
460 }
461
462 first = 1;
463
464 fmpz_zero(C);
465 fmpz_zero(S);
466
467 while (heap_len > 1 && mpoly_monomial_equal(heap[1].exp, exp, N))
468 {
469 exp_list[--exp_next] = heap[1].exp;
470 x = _mpoly_heap_pop(heap, &heap_len, N, cmpmask);
471
472 largest[x->i] |= topbit;
473
474 fmpz_mul(t1, poly2 + x->i, gc + x->j);
475 fmpz_add(S, S, t1);
476
477 if (!mpoly_monomial_gt(finalexp, exp, N, cmpmask))
478 {
479 mpn_sub_n(temp2, fik + x->i*N, ge + x->j*N, N);
480 fmpz_set_mpn_signed(t2, temp2, N);
481 fmpz_addmul(C, t1, t2);
482 }
483
484 if (first)
485 {
486 if (bits <= FLINT_BITS)
487 mpoly_monomial_sub(ge + gnext*N, exp, exp2 + 0, N);
488 else
489 mpoly_monomial_sub_mp(ge + gnext*N, exp, exp2 + 0, N);
490
491 first = 0;
492 }
493
494 Q[Q_len++] = x;
495
496 while ((x = x->next) != NULL)
497 {
498 largest[x->i] |= topbit;
499
500 fmpz_mul(t1, poly2 + x->i, gc + x->j);
501 fmpz_add(S, S, t1);
502
503 if (!mpoly_monomial_gt(finalexp, exp, N, cmpmask))
504 {
505 mpn_sub_n(temp2, fik + x->i*N, ge + x->j*N, N);
506 fmpz_set_mpn_signed(t2, temp2, N);
507 fmpz_addmul(C, t1, t2);
508 }
509
510 Q[Q_len++] = x;
511 }
512 }
513
514 while (Q_len > 0)
515 {
516 slong i, j;
517
518 x = Q[--Q_len];
519 i = x->i;
520 j = x->j;
521
522 if (i < len2 - 1 && largest[i + 1] == (j | topbit))
523 {
524 x->i++;
525 x->next = NULL;
526
527 if (bits <= FLINT_BITS)
528 mpoly_monomial_add(exp_list[exp_next], exp2 + (i + 1)*N, ge + j*N, N);
529 else
530 mpoly_monomial_add_mp(exp_list[exp_next], exp2 + (i + 1)*N, ge + j*N, N);
531
532 if (!_mpoly_heap_insert(heap, exp_list[exp_next++], x,
533 &next_loc, &heap_len, N, cmpmask))
534 exp_next--;
535
536 largest[i + 1] = j + 1;
537 } else
538 reuse[--next_free] = x;
539
540 if (j < gnext - 1 && (largest[i] & mask) < j + 2)
541 {
542 x = reuse[next_free++];
543
544 x->i = i;
545 x->j = j + 1;
546 x->next = NULL;
547
548 if (bits <= FLINT_BITS)
549 mpoly_monomial_add(exp_list[exp_next], exp2 + i*N, ge + (j + 1)*N, N);
550 else
551 mpoly_monomial_add_mp(exp_list[exp_next], exp2 + i*N, ge + (j + 1)*N, N);
552
553 if (!_mpoly_heap_insert(heap, exp_list[exp_next++], x,
554 &next_loc, &heap_len, N, cmpmask))
555 exp_next--;
556
557 largest[i] = j + 2;
558 }
559 }
560
561 if (!fmpz_is_zero(C))
562 {
563 if (bits <= FLINT_BITS)
564 mpoly_monomial_mul_ui(temp2, exp2 + 0, N, k);
565 else
566 mpoly_monomial_mul_ui_mp(temp2, exp2 + 0, N, k);
567
568 mpn_sub_n(temp2, exp_copy, temp2, N);
569 fmpz_set_mpn_signed(t2, temp2, N);
570 fmpz_divexact(temp1, C, t2);
571 fmpz_add(S, S, temp1);
572 fmpz_divexact(gc + gnext, temp1, poly2 + 0);
573
574 if ((largest[1] & topbit) != 0)
575 {
576 x = reuse[next_free++];
577
578 x->i = 1;
579 x->j = gnext;
580 x->next = NULL;
581
582 if (bits <= FLINT_BITS)
583 mpoly_monomial_add(exp_list[exp_next], exp2 + N, ge + gnext*N, N);
584 else
585 mpoly_monomial_add_mp(exp_list[exp_next], exp2 + N, ge + gnext*N, N);
586
587 if (!_mpoly_heap_insert(heap, exp_list[exp_next++], x,
588 &next_loc, &heap_len, N, cmpmask))
589 exp_next--;
590
591 largest[1] = gnext + 1;
592 }
593 }
594
595 if (!fmpz_is_zero(S))
596 {
597 fmpz_set(p1 + rnext, S);
598 if (bits <= FLINT_BITS)
599 mpoly_monomial_add(e1 + rnext*N, ge + gnext*N, exp2 + 0, N);
600 else
601 mpoly_monomial_add_mp(e1 + rnext*N, ge + gnext*N, exp2 + 0, N);
602
603 } else
604 rnext--;
605
606 if (fmpz_is_zero(C))
607 gnext--;
608 }
609
610 rnext++;
611
612 (*poly1) = p1;
613 (*exp1) = e1;
614
615 fmpz_clear(t1);
616 fmpz_clear(t2);
617 fmpz_clear(C);
618 fmpz_clear(S);
619 fmpz_clear(temp1);
620
621 flint_free(ge);
622 for (i = 0; i < g_alloc; i++)
623 fmpz_clear(gc + i);
624 flint_free(gc);
625
626 TMP_END;
627
628 return rnext;
629 }
630
fmpz_mpoly_pow_fps(fmpz_mpoly_t A,const fmpz_mpoly_t B,ulong k,const fmpz_mpoly_ctx_t ctx)631 void fmpz_mpoly_pow_fps(fmpz_mpoly_t A, const fmpz_mpoly_t B,
632 ulong k, const fmpz_mpoly_ctx_t ctx)
633 {
634 slong i, N, len = 0;
635 fmpz * maxBfields;
636 flint_bitcnt_t exp_bits;
637 ulong * cmpmask;
638 ulong * Bexp = B->exps;
639 int freeBexp = 0;
640 TMP_INIT;
641
642 FLINT_ASSERT(k >= 2);
643 FLINT_ASSERT(B->length > 0);
644
645 TMP_START;
646
647 maxBfields = (fmpz *) TMP_ALLOC(ctx->minfo->nfields*sizeof(fmpz));
648 for (i = 0; i < ctx->minfo->nfields; i++)
649 fmpz_init(maxBfields + i);
650
651 mpoly_max_fields_fmpz(maxBfields, B->exps, B->length, B->bits, ctx->minfo);
652 _fmpz_vec_scalar_mul_ui(maxBfields, maxBfields, ctx->minfo->nfields, k);
653
654 exp_bits = _fmpz_vec_max_bits(maxBfields, ctx->minfo->nfields);
655 exp_bits = FLINT_MAX(MPOLY_MIN_BITS, exp_bits + 1);
656 exp_bits = FLINT_MAX(exp_bits, B->bits);
657 exp_bits = mpoly_fix_bits(exp_bits, ctx->minfo);
658
659 for (i = 0; i < ctx->minfo->nfields; i++)
660 fmpz_clear(maxBfields + i);
661
662 N = mpoly_words_per_exp(exp_bits, ctx->minfo);
663 cmpmask = (ulong*) TMP_ALLOC(N*sizeof(ulong));
664 mpoly_get_cmpmask(cmpmask, N, exp_bits, ctx->minfo);
665
666 if (exp_bits > B->bits)
667 {
668 freeBexp = 1;
669 Bexp = (ulong *) flint_malloc(N*B->length*sizeof(ulong));
670 mpoly_repack_monomials(Bexp, exp_bits, B->exps, B->bits,
671 B->length, ctx->minfo);
672 }
673
674 if (B->length == 1)
675 {
676 fmpz_mpoly_fit_length(A, 1, ctx);
677 fmpz_mpoly_fit_bits(A, exp_bits, ctx);
678 A->bits = exp_bits;
679
680 fmpz_pow_ui(A->coeffs + 0, B->coeffs + 0, k);
681
682 if (exp_bits <= FLINT_BITS)
683 mpoly_monomial_mul_ui(A->exps, Bexp, N, k);
684 else
685 mpoly_monomial_mul_ui_mp(A->exps, Bexp, N, k);
686
687 len = 1;
688 goto cleanup;
689 }
690
691 if (A == B)
692 {
693 fmpz_mpoly_t T;
694
695 fmpz_mpoly_init2(T, k*(B->length - 1) + 1, ctx);
696 fmpz_mpoly_fit_bits(T, exp_bits, ctx);
697 T->bits = exp_bits;
698
699 len = _fmpz_mpoly_pow_fps(&T->coeffs, &T->exps, &T->alloc,
700 B->coeffs, Bexp, B->length, k, exp_bits, N, cmpmask);
701
702 fmpz_mpoly_swap(T, A, ctx);
703 fmpz_mpoly_clear(T, ctx);
704 }
705 else
706 {
707 fmpz_mpoly_fit_length(A, k*(B->length - 1) + 1, ctx);
708 fmpz_mpoly_fit_bits(A, exp_bits, ctx);
709 A->bits = exp_bits;
710
711 len = _fmpz_mpoly_pow_fps(&A->coeffs, &A->exps, &A->alloc,
712 B->coeffs, Bexp, B->length, k, exp_bits, N, cmpmask);
713 }
714
715 cleanup:
716
717 if (freeBexp)
718 flint_free(Bexp);
719
720 _fmpz_mpoly_set_length(A, len, ctx);
721
722 TMP_END;
723 }
724