1 /*
2 Copyright (C) 2018 Daniel Schultz
3
4 This file is part of FLINT.
5
6 FLINT is free software: you can redistribute it and/or modify it under
7 the terms of the GNU Lesser General Public License (LGPL) as published
8 by the Free Software Foundation; either version 2.1 of the License, or
9 (at your option) any later version. See <https://www.gnu.org/licenses/>.
10 */
11
12 #include "nmod_poly.h"
13 #include "nmod_mpoly.h"
14
_nmod_mpoly_mul_dense(nmod_mpoly_t P,const nmod_mpoly_t A,fmpz * maxAfields,const nmod_mpoly_t B,fmpz * maxBfields,const nmod_mpoly_ctx_t ctx)15 int _nmod_mpoly_mul_dense(nmod_mpoly_t P,
16 const nmod_mpoly_t A, fmpz * maxAfields,
17 const nmod_mpoly_t B, fmpz * maxBfields,
18 const nmod_mpoly_ctx_t ctx)
19 {
20 int success = 1;
21 slong i;
22 slong nvars = ctx->minfo->nvars;
23 nmod_mpolyd_ctx_t dctx;
24 nmod_mpolyd_t Ad, Bd, Pd;
25 nmod_poly_t Au, Bu, Pu;
26 slong * Abounds, * Bbounds, * Pbounds;
27 TMP_INIT;
28
29 FLINT_ASSERT(A->length != 0);
30 FLINT_ASSERT(B->length != 0);
31
32 FLINT_ASSERT(A->bits <= FLINT_BITS);
33 FLINT_ASSERT(B->bits <= FLINT_BITS);
34
35 TMP_START;
36
37 /* set the ordering of variables for the dense representation */
38 nmod_mpolyd_ctx_init(dctx, nvars);
39
40 /*
41 for each variable v except for the outermost variable,
42 we need to pack to degree deg_v(A) + deg_v(B)
43 */
44 Abounds = (slong *) TMP_ALLOC(ctx->minfo->nvars*sizeof(slong));
45 Bbounds = (slong *) TMP_ALLOC(ctx->minfo->nvars*sizeof(slong));
46 Pbounds = (slong *) TMP_ALLOC(ctx->minfo->nvars*sizeof(slong));
47 mpoly_get_monomial_ui_unpacked_ffmpz((ulong *)Abounds, maxAfields, ctx->minfo);
48 mpoly_get_monomial_ui_unpacked_ffmpz((ulong *)Bbounds, maxBfields, ctx->minfo);
49 for (i = 0; i < ctx->minfo->nvars; i++)
50 {
51 Abounds[i] = Abounds[i] + 1;
52 Bbounds[i] = Bbounds[i] + 1;
53 Pbounds[i] = Abounds[i] + Bbounds[i] - 1;
54 if ((Abounds[i] | Bbounds[i] | Pbounds[i]) < WORD(0))
55 {
56 goto failed_stage1;
57 }
58 if (i != dctx->perm[0])
59 {
60 /* variable of index i is not the outermost */
61 Abounds[i] = Pbounds[i];
62 Bbounds[i] = Pbounds[i];
63 }
64 }
65
66 nmod_mpolyd_init(Ad, nvars);
67 nmod_mpolyd_init(Bd, nvars);
68 nmod_mpolyd_init(Pd, nvars);
69
70 success = 1;
71 success = success && nmod_mpolyd_set_degbounds_perm(Ad, dctx, Abounds);
72 success = success && nmod_mpolyd_set_degbounds_perm(Bd, dctx, Bbounds);
73 success = success && nmod_mpolyd_set_degbounds_perm(Pd, dctx, Pbounds);
74 if (!success)
75 {
76 goto failed_stage2;
77 }
78
79 nmod_mpoly_convert_to_nmod_mpolyd_degbound(Ad, dctx, A, ctx);
80 nmod_mpoly_convert_to_nmod_mpolyd_degbound(Bd, dctx, B, ctx);
81
82 /* let Au and Bu borrow Ad and Bd */
83 Au->alloc = Ad->coeff_alloc;
84 Au->coeffs = Ad->coeffs;
85 Au->length = nmod_mpolyd_length(Ad);
86 Au->mod = ctx->mod;
87
88 Bu->alloc = Bd->coeff_alloc;
89 Bu->coeffs = Bd->coeffs;
90 Bu->length = nmod_mpolyd_length(Bd);
91 Bu->mod = ctx->mod;
92
93 /* manually move P to Pu */
94 Pu->alloc = Pd->coeff_alloc;
95 Pu->coeffs = Pd->coeffs;
96 Pu->length = 0;
97 Pu->mod = ctx->mod;
98
99 nmod_poly_mul(Pu, Au, Bu);
100
101 /* manually move Pu to P */
102 Pd->coeff_alloc = Pu->alloc;
103 Pd->coeffs = Pu->coeffs;
104 for (i = Pu->length; i < Pd->coeff_alloc; i++)
105 Pd->coeffs[i] = UWORD(0);
106
107 nmod_mpolyd_clear(Bd);
108 nmod_mpolyd_clear(Ad);
109 nmod_mpoly_convert_from_nmod_mpolyd(P, ctx, Pd, dctx);
110 nmod_mpolyd_clear(Pd);
111
112 nmod_mpolyd_ctx_clear(dctx);
113
114 done:
115 TMP_END;
116 return success;
117
118 failed_stage2:
119 nmod_mpolyd_clear(Ad);
120 nmod_mpolyd_clear(Bd);
121 nmod_mpolyd_clear(Pd);
122
123 failed_stage1:
124 nmod_mpolyd_ctx_clear(dctx);
125 success = 0;
126 goto done;
127 }
128
129
nmod_mpoly_mul_dense(nmod_mpoly_t A,const nmod_mpoly_t B,const nmod_mpoly_t C,const nmod_mpoly_ctx_t ctx)130 int nmod_mpoly_mul_dense(nmod_mpoly_t A, const nmod_mpoly_t B,
131 const nmod_mpoly_t C, const nmod_mpoly_ctx_t ctx)
132 {
133 slong i;
134 int success;
135 fmpz * maxBfields, * maxCfields;
136 TMP_INIT;
137
138 if (B->length == 0 || C->length == 0)
139 {
140 nmod_mpoly_zero(A, ctx);
141 return 1;
142 }
143
144 if (B->bits > FLINT_BITS || C->bits > FLINT_BITS)
145 {
146 return 0;
147 }
148
149 TMP_START;
150
151 maxBfields = (fmpz *) TMP_ALLOC(ctx->minfo->nfields*sizeof(fmpz));
152 maxCfields = (fmpz *) TMP_ALLOC(ctx->minfo->nfields*sizeof(fmpz));
153 for (i = 0; i < ctx->minfo->nfields; i++)
154 {
155 fmpz_init(maxBfields + i);
156 fmpz_init(maxCfields + i);
157 }
158 mpoly_max_fields_fmpz(maxBfields, B->exps, B->length, B->bits, ctx->minfo);
159 mpoly_max_fields_fmpz(maxCfields, C->exps, C->length, C->bits, ctx->minfo);
160
161 success = _nmod_mpoly_mul_dense(A, B, maxBfields, C, maxCfields, ctx);
162
163 for (i = 0; i < ctx->minfo->nfields; i++)
164 {
165 fmpz_clear(maxBfields + i);
166 fmpz_clear(maxCfields + i);
167 }
168
169 TMP_END;
170 return success;
171 }
172