1 /*
2     Copyright (C) 2020 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 "mpoly.h"
13 
14 /* unfortunate function missing from fmpz_mat */
fmpz_mat_mul_vec(fmpz * v,const fmpz_mat_t M,fmpz * u)15 void fmpz_mat_mul_vec(fmpz * v, const fmpz_mat_t M, fmpz * u)
16 {
17     slong i;
18     slong r = fmpz_mat_nrows(M);
19     slong c = fmpz_mat_ncols(M);
20 
21     for (i = 0; i < r; i++)
22         _fmpz_vec_dot(v + i, M->rows[i], u, c);
23 }
24 
25 /* Fill the compose matrix for the T_mpoly_compose_T_mpoly_gen functions */
mpoly_compose_mat_gen(fmpz_mat_t M,const slong * c,const mpoly_ctx_t mctxB,const mpoly_ctx_t mctxAC)26 void mpoly_compose_mat_gen(fmpz_mat_t M, const slong * c,
27                              const mpoly_ctx_t mctxB, const mpoly_ctx_t mctxAC)
28 {
29     slong i, gi, j;
30     fmpz * t;
31 
32     FLINT_ASSERT(fmpz_mat_nrows(M) == mctxAC->nfields + 1);
33     FLINT_ASSERT(fmpz_mat_ncols(M) == mctxB->nfields);
34 
35     fmpz_mat_zero(M);
36 
37     t = _fmpz_vec_init(mctxAC->nfields);
38 
39     for (i = 0; i < mctxB->nvars; i++)
40     {
41         gi = mpoly_gen_index(i, mctxB);
42         if (0 <= c[i] && c[i] < mctxAC->nfields)
43         {
44             mpoly_gen_fields_fmpz(t, c[i], mctxAC);
45             for (j = 0; j < mctxAC->nfields; j++)
46                 fmpz_swap(fmpz_mat_entry(M, j, gi), t + j);
47         }
48         else
49         {
50             /* c[i] corresponds to zero */
51             fmpz_one(fmpz_mat_entry(M, mctxAC->nfields, gi));
52         }
53     }
54 
55     _fmpz_vec_clear(t, mctxAC->nfields);
56 }
57 
58 /*
59     Fill the column of the compose matrix for the variable Bvar assuming
60     Bvar maps to the monomial in (Cexps, Cbits).
61     NULL may be passed in for Cexp to indicate that Bbar maps to 0.
62 */
mpoly_compose_mat_fill_column(fmpz_mat_t M,const ulong * Cexp,flint_bitcnt_t Cbits,slong Bvar,const mpoly_ctx_t mctxB,const mpoly_ctx_t mctxAC)63 void mpoly_compose_mat_fill_column(fmpz_mat_t M, const ulong * Cexp,
64                      flint_bitcnt_t Cbits, slong Bvar, const mpoly_ctx_t mctxB,
65                                                       const mpoly_ctx_t mctxAC)
66 {
67     slong Bidx = mpoly_gen_index(Bvar, mctxB);
68     slong j;
69     fmpz * t;
70 
71     if (Cexp == NULL)
72     {
73         j = mctxAC->nfields;
74         fmpz_one(fmpz_mat_entry(M, j, Bidx));
75 
76         for (j--; j >= 0; j--)
77             fmpz_zero(fmpz_mat_entry(M, j, Bidx));
78 
79         return;
80     }
81 
82     t = _fmpz_vec_init(mctxAC->nfields);
83 
84     mpoly_unpack_vec_fmpz(t, Cexp, Cbits, mctxAC->nfields, 1);
85 
86     j = mctxAC->nfields;
87     fmpz_zero(fmpz_mat_entry(M, j, Bidx));
88     for (j--; j >= 0; j--)
89         fmpz_swap(fmpz_mat_entry(M, j, Bidx), t + j);
90 
91     _fmpz_vec_clear(t, mctxAC->nfields);
92 }
93