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 <http://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.n = ctx->ffinfo->mod.n;
87 Au->mod.ninv = ctx->ffinfo->mod.ninv;
88 Au->mod.norm = ctx->ffinfo->mod.norm;
89
90 Bu->alloc = Bd->coeff_alloc;
91 Bu->coeffs = Bd->coeffs;
92 Bu->length = nmod_mpolyd_length(Bd);
93 Bu->mod.n = ctx->ffinfo->mod.n;
94 Bu->mod.ninv = ctx->ffinfo->mod.ninv;
95 Bu->mod.norm = ctx->ffinfo->mod.norm;
96
97 /* manually move P to Pu */
98 Pu->alloc = Pd->coeff_alloc;
99 Pu->coeffs = Pd->coeffs;
100 Pu->length = 0;
101 Pu->mod.n = ctx->ffinfo->mod.n;
102 Pu->mod.ninv = ctx->ffinfo->mod.ninv;
103 Pu->mod.norm = ctx->ffinfo->mod.norm;
104
105 nmod_poly_mul(Pu, Au, Bu);
106
107 /* manually move Pu to P */
108 Pd->coeff_alloc = Pu->alloc;
109 Pd->coeffs = Pu->coeffs;
110 for (i = Pu->length; i < Pd->coeff_alloc; i++)
111 Pd->coeffs[i] = UWORD(0);
112
113 nmod_mpolyd_clear(Bd);
114 nmod_mpolyd_clear(Ad);
115 nmod_mpoly_convert_from_nmod_mpolyd(P, ctx, Pd, dctx);
116 nmod_mpolyd_clear(Pd);
117
118 nmod_mpolyd_ctx_clear(dctx);
119
120 done:
121 TMP_END;
122 return success;
123
124 failed_stage2:
125 nmod_mpolyd_clear(Ad);
126 nmod_mpolyd_clear(Bd);
127 nmod_mpolyd_clear(Pd);
128
129 failed_stage1:
130 nmod_mpolyd_ctx_clear(dctx);
131 success = 0;
132 goto done;
133 }
134
135
nmod_mpoly_mul_dense(nmod_mpoly_t A,const nmod_mpoly_t B,const nmod_mpoly_t C,const nmod_mpoly_ctx_t ctx)136 int nmod_mpoly_mul_dense(nmod_mpoly_t A, const nmod_mpoly_t B,
137 const nmod_mpoly_t C, const nmod_mpoly_ctx_t ctx)
138 {
139 slong i;
140 int success;
141 fmpz * maxBfields, * maxCfields;
142 TMP_INIT;
143
144 if (B->length == 0 || C->length == 0)
145 {
146 nmod_mpoly_zero(A, ctx);
147 return 1;
148 }
149
150 if (B->bits > FLINT_BITS || C->bits > FLINT_BITS)
151 {
152 return 0;
153 }
154
155 TMP_START;
156
157 maxBfields = (fmpz *) TMP_ALLOC(ctx->minfo->nfields*sizeof(fmpz));
158 maxCfields = (fmpz *) TMP_ALLOC(ctx->minfo->nfields*sizeof(fmpz));
159 for (i = 0; i < ctx->minfo->nfields; i++)
160 {
161 fmpz_init(maxBfields + i);
162 fmpz_init(maxCfields + i);
163 }
164 mpoly_max_fields_fmpz(maxBfields, B->exps, B->length, B->bits, ctx->minfo);
165 mpoly_max_fields_fmpz(maxCfields, C->exps, C->length, C->bits, ctx->minfo);
166
167 success = _nmod_mpoly_mul_dense(A, B, maxBfields, C, maxCfields, ctx);
168
169 for (i = 0; i < ctx->minfo->nfields; i++)
170 {
171 fmpz_clear(maxBfields + i);
172 fmpz_clear(maxCfields + i);
173 }
174
175 TMP_END;
176 return success;
177 }
178