| Peter Bruin on Thu, 16 Jan 2014 14:25:11 +0100 |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
| Re: new FFM_mul (and FlxqM_mul, FqM_mul, ...) |
Hello, Bill Allombert <Bill.Allombert@math.u-bordeaux.fr> wrote: >> The more I think about it, the more I actually prefer this approach to >> the code duplication in my patch, so if you want to do it in this way, >> I'll certainly be happy with that. > > Yes. Would you provide a patch for gen_matmul? Here is a new patch using the bb_field interface. In addition to the previous patch, it implements FFM_FFC_mul, gen_matmul and gen_matcolmul; I think we now have all matrix * matrix and matrix * column vector functions over finite fields. The matrix * matrix functions are covered by the current test suite, but the matrix * column vector functions are not. Best wishes, Peter PS: I just reported a bug (matrix multiplication 0x0 * 0xN -> 0x0 instead of 0xN); this patch is based on the one attached to that bug report.
diff --git a/doc/develop.tex b/doc/develop.tex
index e5cfcae4..092c392 100644
--- a/doc/develop.tex
+++ b/doc/develop.tex
@@ -434,6 +434,10 @@ fields:
\fun{GEN}{gen_matid}{long n, void *E, const struct bb_field *S}
+\fun{GEN}{gen_matcolmul}{GEN a, GEN b, void *E, const struct bb_field *ff}
+
+\fun{GEN}{gen_matmul}{GEN a, GEN b, void *E, const struct bb_field *ff}
+
\subsec{Functions returning black box fields}
\fun{const struct bb_field *}{get_Fp_field}{void **pt_E, GEN p}
diff --git a/doc/usersch5.tex b/doc/usersch5.tex
index 42680b3..06787cc 100644
--- a/doc/usersch5.tex
+++ b/doc/usersch5.tex
@@ -3757,6 +3757,8 @@ as \kbd{gauss}, where $b$ is a \kbd{FqM}.
\fun{GEN}{FqM_FqC_gauss}{GEN a, GEN b, GEN T, GEN p}
as \kbd{gauss}, where $b$ is a \kbd{FqC}.
+\fun{GEN}{FqM_FqC_mul}{GEN a, GEN b, GEN T, GEN p}
+
\fun{GEN}{FqM_ker}{GEN x, GEN T, GEN p} as \kbd{ker}
\fun{GEN}{FqM_image}{GEN x, GEN T, GEN p} as \kbd{image}
@@ -3764,6 +3766,8 @@ as \kbd{gauss}, where $b$ is a \kbd{FqC}.
\fun{GEN}{FqM_inv}{GEN x, GEN T, GEN p} returns the inverse of \kbd{x}, or
\kbd{NULL} if \kbd{x} is not invertible.
+\fun{GEN}{FqM_FqC_mul}{GEN a, GEN b, GEN T, GEN p}
+
\fun{long}{FqM_rank}{GEN x, GEN T, GEN p} as \kbd{rank}
\fun{GEN}{FqM_suppl}{GEN x, GEN T, GEN p} as \kbd{suppl}
@@ -4033,6 +4037,8 @@ and \kbd{y}
\fun{GEN}{FlxqM_FlxqC_gauss}{GEN a, GEN b, GEN T, ulong p}
+\fun{GEN}{FlxqM_FlxqC_mul}{GEN a, GEN b, GEN T, ulong p}
+
\fun{GEN}{FlxqM_ker}{GEN x, GEN T, ulong p}
\fun{GEN}{FlxqM_image}{GEN x, GEN T, ulong p}
@@ -4041,6 +4047,8 @@ and \kbd{y}
\fun{GEN}{FlxqM_inv}{GEN x, GEN T, ulong p}
+\fun{GEN}{FlxqM_mul}{GEN a, GEN b, GEN T, ulong p}
+
\fun{long}{FlxqM_rank}{GEN x, GEN T, ulong p}
\fun{GEN}{matid_FlxqM}{long n, GEN T, ulong p}
@@ -5395,6 +5403,8 @@ $a=\sigma(X)$ where $\sigma$ is an automorphism of the algebra $\F_2[X]/T(X)$.
\subsec{\kbd{F2xqV}, \kbd{F2xqM}}. See \kbd{FqV}, \kbd{FqM} operations.
+\fun{GEN}{F2xqM_F2xqC_mul}{GEN a, GEN b, GEN T}
+
\fun{GEN}{F2xqM_ker}{GEN x, GEN T}
\fun{GEN}{F2xqM_det}{GEN a, GEN T}
@@ -5403,6 +5413,8 @@ $a=\sigma(X)$ where $\sigma$ is an automorphism of the algebra $\F_2[X]/T(X)$.
\fun{GEN}{F2xqM_inv}{GEN a, GEN T}
+\fun{GEN}{F2xqM_mul}{GEN a, GEN b, GEN T}
+
\fun{long}{F2xqM_rank}{GEN x, GEN T}
\fun{GEN}{matid_F2xqM}{long n, GEN T}
@@ -9090,9 +9102,13 @@ of the univariate polynomial \kbd{f} over the definition field of the
\typ{FFELT} \kbd{a}. The coefficients of \kbd{f} must be of type \typ{INT},
\typ{INTMOD} or \typ{FFELT} and compatible with \kbd{a}.
-\fun{GEN}{FFM_ker}{GEN M, GEN ff} return the kernel of the \typ{MAT} \kbd{M}
+\fun{GEN}{FFM_FFC_mul}{GEN M, GEN C, GEN ff} returns the product of
+the matrix~\kbd{M} (\typ{MAT}) and the column vector~\kbd{C}
+(\typ{COL}) over the finite field given by \kbd{ff} (\typ{FFELT}).
+
+\fun{GEN}{FFM_ker}{GEN M, GEN ff} returns the kernel of the \typ{MAT} \kbd{M}
defined over the finite field given by the \typ{FFELT} \kbd{ff} (obtained
-by \tet{RgM_is_FFM(M,\&ff)}.
+by \tet{RgM_is_FFM(M,\&ff)}).
\fun{GEN}{FFM_det}{GEN M, GEN ff}
@@ -9100,6 +9116,10 @@ by \tet{RgM_is_FFM(M,\&ff)}.
\fun{GEN}{FFM_inv}{GEN M, GEN ff}
+\fun{GEN}{FFM_mul}{GEN M, GEN N, GEN ff} returns the product of the
+matrices \kbd{M} and~\kbd{N} (\typ{MAT}) over the finite field given
+by \kbd{ff} (\typ{FFELT}).
+
\fun{long}{FFM_rank}{GEN M, GEN ff}
\section{Transcendental functions}
diff --git a/src/basemath/FF.c b/src/basemath/FF.c
index 00aee57..c5116b6 100644
--- a/src/basemath/FF.c
+++ b/src/basemath/FF.c
@@ -1738,3 +1738,39 @@ FFM_det(GEN M, GEN ff)
}
return gerepilecopy(av, mkFF_i(ff, d));
}
+
+GEN
+FFM_FFC_mul(GEN M, GEN C, GEN ff)
+{
+ pari_sp av = avma;
+ ulong pp;
+ GEN P, T, p;
+ _getFF(ff, &T, &p, &pp);
+ M = FFM_to_raw(M);
+ C = FFC_to_raw(C);
+ switch (ff[1])
+ {
+ case t_FF_FpXQ: P = FqM_FqC_mul(M, C, T, p); break;
+ case t_FF_F2xq: P = F2xqM_F2xqC_mul(M, C, T); break;
+ default: P = FlxqM_FlxqC_mul(M, C, T, pp); break;
+ }
+ return gerepilecopy(av, raw_to_FFM(P, ff));
+}
+
+GEN
+FFM_mul(GEN M, GEN N, GEN ff)
+{
+ pari_sp av = avma;
+ ulong pp;
+ GEN P, T, p;
+ _getFF(ff, &T, &p, &pp);
+ M = FFM_to_raw(M);
+ N = FFM_to_raw(N);
+ switch (ff[1])
+ {
+ case t_FF_FpXQ: P = FqM_mul(M, N, T, p); break;
+ case t_FF_F2xq: P = F2xqM_mul(M, N, T); break;
+ default: P = FlxqM_mul(M, N, T, pp); break;
+ }
+ return gerepilecopy(av, raw_to_FFM(P, ff));
+}
diff --git a/src/basemath/RgV.c b/src/basemath/RgV.c
index 6ffa39a..efde7c2 100644
--- a/src/basemath/RgV.c
+++ b/src/basemath/RgV.c
@@ -408,14 +408,22 @@ RgM_RgC_mul_i(GEN x, GEN y, long lx, long l)
for (i=1; i<l; i++) gel(z,i) = RgMrow_RgC_mul_i(x,y,i,lx);
return z;
}
+
GEN
RgM_RgC_mul(GEN x, GEN y)
{
long lx = lg(x);
+ GEN ffx = NULL, ffy = NULL;
if (lx != lg(y)) pari_err_OP("operation 'RgM_RgC_mul'", x,y);
if (lx == 1) return cgetg(1,t_COL);
+ if (RgM_is_FFM(x, &ffx) && RgC_is_FFC(y, &ffy)) {
+ if (!FF_samefield(ffx, ffy))
+ pari_err_OP("*", ffx, ffy);
+ return FFM_FFC_mul(x, y, ffx);
+ }
return RgM_RgC_mul_i(x, y, lx, lgcols(x));
}
+
GEN
RgV_RgM_mul(GEN x, GEN y)
{
@@ -478,12 +486,17 @@ RgM_mul(GEN x, GEN y)
{
pari_sp av = avma;
long j, l, lx, ly = lg(y);
- GEN z;
+ GEN z, ffx = NULL, ffy = NULL;
if (ly == 1) return cgetg(1,t_MAT);
lx = lg(x);
if (lx != lgcols(y)) pari_err_OP("operation 'RgM_mul'", x,y);
if (lx == 1) return zeromatcopy(0, ly-1);
if (is_modular_mul(x,y,&z)) return gerepileupto(av, z);
+ if (RgM_is_FFM(x, &ffx) && RgM_is_FFM(y, &ffy)) {
+ if (!FF_samefield(ffx, ffy))
+ pari_err_OP("*", ffx, ffy);
+ return FFM_mul(x, y, ffx);
+ }
z = cgetg(ly, t_MAT);
l = lgcols(x);
for (j=1; j<ly; j++) gel(z,j) = RgM_RgC_mul_i(x, gel(y,j), lx, l);
diff --git a/src/basemath/alglin1.c b/src/basemath/alglin1.c
index 43ccc91..f2f9b95 100644
--- a/src/basemath/alglin1.c
+++ b/src/basemath/alglin1.c
@@ -360,6 +360,54 @@ gen_Gauss(GEN a, GEN b, void *E, const struct bb_field *ff)
return u;
}
+/* compatible t_MAT * t_COL, lgA = lg(A) = lg(B) > 1, l = lgcols(A) */
+static GEN
+gen_matcolmul_i(GEN A, GEN B, ulong lgA, ulong l,
+ void *E, const struct bb_field *ff)
+{
+ GEN C = cgetg(l, t_COL);
+ ulong i;
+ for (i = 1; i < l; i++) {
+ pari_sp av = avma;
+ GEN e = ff->mul(E, gcoeff(A, i, 1), gel(B, 1));
+ long k;
+ for(k = 2; k < lgA; k++)
+ e = ff->add(E, e, ff->mul(E, gcoeff(A, i, k), gel(B, k)));
+ gel(C, i) = gerepileupto(av, ff->red(E, e));
+ }
+ return C;
+}
+
+GEN
+gen_matcolmul(GEN A, GEN B, void *E, const struct bb_field *ff)
+{
+ ulong lgA = lg(A);
+ if (lgA != lg(B))
+ pari_err_OP("operation 'gen_matcolmul'", A, B);
+ if (lgA == 1)
+ return cgetg(1, t_COL);
+ return gen_matcolmul_i(A, B, lgA, lgcols(A), E, ff);
+}
+
+GEN
+gen_matmul(GEN A, GEN B, void *E, const struct bb_field *ff)
+{
+ ulong j, l, lgA, lgB = lg(B);
+ GEN C;
+ if (lgB == 1)
+ return cgetg(1, t_MAT);
+ lgA = lg(A);
+ if (lgA != lgcols(B))
+ pari_err_OP("operation 'gen_matmul'", A, B);
+ if (lgA == 1)
+ return zeromatcopy(0, lgB - 1);
+ l = lgcols(A);
+ C = cgetg(lgB, t_MAT);
+ for(j = 1; j < lgB; j++)
+ gel(C, j) = gen_matcolmul_i(A, gel(B, j), lgA, l, E, ff);
+ return C;
+}
+
static GEN
image_from_pivot(GEN x, GEN d, long r)
{
@@ -917,6 +965,21 @@ FlxqM_det(GEN a, GEN T, ulong p)
const struct bb_field *S = get_Flxq_field(&E, T, p);
return gen_det(a, E, S);
}
+
+GEN
+FlxqM_FlxqC_mul(GEN A, GEN B, GEN T, unsigned long p) {
+ void *E;
+ const struct bb_field *ff = get_Flxq_field(&E, T, p);
+ return gen_matcolmul(A, B, E, ff);
+}
+
+GEN
+FlxqM_mul(GEN A, GEN B, GEN T, unsigned long p) {
+ void *E;
+ const struct bb_field *ff = get_Flxq_field(&E, T, p);
+ return gen_matmul(A, B, E, ff);
+}
+
GEN
F2xqM_det(GEN a, GEN T)
{
@@ -944,6 +1007,20 @@ F2xqM_inv(GEN a, GEN T)
return gerepilecopy(av, u);
}
+GEN
+F2xqM_F2xqC_mul(GEN A, GEN B, GEN T) {
+ void *E;
+ const struct bb_field *ff = get_F2xq_field(&E, T);
+ return gen_matcolmul(A, B, E, ff);
+}
+
+GEN
+F2xqM_mul(GEN A, GEN B, GEN T) {
+ void *E;
+ const struct bb_field *ff = get_F2xq_field(&E, T);
+ return gen_matmul(A, B, E, ff);
+}
+
static GEN
FqM_gauss_pivot_gen(GEN x, GEN T, GEN p, long *rr)
{
@@ -992,6 +1069,20 @@ FqM_det(GEN x, GEN T, GEN p)
return gen_det(x, E, S);
}
+GEN
+FqM_FqC_mul(GEN A, GEN B, GEN T, GEN p) {
+ void *E;
+ const struct bb_field *ff = get_Fq_field(&E, T, p);
+ return gen_matcolmul(A, B, E, ff);
+}
+
+GEN
+FqM_mul(GEN A, GEN B, GEN T, GEN p) {
+ void *E;
+ const struct bb_field *ff = get_Fq_field(&E, T, p);
+ return gen_matmul(A, B, E, ff);
+}
+
static GEN
FpM_ker_gen(GEN x, GEN p, long deplin)
{
diff --git a/src/headers/paridecl.h b/src/headers/paridecl.h
index 5ffb570..ef9d85d 100644
--- a/src/headers/paridecl.h
+++ b/src/headers/paridecl.h
@@ -954,10 +954,12 @@ GEN F2m_ker(GEN x);
GEN F2m_ker_sp(GEN x, long deplin);
long F2m_rank(GEN x);
GEN F2m_suppl(GEN x);
+GEN F2xqM_F2xqC_mul(GEN a, GEN b, GEN T);
GEN F2xqM_det(GEN a, GEN T);
GEN F2xqM_ker(GEN x, GEN T);
GEN F2xqM_image(GEN x, GEN T);
GEN F2xqM_inv(GEN a, GEN T);
+GEN F2xqM_mul(GEN a, GEN b, GEN T);
long F2xqM_rank(GEN x, GEN T);
GEN Flm_Flc_gauss(GEN a, GEN b, ulong p);
GEN Flm_Flc_invimage(GEN mat, GEN y, ulong p);
@@ -974,11 +976,13 @@ GEN Flm_ker_sp(GEN x, ulong p, long deplin);
long Flm_rank(GEN x, ulong p);
GEN Flm_suppl(GEN x, ulong p);
GEN FlxqM_FlxqC_gauss(GEN a, GEN b, GEN T, ulong p);
+GEN FlxqM_FlxqC_mul(GEN a, GEN b, GEN T, ulong p);
GEN FlxqM_det(GEN a, GEN T, ulong p);
GEN FlxqM_gauss(GEN a, GEN b, GEN T, ulong p);
GEN FlxqM_ker(GEN x, GEN T, ulong p);
GEN FlxqM_image(GEN x, GEN T, ulong p);
GEN FlxqM_inv(GEN x, GEN T, ulong p);
+GEN FlxqM_mul(GEN a, GEN b, GEN T, ulong p);
long FlxqM_rank(GEN x, GEN T, ulong p);
GEN FpM_FpC_gauss(GEN a, GEN b, GEN p);
GEN FpM_FpC_invimage(GEN m, GEN v, GEN p);
@@ -994,12 +998,14 @@ GEN FpM_ker(GEN x, GEN p);
long FpM_rank(GEN x, GEN p);
GEN FpM_suppl(GEN x, GEN p);
GEN FqM_FqC_gauss(GEN a, GEN b, GEN T, GEN p);
+GEN FqM_FqC_mul(GEN a, GEN b, GEN T, GEN p);
GEN FqM_deplin(GEN x, GEN T, GEN p);
GEN FqM_det(GEN x, GEN T, GEN p);
GEN FqM_gauss(GEN a, GEN b, GEN T, GEN p);
GEN FqM_ker(GEN x, GEN T, GEN p);
GEN FqM_image(GEN x, GEN T, GEN p);
GEN FqM_inv(GEN x, GEN T, GEN p);
+GEN FqM_mul(GEN a, GEN b, GEN T, GEN p);
long FqM_rank(GEN a, GEN T, GEN p);
GEN FqM_suppl(GEN x, GEN T, GEN p);
GEN QM_inv(GEN M, GEN dM);
@@ -1035,6 +1041,8 @@ GEN gen_Gauss(GEN a, GEN b, void *E, const struct bb_field *ff);
GEN gen_Gauss_pivot(GEN x, long *rr, void *E, const struct bb_field *ff);
GEN gen_det(GEN a, void *E, const struct bb_field *ff);
GEN gen_ker(GEN x, long deplin, void *E, const struct bb_field *ff);
+GEN gen_matcolmul(GEN a, GEN b, void *E, const struct bb_field *ff);
+GEN gen_matmul(GEN a, GEN b, void *E, const struct bb_field *ff);
GEN image(GEN x);
GEN image2(GEN x);
GEN imagecompl(GEN x);
@@ -2190,10 +2198,12 @@ GEN FF_to_FpXQ(GEN x);
GEN FF_to_FpXQ_i(GEN x);
GEN FF_trace(GEN x);
GEN FF_zero(GEN a);
+GEN FFM_FFC_mul(GEN M, GEN C, GEN ff);
GEN FFM_det(GEN M, GEN ff);
GEN FFM_image(GEN M, GEN ff);
GEN FFM_inv(GEN M, GEN ff);
GEN FFM_ker(GEN M, GEN ff);
+GEN FFM_mul(GEN M, GEN N, GEN ff);
long FFM_rank(GEN M, GEN ff);
GEN FFX_factor(GEN f, GEN x);
GEN FFX_roots(GEN f, GEN x);