Code coverage tests

This page documents the degree to which the PARI/GP source code is tested by our public test suite, distributed with the source distribution in directory src/test/. This is measured by the gcov utility; we then process gcov output using the lcov frond-end.

We test a few variants depending on Configure flags on the pari.math.u-bordeaux.fr machine (x86_64 architecture), and agregate them in the final report:

The target is 90% coverage for all mathematical modules (given that branches depending on DEBUGLEVEL or DEBUGMEM are not covered). This script is run to produce the results below.

LCOV - code coverage report
Current view: top level - basemath - FpV.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 19226-b907b8d) Lines: 611 815 75.0 %
Date: 2016-07-29 07:10:27 Functions: 86 113 76.1 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8             : ANY WARRANTY WHATSOEVER.
       9             : 
      10             : Check the License for details. You should have received a copy of it, along
      11             : with the package; see the file 'COPYING'. If not, write to the Free Software
      12             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13             : 
      14             : #include "pari.h"
      15             : #include "paripriv.h"
      16             : 
      17             : /********************************************************************/
      18             : /**                                                                **/
      19             : /**                           REDUCTION                            **/
      20             : /**                                                                **/
      21             : /********************************************************************/
      22             : /* z in Z^n, return lift(Col(z) * Mod(1,p)) */
      23             : GEN
      24     4057635 : FpC_red(GEN z, GEN p)
      25             : {
      26     4057635 :   long i,l = lg(z);
      27     4057635 :   GEN x = cgetg(l, t_COL);
      28     4057635 :   for (i=1; i<l; i++) gel(x,i) = modii(gel(z,i),p);
      29     4057635 :   return x;
      30             : }
      31             : 
      32             : /* z in Z^n, return lift(Vec(z) * Mod(1,p)) */
      33             : GEN
      34      202006 : FpV_red(GEN z, GEN p)
      35             : {
      36      202006 :   long i,l = lg(z);
      37      202006 :   GEN x = cgetg(l, t_VEC);
      38      202006 :   for (i=1; i<l; i++) gel(x,i) = modii(gel(z,i),p);
      39      202006 :   return x;
      40             : }
      41             : GEN
      42      181874 : FpC_center(GEN z, GEN p, GEN pov2)
      43             : {
      44      181874 :   long i,l = lg(z);
      45      181874 :   GEN x = cgetg(l, t_COL);
      46      181874 :   for (i=1; i<l; i++) gel(x,i) = Fp_center(gel(z,i),p, pov2);
      47      181874 :   return x;
      48             : }
      49             : 
      50             : /* assume 0 <= u < p and ps2 = p>>1 */
      51             : INLINE void
      52          28 : Fp_center_inplace(GEN u, GEN p, GEN ps2)
      53          28 : { if (absi_cmp(u,ps2) > 0) subiiz(u,p,u); }
      54             : 
      55             : void
      56          14 : FpC_center_inplace(GEN z, GEN p, GEN ps2)
      57             : {
      58          14 :   long i,l = lg(z);
      59          42 :   for (i=1; i<l; i++)
      60          28 :     Fp_center_inplace(gel(z,i), p, ps2);
      61          14 : }
      62             : 
      63             : GEN
      64      111678 : Flv_center(GEN z, ulong p, ulong ps2)
      65             : {
      66      111678 :   long i, l = lg(z);
      67      111678 :   GEN x = cgetg(l,t_VECSMALL);
      68      111678 :   for (i=1; i<l; i++) x[i] = Fl_center(z[i],p,ps2);
      69      111678 :   return x;
      70             : }
      71             : 
      72             : /* z in Mat m,n(Z), return lift(z * Mod(1,p)) */
      73             : GEN
      74      764777 : FpM_red(GEN z, GEN p)
      75             : {
      76      764777 :   long i, l = lg(z);
      77      764777 :   GEN x = cgetg(l,t_MAT);
      78      764777 :   for (i=1; i<l; i++) gel(x,i) = FpC_red(gel(z,i), p);
      79      764777 :   return x;
      80             : }
      81             : GEN
      82       49805 : FpM_center(GEN z, GEN p, GEN pov2)
      83             : {
      84       49805 :   long i, l = lg(z);
      85       49805 :   GEN x = cgetg(l,t_MAT);
      86       49805 :   for (i=1; i<l; i++) gel(x,i) = FpC_center(gel(z,i), p, pov2);
      87       49805 :   return x;
      88             : }
      89             : 
      90             : void
      91          14 : FpM_center_inplace(GEN z, GEN p, GEN pov2)
      92             : {
      93          14 :   long i, l = lg(z);
      94          14 :   for (i=1; i<l; i++) FpC_center_inplace(gel(z,i), p, pov2);
      95          14 : }
      96             : GEN
      97        7896 : Flm_center(GEN z, ulong p, ulong ps2)
      98             : {
      99        7896 :   long i, l = lg(z);
     100        7896 :   GEN x = cgetg(l,t_MAT);
     101        7896 :   for (i=1; i<l; i++) gel(x,i) = Flv_center(gel(z,i),p,ps2);
     102        7896 :   return x;
     103             : }
     104             : 
     105             : /********************************************************************/
     106             : /**                                                                **/
     107             : /**                           ADD, SUB                             **/
     108             : /**                                                                **/
     109             : /********************************************************************/
     110             : GEN
     111     7348635 : FpC_add(GEN x, GEN y, GEN p)
     112             : {
     113     7348635 :   long i, lx = lg(x);
     114     7348635 :   GEN z = cgetg(lx, t_COL);
     115     7348635 :   for (i = 1; i < lx; i++) gel(z, i) = Fp_add(gel(x, i), gel(y, i), p);
     116     7348635 :   return z;
     117             : }
     118             : GEN
     119           0 : FpV_add(GEN x, GEN y, GEN p)
     120             : {
     121           0 :   long i, lx = lg(x);
     122           0 :   GEN z = cgetg(lx, t_VEC);
     123           0 :   for (i = 1; i < lx; i++) gel(z, i) = Fp_add(gel(x, i), gel(y, i), p);
     124           0 :   return z;
     125             : }
     126             : GEN
     127      267547 : FpM_add(GEN x, GEN y, GEN p)
     128             : {
     129      267547 :   long lx = lg(x), j;
     130             :   GEN z;
     131      267547 :   if (lx == 1) return cgetg(1, t_MAT);
     132      267547 :   z = cgetg(lx, t_MAT);
     133      267547 :   for (j = 1; j < lx; j++) gel(z,j) = FpC_add(gel(x,j), gel(y,j), p);
     134      267547 :   return z;
     135             : }
     136             : 
     137             : GEN
     138   130418190 : Flv_add(GEN x, GEN y, ulong p)
     139             : {
     140   130418190 :   long i, l = lg(x);
     141   130418190 :   GEN z = cgetg(l, t_VECSMALL);
     142   130418190 :   if (p==2)
     143       29554 :     for (i = 1; i < l; i++) z[i] = x[i]^y[i];
     144             :   else
     145   130388636 :     for (i = 1; i < l; i++) z[i] = Fl_add(x[i], y[i], p);
     146   130418190 :   return z;
     147             : }
     148             : 
     149             : void
     150      460545 : Flv_add_inplace(GEN x, GEN y, ulong p)
     151             : {
     152      460545 :   long i, l = lg(x);
     153      460545 :   if (p==2)
     154      458445 :     for (i = 1; i < l; i++) x[i] ^= y[i];
     155             :   else
     156        2100 :     for (i = 1; i < l; i++) x[i] = Fl_add(x[i], y[i], p);
     157      460545 : }
     158             : 
     159             : ulong
     160        3822 : Flv_sum(GEN x, ulong p)
     161             : {
     162        3822 :   long i, l = lg(x);
     163        3822 :   ulong s = 0;
     164        3822 :   if (p==2)
     165        3822 :     for (i = 1; i < l; i++) s ^= x[i];
     166             :   else
     167           0 :     for (i = 1; i < l; i++) s = Fl_add(s, x[i], p);
     168        3822 :   return s;
     169             : }
     170             : 
     171             : GEN
     172      271726 : FpC_sub(GEN x, GEN y, GEN p)
     173             : {
     174      271726 :   long i, lx = lg(x);
     175      271726 :   GEN z = cgetg(lx, t_COL);
     176      271726 :   for (i = 1; i < lx; i++) gel(z, i) = Fp_sub(gel(x, i), gel(y, i), p);
     177      271726 :   return z;
     178             : }
     179             : GEN
     180           0 : FpV_sub(GEN x, GEN y, GEN p)
     181             : {
     182           0 :   long i, lx = lg(x);
     183           0 :   GEN z = cgetg(lx, t_VEC);
     184           0 :   for (i = 1; i < lx; i++) gel(z, i) = Fp_sub(gel(x, i), gel(y, i), p);
     185           0 :   return z;
     186             : }
     187             : 
     188             : GEN
     189           0 : FpM_sub(GEN x, GEN y, GEN p)
     190             : {
     191           0 :   long i, l = lg(x);
     192           0 :   GEN z = cgetg(l, t_MAT);
     193           0 :   for (i = 1; i < l; i++) gel(z, i) = FpC_sub(gel(x, i), gel(y, i), p);
     194           0 :   return z;
     195             : }
     196             : 
     197             : GEN
     198           0 : Flv_sub(GEN x, GEN y, ulong p)
     199             : {
     200           0 :   long i, l = lg(x);
     201           0 :   GEN z = cgetg(l, t_VECSMALL);
     202           0 :   for (i = 1; i < l; i++) z[i] = Fl_sub(x[i], y[i], p);
     203           0 :   return z;
     204             : }
     205             : 
     206             : void
     207           0 : Flv_sub_inplace(GEN x, GEN y, ulong p)
     208             : {
     209           0 :   long i, l = lg(x);
     210           0 :   for (i = 1; i < l; i++) x[i] = Fl_sub(x[i], y[i], p);
     211           0 : }
     212             : 
     213             : GEN
     214        5866 : Flm_Fl_add(GEN x, ulong y, ulong p)
     215             : {
     216        5866 :   long l = lg(x), i, j;
     217        5866 :   GEN z = cgetg(l,t_MAT);
     218             : 
     219        5866 :   if (l==1) return z;
     220        5866 :   if (l != lgcols(x)) pari_err_OP( "+", x, utoi(y));
     221       24465 :   for (i=1; i<l; i++)
     222             :   {
     223       18599 :     GEN zi = cgetg(l,t_VECSMALL), xi = gel(x,i);
     224       18599 :     gel(z,i) = zi;
     225       18599 :     for (j=1; j<l; j++) zi[j] = xi[j];
     226       18599 :     zi[i] = Fl_add(zi[i], y, p);
     227             :   }
     228        5866 :   return z;
     229             : }
     230             : 
     231             : GEN
     232       13622 : Flm_add(GEN x, GEN y, ulong p)
     233             : {
     234       13622 :   long i, l = lg(x);
     235       13622 :   GEN z = cgetg(l,t_MAT);
     236       13622 :   for (i = 1; i < l; i++) gel(z,i) = Flv_add(gel(x,i),gel(y,i),p);
     237       13622 :   return z;
     238             : }
     239             : 
     240             : GEN
     241           0 : Flm_sub(GEN x, GEN y, ulong p)
     242             : {
     243           0 :   long i, l = lg(x);
     244           0 :   GEN z = cgetg(l, t_MAT);
     245           0 :   for (i = 1; i < l; i++) gel(z, i) = Flv_sub(gel(x, i), gel(y, i), p);
     246           0 :   return z;
     247             : }
     248             : 
     249             : /********************************************************************/
     250             : /**                                                                **/
     251             : /**                           MULTIPLICATION                       **/
     252             : /**                                                                **/
     253             : /********************************************************************/
     254             : GEN
     255      353724 : FpC_Fp_mul(GEN x, GEN y, GEN p)
     256             : {
     257      353724 :   long i, l = lg(x);
     258      353724 :   GEN z = cgetg(l, t_COL);
     259      353724 :   for (i=1;i<l;i++) gel(z,i) = Fp_mul(gel(x,i),y,p);
     260      353724 :   return z;
     261             : }
     262             : GEN
     263      728081 : Flv_Fl_mul(GEN x, ulong y, ulong p)
     264             : {
     265      728081 :   long i, l = lg(x);
     266      728081 :   GEN z = cgetg(l, t_VECSMALL);
     267      728081 :   for (i=1;i<l;i++) z[i] = Fl_mul(x[i], y, p);
     268      728081 :   return z;
     269             : }
     270             : GEN
     271           0 : Flv_Fl_div(GEN x, ulong y, ulong p)
     272             : {
     273           0 :   return Flv_Fl_mul(x, Fl_inv(y, p), p);
     274             : }
     275             : void
     276           0 : Flv_Fl_div_inplace(GEN x, ulong y, ulong p)
     277             : {
     278           0 :   Flv_Fl_mul_inplace(x, Fl_inv(y, p), p);
     279           0 : }
     280             : GEN
     281      267547 : FpM_Fp_mul(GEN X, GEN c, GEN p) {
     282      267547 :   long i, j, h, l = lg(X);
     283      267547 :   GEN A = cgetg(l, t_MAT);
     284      267547 :   if (l == 1) return A;
     285      267547 :   h = lgcols(X);
     286     7488992 :   for (j=1; j<l; j++)
     287             :   {
     288     7221445 :     GEN a = cgetg(h, t_COL), x = gel(X, j);
     289     7221445 :     for (i = 1; i < h; i++) gel(a,i) = Fp_mul(gel(x,i), c, p);
     290     7221445 :     gel(A,j) = a;
     291             :   }
     292      267547 :   return A;
     293             : }
     294             : 
     295             : /* x *= y */
     296             : void
     297      993870 : Flv_Fl_mul_part_inplace(GEN x, ulong y, ulong p, long l)
     298             : {
     299             :   long i;
     300      993870 :   for (i=1;i<=l;i++) x[i] = Fl_mul(x[i], y, p);
     301      993870 : }
     302             : void
     303           0 : Flv_Fl_mul_inplace(GEN x, ulong y, ulong p)
     304           0 : { Flv_Fl_mul_part_inplace(x, y, p, lg(x)-1); }
     305             : 
     306             : /* set y *= x */
     307             : void
     308       27819 : Flm_Fl_mul_inplace(GEN y, ulong x, ulong p)
     309             : {
     310       27819 :   long i, j, m = lgcols(y), l = lg(y);
     311       27819 :   if (HIGHWORD(x | p))
     312      292495 :     for(j=1; j<l; j++)
     313      264676 :       for(i=1; i<m; i++) ucoeff(y,i,j) = Fl_mul(ucoeff(y,i,j), x, p);
     314             :   else
     315           0 :     for(j=1; j<l; j++)
     316           0 :       for(i=1; i<m; i++) ucoeff(y,i,j) = (ucoeff(y,i,j) * x) % p;
     317       27819 : }
     318             : /* return x * y */
     319             : GEN
     320      373303 : Flm_Fl_mul(GEN y, ulong x, ulong p)
     321             : {
     322      373303 :   long i, j, m = lgcols(y), l = lg(y);
     323      373303 :   GEN z = cgetg(l, t_MAT);
     324      373303 :   if (HIGHWORD(x | p))
     325           0 :     for(j=1; j<l; j++) {
     326           0 :       GEN c = cgetg(m, t_VECSMALL); gel(z,j) = c;
     327           0 :       for(i=1; i<m; i++) c[i] = Fl_mul(ucoeff(y,i,j), x, p);
     328             :     }
     329             :   else
     330    10729656 :     for(j=1; j<l; j++) {
     331    10356353 :       GEN c = cgetg(m, t_VECSMALL); gel(z,j) = c;
     332    10356353 :       for(i=1; i<m; i++) c[i] = (ucoeff(y,i,j) * x) % p;
     333             :     }
     334      373303 :   return z;
     335             : }
     336             : 
     337             : GEN
     338          42 : Flv_neg(GEN v, ulong p)
     339             : {
     340          42 :   long i, m = lg(v);
     341          42 :   GEN c = cgetg(m, t_VECSMALL);
     342          42 :   for(i=1; i<m; i++) uel(c,i) = Fl_neg(uel(v,i), p);
     343          42 :   return c;
     344             : }
     345             : 
     346             : void
     347        3347 : Flv_neg_inplace(GEN v, ulong p)
     348             : {
     349             :   long i;
     350      148441 :   for (i = 1; i < lg(v); ++i)
     351      145094 :     v[i] = Fl_neg(v[i], p);
     352        3347 : }
     353             : 
     354             : GEN
     355          14 : Flm_neg(GEN y, ulong p)
     356             : {
     357          14 :   long j, l = lg(y);
     358          14 :   GEN z = cgetg(l, t_MAT);
     359          56 :   for(j=1; j<l; j++)
     360          42 :     gel(z,j) = Flv_neg(gel(y,j), p);
     361          14 :   return z;
     362             : }
     363             : 
     364             : /* x[i,]*y. Assume lx > 1 and 0 < i < lgcols(x) */
     365             : static GEN
     366     9468637 : ZMrow_ZC_mul_i(GEN x, GEN y, long lx, long i)
     367             : {
     368     9468637 :   GEN c = mulii(gcoeff(x,i,1), gel(y,1));
     369             :   long k;
     370   226451536 :   for (k = 2; k < lx; k++)
     371             :   {
     372   216982899 :     GEN t = mulii(gcoeff(x,i,k), gel(y,k));
     373   216982899 :     if (signe(t)) c = addii(c, t);
     374             :   }
     375     9468637 :   return c;
     376             : }
     377             : 
     378             : static long
     379    21635677 : zmrow_zc_mul(GEN x, GEN y, long lx, long i)
     380             : {
     381             :   long k;
     382    21635677 :   long c = coeff(x,i,1) * y[1];
     383   289228198 :   for (k = 2; k < lx; k++)
     384   267592521 :     c += coeff(x,i,k) * y[k];
     385    21635677 :   return c;
     386             : }
     387             : 
     388             : GEN
     389     1653239 : zm_zc_mul(GEN x, GEN y)
     390             : {
     391     1653239 :   long lx = lg(x), l, i;
     392             :   GEN z;
     393     1653239 :   if (lx == 1) return cgetg(1, t_VECSMALL);
     394     1653239 :   l = lg(gel(x,1));
     395     1653239 :   z = cgetg(l,t_VECSMALL);
     396     1653239 :   for (i=1; i<l; i++) z[i] = zmrow_zc_mul(x, y, lx, i);
     397     1653239 :   return z;
     398             : }
     399             : 
     400             : GEN
     401         252 : zm_mul(GEN x, GEN y)
     402             : {
     403         252 :   long i,j,lx=lg(x), ly=lg(y);
     404             :   GEN z;
     405         252 :   if (ly==1) return cgetg(1,t_MAT);
     406         252 :   z = cgetg(ly,t_MAT);
     407         252 :   if (lx==1)
     408             :   {
     409           0 :     for (i=1; i<ly; i++) gel(z,i) = cgetg(1,t_VECSMALL);
     410           0 :     return z;
     411             :   }
     412        2128 :   for (j=1; j<ly; j++)
     413        1876 :     gel(z,j) = zm_zc_mul(x, gel(y,j));
     414         252 :   return z;
     415             : }
     416             : 
     417             : static ulong
     418   178229416 : Flmrow_Flc_mul_SMALL(GEN x, GEN y, ulong p, long lx, long i)
     419             : {
     420   178229416 :   ulong c = ucoeff(x,i,1) * uel(y,1);
     421             :   long k;
     422  5820372688 :   for (k = 2; k < lx; k++) {
     423  5642143272 :     c += ucoeff(x,i,k) * uel(y,k);
     424  5642143272 :     if (c & HIGHBIT) c %= p;
     425             :   }
     426   178229416 :   return c % p;
     427             : }
     428             : 
     429             : static ulong
     430    44020131 : Flmrow_Flc_mul_i(GEN x, GEN y, ulong p, ulong pi, long lx, long i)
     431             : {
     432             :   ulong l0, l1, v1, h0, h1;
     433    44020131 :   long k = 1;
     434             :   LOCAL_OVERFLOW;
     435             :   LOCAL_HIREMAINDER;
     436    44020131 :   l1 = mulll(ucoeff(x,i,k), uel(y,k)); h1 = hiremainder; v1 = 0;
     437   706900926 :   while (++k < lx) {
     438   618860664 :     l0 = mulll(ucoeff(x,i,k), uel(y,k)); h0 = hiremainder;
     439   618860664 :     l1 = addll(l0, l1); h1 = addllx(h0, h1); v1 += overflow;
     440             :   }
     441    44020131 :   if (v1 == 0) return remll_pre(h1, l1, p, pi);
     442     1475646 :   else return remlll_pre(v1, h1, l1, p, pi);
     443             : }
     444             : 
     445             : static GEN
     446      396016 : Flm_Flc_mul_i_2(GEN x, GEN y, long lx, long l)
     447             : {
     448             :   long i,j;
     449      396016 :   GEN z = NULL;
     450             : 
     451     4339982 :   for (j=1; j<lx; j++)
     452             :   {
     453     3943966 :     if (!y[j]) continue;
     454      425965 :     if (!z) z = Flv_copy(gel(x,j));
     455      158179 :     else for (i = 1; i < l; i++) z[i] ^= coeff(x,i,j);
     456             :   }
     457      396016 :   if (!z) z = zero_zv(l-1);
     458      396016 :   return z;
     459             : }
     460             : 
     461             : static GEN
     462      446469 : FpM_FpC_mul_i(GEN x, GEN y, long lx, long l, GEN p)
     463             : {
     464      446469 :   GEN z = cgetg(l,t_COL);
     465             :   long i;
     466     8943843 :   for (i = 1; i < l; i++)
     467             :   {
     468     8497374 :     pari_sp av = avma;
     469     8497374 :     GEN c = ZMrow_ZC_mul_i(x, y, lx, i);
     470     8497374 :     gel(z,i) = gerepileuptoint(av, modii(c,p));
     471             :   }
     472      446469 :   return z;
     473             : }
     474             : static GEN
     475     7773734 : Flm_Flc_mul_i_SMALL(GEN x, GEN y, long lx, long l, ulong p)
     476             : {
     477     7773734 :   GEN z = cgetg(l,t_VECSMALL);
     478             :   long i;
     479     7789596 :   for (i = 1; i < l; i++) z[i] = Flmrow_Flc_mul_SMALL(x, y, p, lx, i);
     480     7778756 :   return z;
     481             : }
     482             : static GEN
     483     3908036 : Flm_Flc_mul_i(GEN x, GEN y, long lx, long l, ulong p, ulong pi)
     484             : {
     485     3908036 :   GEN z = cgetg(l,t_VECSMALL);
     486             :   long i;
     487     3915783 :   for (i = 1; i < l; i++) z[i] = Flmrow_Flc_mul_i(x, y, p, pi, lx, i);
     488     3914802 :   return z;
     489             : }
     490             : INLINE GEN
     491      811927 : F2m_F2c_mul_i(GEN x, GEN y, long lx, long l)
     492             : {
     493             :   long j;
     494      811927 :   GEN z = NULL;
     495             : 
     496     5944234 :   for (j=1; j<lx; j++)
     497             :   {
     498     5132307 :     if (!F2v_coeff(y,j)) continue;
     499     1409582 :     if (!z) z = vecsmall_copy(gel(x,j));
     500      618205 :     else F2v_add_inplace(z,gel(x,j));
     501             :   }
     502      811927 :   if (!z) z = zero_F2v(l);
     503      811927 :   return z;
     504             : }
     505             : 
     506             : GEN
     507      415070 : FpM_mul(GEN x, GEN y, GEN p)
     508             : {
     509      415070 :   long j, l, lx=lg(x), ly=lg(y);
     510             :   GEN z;
     511      415070 :   if (ly==1) return cgetg(1,t_MAT);
     512      415070 :   if (lx==1) return zeromat(0, ly-1);
     513      415070 :   if (lgefint(p) == 3)
     514             :   {
     515      414345 :     pari_sp av = avma;
     516      414345 :     ulong pp = uel(p,2);
     517      414345 :     if (pp == 2)
     518             :     {
     519      145188 :       x = ZM_to_F2m(x);
     520      145188 :       y = ZM_to_F2m(y);
     521      145188 :       z = F2m_to_ZM(F2m_mul(x,y));
     522             :     }
     523             :     else
     524             :     {
     525      269157 :       x = ZM_to_Flm(x, pp);
     526      269157 :       y = ZM_to_Flm(y, pp);
     527      269157 :       z = Flm_to_ZM(Flm_mul(x,y, pp));
     528             :     }
     529      414345 :     return gerepileupto(av, z);
     530             :   }
     531         725 :   l = lgcols(x); z = cgetg(ly,t_MAT);
     532         725 :   for (j=1; j<ly; j++) gel(z,j) = FpM_FpC_mul_i(x, gel(y,j), lx, l, p);
     533         725 :   return z;
     534             : }
     535             : GEN
     536      356527 : Flm_mul(GEN x, GEN y, ulong p)
     537             : {
     538      356527 :   long i,j,l,lx=lg(x), ly=lg(y);
     539             :   GEN z;
     540      356527 :   if (ly==1) return cgetg(1,t_MAT);
     541      356527 :   z = cgetg(ly,t_MAT);
     542      356527 :   if (lx==1)
     543             :   {
     544           0 :     for (i=1; i<ly; i++) gel(z,i) = cgetg(1,t_VECSMALL);
     545           0 :     return z;
     546             :   }
     547      356527 :   l = lgcols(x);
     548      356527 :   if (SMALL_ULONG(p)) {
     549     4968429 :     for (j=1; j<ly; j++)
     550     4618706 :       gel(z,j) = Flm_Flc_mul_i_SMALL(x, gel(y,j), lx, l, p);
     551             :   } else {
     552        6804 :     ulong pi = get_Fl_red(p);
     553       67562 :     for (j=1; j<ly; j++)
     554       60758 :       gel(z,j) = Flm_Flc_mul_i(x, gel(y,j), lx, l, p, pi);
     555             :   }
     556      356527 :   return z;
     557             : }
     558             : GEN
     559      145237 : F2m_mul(GEN x, GEN y)
     560             : {
     561      145237 :   long i,j,l,lx=lg(x), ly=lg(y);
     562             :   GEN z;
     563      145237 :   if (ly==1) return cgetg(1,t_MAT);
     564      145237 :   z = cgetg(ly,t_MAT);
     565      145237 :   if (lx==1)
     566             :   {
     567           0 :     for (i=1; i<ly; i++) gel(z,i) = mkvecsmall(0);
     568           0 :     return z;
     569             :   }
     570      145237 :   l = coeff(x,1,1);
     571      145237 :   for (j=1; j<ly; j++) gel(z,j) = F2m_F2c_mul_i(x, gel(y,j), lx, l);
     572      145237 :   return z;
     573             : }
     574             : 
     575             : static GEN
     576       14203 : _Flm_mul(void *p , GEN x, GEN y)
     577       14203 : { return Flm_mul(x,y,*(ulong*)p); }
     578             : static GEN
     579       42336 : _Flm_sqr(void *p, GEN x)
     580       42336 : { return Flm_mul(x,x,*(ulong*)p); }
     581             : GEN
     582       19628 : Flm_powu(GEN x, ulong n, ulong p)
     583             : {
     584       19628 :   pari_sp av = avma;
     585       19628 :   if (!n) return matid(lg(x)-1);
     586       19628 :   return gerepileupto(av, gen_powu(x, n, (void*)&p, &_Flm_sqr, &_Flm_mul));
     587             : }
     588             : static GEN
     589           0 : _F2m_mul(void *data, GEN x, GEN y)
     590           0 : { (void) data; return F2m_mul(x,y); }
     591             : static GEN
     592           0 : _F2m_sqr(void *data, GEN x)
     593           0 : { (void) data; return F2m_mul(x,x); }
     594             : GEN
     595           0 : F2m_powu(GEN x, ulong n)
     596             : {
     597           0 :   pari_sp av = avma;
     598           0 :   if (!n) return matid(lg(x)-1);
     599           0 :   return gerepileupto(av, gen_powu(x, n,NULL, &_F2m_sqr, &_F2m_mul));
     600             : }
     601             : static GEN
     602           0 : _FpM_mul(void *p , GEN x, GEN y)
     603           0 : { return FpM_mul(x,y,(GEN)p); }
     604             : static GEN
     605           0 : _FpM_sqr(void *p, GEN x)
     606           0 : { return FpM_mul(x,x,(GEN)p); }
     607             : GEN
     608           0 : FpM_powu(GEN x, ulong n, GEN p)
     609             : {
     610           0 :   pari_sp av = avma;
     611           0 :   if (!n) return matid(lg(x)-1);
     612           0 :   if (lgefint(p) == 3)
     613             :   {
     614           0 :     pari_sp av = avma;
     615           0 :     ulong pp = uel(p,2);
     616             :     GEN z;
     617           0 :     if (pp == 2)
     618           0 :       z = F2m_to_ZM(F2m_powu(ZM_to_F2m(x),n));
     619             :     else
     620           0 :       z = Flm_to_ZM(Flm_powu(ZM_to_Flm(x, pp), n, pp));
     621           0 :     return gerepileupto(av, z);
     622             :   }
     623           0 :   return gerepileupto(av, gen_powu(x, n, (void*)p, &_FpM_sqr, &_FpM_mul));
     624             : }
     625             : 
     626             : /*Multiple a column vector by a line vector to make a matrix*/
     627             : GEN
     628        1043 : FpC_FpV_mul(GEN x, GEN y, GEN p)
     629             : {
     630        1043 :   long i,j, lx=lg(x), ly=lg(y);
     631             :   GEN z;
     632        1043 :   if (ly==1) return cgetg(1,t_MAT);
     633        1043 :   z = cgetg(ly,t_MAT);
     634       23037 :   for (j=1; j < ly; j++)
     635             :   {
     636       21994 :     gel(z,j) = cgetg(lx,t_COL);
     637       21994 :     for (i=1; i<lx; i++) gcoeff(z,i,j) = Fp_mul(gel(x,i),gel(y,j), p);
     638             :   }
     639        1043 :   return z;
     640             : }
     641             : 
     642             : /* Multiply a line vector by a column and return a scalar (t_INT) */
     643             : GEN
     644     2262502 : FpV_dotproduct(GEN x, GEN y, GEN p)
     645             : {
     646     2262502 :   long i, lx = lg(x);
     647             :   pari_sp av;
     648             :   GEN c;
     649     2262502 :   if (lx == 1) return gen_0;
     650     2262467 :   av = avma; c = mulii(gel(x,1),gel(y,1));
     651     2262467 :   for (i=2; i<lx; i++) c = addii(c, mulii(gel(x,i),gel(y,i)));
     652     2262467 :   return gerepileuptoint(av, modii(c,p));
     653             : }
     654             : GEN
     655           0 : FpV_dotsquare(GEN x, GEN p)
     656             : {
     657           0 :   long i, lx = lg(x);
     658             :   pari_sp av;
     659             :   GEN c;
     660           0 :   if (lx == 1) return gen_0;
     661           0 :   av = avma; c = sqri(gel(x,1));
     662           0 :   for (i=2; i<lx; i++) c = addii(c, sqri(gel(x,i)));
     663           0 :   return gerepileuptoint(av, modii(c,p));
     664             : }
     665             : 
     666             : INLINE ulong
     667      812692 : Flv_dotproductspec_SMALL(GEN x, GEN y, ulong p, long lx)
     668             : {
     669      812692 :   ulong c = uel(x,0) * uel(y,0);
     670             :   long k;
     671    12975916 :   for (k = 1; k < lx; k++) {
     672    12163224 :     c += uel(x,k) * uel(y,k);
     673    12163224 :     if (c & HIGHBIT) c %= p;
     674             :   }
     675      812692 :   return c % p;
     676             : }
     677             : 
     678             : INLINE ulong
     679        4790 : Flv_dotproductspec_i(GEN x, GEN y, ulong p, ulong pi, long lx)
     680             : {
     681             :   ulong l0, l1, v1, h0, h1;
     682        4790 :   long i = 0;
     683             :   LOCAL_OVERFLOW;
     684             :   LOCAL_HIREMAINDER;
     685        4790 :   l1 = mulll(uel(x,i), uel(y,i)); h1 = hiremainder; v1 = 0;
     686      147856 :   while (++i < lx) {
     687      138276 :     l0 = mulll(uel(x,i), uel(y,i)); h0 = hiremainder;
     688      138276 :     l1 = addll(l0, l1); h1 = addllx(h0, h1); v1 += overflow;
     689             :   }
     690        4790 :   if (v1 == 0) return remll_pre(h1, l1, p, pi);
     691           0 :   else return remlll_pre(v1, h1, l1, p, pi);
     692             : }
     693             : 
     694             : ulong
     695      373303 : Flv_dotproduct(GEN x, GEN y, ulong p)
     696             : {
     697      373303 :   long lx = lg(x)-1;
     698      373303 :   if (lx == 1) return 0;
     699      373303 :   if (SMALL_ULONG(p))
     700      373303 :     return Flv_dotproductspec_SMALL(x+1, y+1, p, lx);
     701             :   else
     702           0 :     return Flv_dotproductspec_i(x+1, y+1, p, get_Fl_red(p), lx);
     703             : }
     704             : 
     705             : ulong
     706       44458 : Flv_dotproduct_pre(GEN x, GEN y, ulong p, ulong pi)
     707             : {
     708       44458 :   long lx = lg(x)-1;
     709       44458 :   if (lx == 1) return 0;
     710       44458 :   if (SMALL_ULONG(p))
     711       39752 :     return Flv_dotproductspec_SMALL(x+1, y+1, p, lx);
     712             :   else
     713        4706 :     return Flv_dotproductspec_i(x+1, y+1, p, pi, lx);
     714             : }
     715             : 
     716             : ulong
     717      408052 : Flx_dotproduct(GEN x, GEN y, ulong p)
     718             : {
     719      408052 :   long lx = minss(lgpol(x), lgpol(y));
     720      408052 :   if (lx == 0) return 0;
     721      399722 :   if (SMALL_ULONG(p))
     722      399638 :     return Flv_dotproductspec_SMALL(x+2, y+2, p, lx);
     723             :   else
     724          84 :     return Flv_dotproductspec_i(x+2, y+2, p, get_Fl_red(p), lx);
     725             : }
     726             : 
     727             : ulong
     728           0 : F2v_dotproduct(GEN x, GEN y)
     729             : {
     730           0 :   long i, lx = lg(x);
     731             :   ulong c;
     732           0 :   if (lx == 2) return 0;
     733           0 :   c = uel(x,2) & uel(y,2);
     734           0 :   for (i=3; i<lx; i++) c ^= uel(x,i) & uel(y,i);
     735             : #ifdef LONG_IS_64BIT
     736           0 :   c ^= c >> 32;
     737             : #endif
     738           0 :   c ^= c >> 16;
     739           0 :   c ^= c >>  8;
     740           0 :   c ^= c >>  4;
     741           0 :   c ^= c >>  2;
     742           0 :   c ^= c >>  1;
     743           0 :   return c & 1;
     744             : }
     745             : 
     746             : GEN
     747      439777 : FpM_FpC_mul(GEN x, GEN y, GEN p)
     748             : {
     749      439777 :   long lx = lg(x);
     750      439777 :   return lx==1? cgetg(1,t_COL): FpM_FpC_mul_i(x, y, lx, lgcols(x), p);
     751             : }
     752             : GEN
     753      743860 : Flm_Flc_mul(GEN x, GEN y, ulong p)
     754             : {
     755      743860 :   long l, lx = lg(x);
     756      743860 :   if (lx==1) return cgetg(1,t_VECSMALL);
     757      638699 :   l = lgcols(x);
     758      638699 :   if (p==2)
     759      396016 :     return Flm_Flc_mul_i_2(x, y, lx, l);
     760      242683 :   else if (SMALL_ULONG(p))
     761      242571 :     return Flm_Flc_mul_i_SMALL(x, y, lx, l, p);
     762             :   else
     763         112 :     return Flm_Flc_mul_i(x, y, lx, l, p, get_Fl_red(p));
     764             : }
     765             : 
     766             : GEN
     767     6760555 : Flm_Flc_mul_pre(GEN x, GEN y, ulong p, ulong pi)
     768             : {
     769     6760555 :   long l, lx = lg(x);
     770     6760555 :   if (lx==1) return cgetg(1,t_VECSMALL);
     771     6760555 :   l = lgcols(x);
     772     6763568 :   if (SMALL_ULONG(p))
     773     2913957 :     return Flm_Flc_mul_i_SMALL(x, y, lx, l, p);
     774             :   else
     775     3849611 :     return Flm_Flc_mul_i(x, y, lx, l, p, pi);
     776             : }
     777             : 
     778             : GEN
     779           0 : F2m_F2c_mul(GEN x, GEN y)
     780             : {
     781           0 :   long l, lx = lg(x);
     782           0 :   if (lx==1) return cgetg(1,t_VECSMALL);
     783           0 :   l = coeff(x,1,1);
     784           0 :   return F2m_F2c_mul_i(x, y, lx, l);
     785             : }
     786             : /* RgV_to_RgX(FpM_FpC_mul(x,y,p), v), p != NULL, memory clean */
     787             : GEN
     788      246444 : FpM_FpC_mul_FpX(GEN x, GEN y, GEN p, long v)
     789             : {
     790      246444 :   long i, l, lx = lg(x);
     791             :   GEN z;
     792      246444 :   if (lx==1) return pol_0(v);
     793      246444 :   l = lgcols(x);
     794      246444 :   z = new_chunk(l+1);
     795      396335 :   for (i=l-1; i; i--)
     796             :   {
     797      380501 :     pari_sp av = avma;
     798      380501 :     GEN p1 = ZMrow_ZC_mul_i(x,y,lx,i);
     799      380501 :     p1 = modii(p1, p);
     800      380501 :     if (signe(p1))
     801             :     {
     802      230610 :       if (i != l-1) stackdummy((pari_sp)(z + l+1), (pari_sp)(z + i+2));
     803      230610 :       gel(z,i+1) = gerepileuptoint(av, p1);
     804      230610 :       break;
     805             :     }
     806      149891 :     avma = av;
     807             :   }
     808      246444 :   if (!i) { avma = (pari_sp)(z + l+1); return pol_0(v); }
     809      230610 :   z[0] = evaltyp(t_POL) | evallg(i+2);
     810      230610 :   z[1] = evalsigne(1) | evalvarn(v);
     811      821372 :   for (; i; i--)
     812             :   {
     813      590762 :     pari_sp av = avma;
     814      590762 :     GEN p1 = ZMrow_ZC_mul_i(x,y,lx,i);
     815      590762 :     gel(z,i+1) = gerepileuptoint(av, modii(p1,p));
     816             :   }
     817      230610 :   return z;
     818             : }
     819             : 
     820             : /********************************************************************/
     821             : /**                                                                **/
     822             : /**                           TRANSPOSITION                        **/
     823             : /**                                                                **/
     824             : /********************************************************************/
     825             : 
     826             : /* == zm_transpose */
     827             : GEN
     828        3952 : Flm_transpose(GEN x)
     829             : {
     830        3952 :   long i, dx, lx = lg(x);
     831             :   GEN y;
     832        3952 :   if (lx == 1) return cgetg(1,t_MAT);
     833        3952 :   dx = lgcols(x); y = cgetg(dx,t_MAT);
     834        3954 :   for (i=1; i<dx; i++) gel(y,i) = Flm_row(x,i);
     835        3951 :   return y;
     836             : }
     837             : 
     838             : /********************************************************************/
     839             : /**                                                                **/
     840             : /**                           SCALAR MATRICES                      **/
     841             : /**                                                                **/
     842             : /********************************************************************/
     843             : 
     844             : GEN
     845          28 : gen_matid(long n, void *E, const struct bb_field *S)
     846             : {
     847          28 :   GEN y = cgetg(n+1,t_MAT), _0, _1;
     848             :   long i;
     849          28 :   if (n < 0) pari_err_DOMAIN("gen_matid", "dimension","<",gen_0,stoi(n));
     850          28 :   _0 = S->s(E,0);
     851          28 :   _1 = S->s(E,1);
     852         112 :   for (i=1; i<=n; i++)
     853             :   {
     854          84 :     GEN z = const_col(n, _0); gel(z,i) = _1;
     855          84 :     gel(y, i) = z;
     856             :   }
     857          28 :   return y;
     858             : }
     859             : 
     860             : GEN
     861       12864 : matid_F2m(long n)
     862             : {
     863       12864 :   GEN y = cgetg(n+1,t_MAT);
     864             :   long i;
     865       12864 :   if (n < 0) pari_err_DOMAIN("matid_F2m", "dimension","<",gen_0,stoi(n));
     866       12864 :   for (i=1; i<=n; i++) { gel(y,i) = zero_F2v(n); F2v_set(gel(y,i),i); }
     867       12864 :   return y;
     868             : }
     869             : 
     870             : GEN
     871          14 : matid_F2xqM(long n, GEN T)
     872             : {
     873             :   void *E;
     874          14 :   const struct bb_field *S = get_F2xq_field(&E, T);
     875          14 :   return gen_matid(n, E, S);
     876             : }
     877             : GEN
     878          14 : matid_FlxqM(long n, GEN T, ulong p)
     879             : {
     880             :   void *E;
     881          14 :   const struct bb_field *S = get_Flxq_field(&E, T, p);
     882          14 :   return gen_matid(n, E, S);
     883             : }
     884             : 
     885             : GEN
     886      250660 : matid_Flm(long n)
     887             : {
     888      250660 :   GEN y = cgetg(n+1,t_MAT);
     889             :   long i;
     890      250660 :   if (n < 0) pari_err_DOMAIN("matid_Flm", "dimension","<",gen_0,stoi(n));
     891      250660 :   for (i=1; i<=n; i++) { gel(y,i) = zero_zv(n); ucoeff(y, i,i) = 1; }
     892      250660 :   return y;
     893             : }
     894             : 
     895             : GEN
     896          42 : scalar_Flm(long s, long n)
     897             : {
     898             :   long i;
     899          42 :   GEN y = cgetg(n+1,t_MAT);
     900          42 :   for (i=1; i<=n; i++) { gel(y,i) = zero_Flv(n); coeff(y, i,i) = s; }
     901          42 :   return y;
     902             : }
     903             : 
     904             : /********************************************************************/
     905             : /**                                                                **/
     906             : /**                           CONVERSIONS                          **/
     907             : /**                                                                **/
     908             : /********************************************************************/
     909             : GEN
     910    15471759 : ZV_to_Flv(GEN x, ulong p)
     911             : {
     912    15471759 :   long i, n = lg(x);
     913    15471759 :   GEN y = cgetg(n,t_VECSMALL);
     914    15471854 :   for (i=1; i<n; i++) y[i] = umodiu(gel(x,i), p);
     915    15471748 :   return y;
     916             : }
     917             : GEN
     918     1541365 : ZM_to_Flm(GEN x, ulong p)
     919             : {
     920     1541365 :   long j,n = lg(x);
     921     1541365 :   GEN y = cgetg(n,t_MAT);
     922     1541376 :   if (n == 1) return y;
     923     1541250 :   for (j=1; j<n; j++) gel(y,j) = ZV_to_Flv(gel(x,j), p);
     924     1541240 :   return y;
     925             : }
     926             : GEN
     927        1365 : ZMV_to_FlmV(GEN z, ulong m)
     928             : {
     929        1365 :   long i, l = lg(z);
     930        1365 :   GEN x = cgetg(l,t_VEC);
     931        1365 :   for (i=1; i<l; i++) gel(x,i) = ZM_to_Flm(gel(z,i), m);
     932        1365 :   return x;
     933             : }
     934             : 
     935             : /*                          TO INTMOD                        */
     936             : static GEN
     937     1269507 : to_intmod(GEN x, GEN p) { retmkintmod(modii(x, p), p); }
     938             : static GEN
     939        1154 : Fl_to_intmod(ulong x, GEN p) { retmkintmod(utoi(x), p); }
     940             : 
     941             : GEN
     942      281628 : Fp_to_mod(GEN z, GEN p)
     943             : {
     944      281628 :   retmkintmod(modii(z, p), icopy(p));
     945             : }
     946             : 
     947             : /* z in Z[X], return z * Mod(1,p), normalized*/
     948             : GEN
     949       53508 : FpX_to_mod(GEN z, GEN p)
     950             : {
     951       53508 :   long i,l = lg(z);
     952       53508 :   GEN x = cgetg(l,t_POL);
     953       53508 :   if (l >2) p = icopy(p);
     954       53508 :   for (i=2; i<l; i++) gel(x,i) = to_intmod(gel(z,i), p);
     955       53508 :   x[1] = z[1]; return normalizepol_lg(x,l);
     956             : }
     957             : 
     958             : /* z in Z^n, return z * Mod(1,p), normalized*/
     959             : GEN
     960        4032 : FpV_to_mod(GEN z, GEN p)
     961             : {
     962        4032 :   long i,l = lg(z);
     963        4032 :   GEN x = cgetg(l, t_VEC);
     964        4032 :   if (l == 1) return x;
     965        4032 :   p = icopy(p);
     966        4032 :   for (i=1; i<l; i++) gel(x,i) = to_intmod(gel(z,i), p);
     967        4032 :   return x;
     968             : }
     969             : /* z in Z^n, return z * Mod(1,p), normalized*/
     970             : GEN
     971          63 : FpC_to_mod(GEN z, GEN p)
     972             : {
     973          63 :   long i, l = lg(z);
     974          63 :   GEN x = cgetg(l, t_COL);
     975          63 :   if (l == 1) return x;
     976          63 :   p = icopy(p);
     977          63 :   for (i=1; i<l; i++) gel(x,i) = to_intmod(gel(z,i), p);
     978          63 :   return x;
     979             : }
     980             : /* z in Mat m,n(Z), return z * Mod(1,p), normalized*/
     981             : GEN
     982         100 : FpM_to_mod(GEN z, GEN p)
     983             : {
     984         100 :   long i, j, m, l = lg(z);
     985         100 :   GEN  x = cgetg(l,t_MAT), y, zi;
     986         100 :   if (l == 1) return x;
     987          93 :   m = lgcols(z);
     988          93 :   p = icopy(p);
     989         790 :   for (i=1; i<l; i++)
     990             :   {
     991         697 :     gel(x,i) = cgetg(m,t_COL);
     992         697 :     y = gel(x,i); zi= gel(z,i);
     993         697 :     for (j=1; j<m; j++) gel(y,j) = to_intmod(gel(zi,j), p);
     994             :   }
     995          93 :   return x;
     996             : }
     997             : GEN
     998          28 : Flc_to_mod(GEN z, ulong pp)
     999             : {
    1000          28 :   long i, l = lg(z);
    1001          28 :   GEN p, x = cgetg(l, t_COL);
    1002          28 :   if (l == 1) return x;
    1003          28 :   p = utoipos(pp);
    1004          28 :   for (i=1; i<l; i++) gel(x,i) = Fl_to_intmod(z[i], p);
    1005          28 :   return x;
    1006             : }
    1007             : GEN
    1008         166 : Flm_to_mod(GEN z, ulong pp)
    1009             : {
    1010         166 :   long i, j, m, l = lg(z);
    1011         166 :   GEN p, x = cgetg(l,t_MAT), y, zi;
    1012         166 :   if (l == 1) return x;
    1013         152 :   m = lgcols(z);
    1014         152 :   p = utoipos(pp);
    1015         547 :   for (i=1; i<l; i++)
    1016             :   {
    1017         395 :     gel(x,i) = cgetg(m,t_COL);
    1018         395 :     y = gel(x,i); zi= gel(z,i);
    1019         395 :     for (j=1; j<m; j++) gel(y,j) = Fl_to_intmod(zi[j], p);
    1020             :   }
    1021         152 :   return x;
    1022             : }
    1023             : 
    1024             : GEN
    1025        1820 : FpVV_to_mod(GEN z, GEN p)
    1026             : {
    1027        1820 :   long i, j, m, l = lg(z);
    1028        1820 :   GEN  x = cgetg(l,t_VEC), y, zi;
    1029        1820 :   if (l == 1) return x;
    1030        1743 :   m = lgcols(z);
    1031        1743 :   p = icopy(p);
    1032        3647 :   for (i=1; i<l; i++)
    1033             :   {
    1034        1904 :     gel(x,i) = cgetg(m,t_VEC);
    1035        1904 :     y = gel(x,i); zi= gel(z,i);
    1036        1904 :     for (j=1; j<m; j++) gel(y,j) = to_intmod(gel(zi,j), p);
    1037             :   }
    1038        1743 :   return x;
    1039             : }
    1040             : 
    1041             : /* z in Z^n, return z * Mod(1,p), normalized*/
    1042             : GEN
    1043           7 : FpXQC_to_mod(GEN z, GEN T, GEN p)
    1044             : {
    1045           7 :   long i,l = lg(z);
    1046           7 :   GEN x = cgetg(l, t_COL); T = FpX_to_mod(T, p);
    1047          21 :   for (i=1; i<l; i++)
    1048          14 :     gel(x,i) = mkpolmod(FpX_to_mod(gel(z,i), p), T);
    1049           7 :   return x;
    1050             : }
    1051             : 
    1052             : /********************************************************************/
    1053             : /*                                                                  */
    1054             : /*                     Blackbox linear algebra                      */
    1055             : /*                                                                  */
    1056             : /********************************************************************/
    1057             : 
    1058             : /* A sparse column (zCs) v is a t_COL with two components C and E which are
    1059             :  * t_VECSMALL of the same length, representing sum_i E[i]*e_{C[i]}, where
    1060             :  * (e_j) is the canonical basis.
    1061             :  * A sparse matrix (zMs) is a t_VEC of zCs */
    1062             : 
    1063             : /* FpCs and FpMs are identical but E[i] is interpreted as a _signed_ C long
    1064             :  * integer representing an element of Fp. This is important since p can be
    1065             :  * large and p+E[i] would not fit in a C long.  */
    1066             : 
    1067             : /* RgCs and RgMs are similar, except that the type of the component is
    1068             :  * unspecified. Functions handling RgCs/RgMs must be independent of the type
    1069             :  * of E. */
    1070             : 
    1071             : /* Most functions take an argument nbrow which is the number of lines of the
    1072             :  * column/matrix, which cannot be derived from the data. */
    1073             : 
    1074             : GEN
    1075           0 : zCs_to_ZC(GEN R, long nbrow)
    1076             : {
    1077           0 :   GEN C = gel(R,1), E = gel(R,2), c = zerocol(nbrow);
    1078           0 :   long j, l = lg(C);
    1079           0 :   for (j = 1; j < l; ++j) gel(c, C[j]) = stoi(E[j]);
    1080           0 :   return c;
    1081             : }
    1082             : 
    1083             : GEN
    1084           0 : zMs_to_ZM(GEN M, long nbrow)
    1085             : {
    1086           0 :   long i, l = lg(M);
    1087           0 :   GEN m = cgetg(l, t_MAT);
    1088           0 :   for (i = 1; i < l; ++i) gel(m,i) = zCs_to_ZC(gel(M,i), nbrow);
    1089           0 :   return m;
    1090             : }
    1091             : 
    1092             : /* Solve equation f(X) = B (mod p) where B is a FpV, and f is an endomorphism.
    1093             :  * Return either a solution as a t_COL, or a kernel vector as a t_VEC. */
    1094             : GEN
    1095         112 : gen_FpM_Wiedemann(void *E, GEN (*f)(void*, GEN), GEN B, GEN p)
    1096             : {
    1097         112 :   pari_sp ltop = avma;
    1098         112 :   long col = 0, n = lg(B)-1, m = 2*n+1;
    1099         112 :   if (ZV_equal0(B)) return zerocol(n);
    1100         224 :   while (++col <= n)
    1101             :   {
    1102         112 :     pari_sp btop = avma, av;
    1103             :     long i, lQ;
    1104         112 :     GEN V, Q, M, W = B;
    1105         112 :     GEN b = cgetg(m+2, t_POL);
    1106         112 :     b[1] = evalsigne(1)|evalvarn(0);
    1107         112 :     gel(b, 2) = gel(W, col);
    1108       49826 :     for (i = 3; i<m+2; i++)
    1109       49714 :       gel(b, i) = cgeti(lgefint(p));
    1110         112 :     av = avma;
    1111       49826 :     for (i = 3; i<m+2; i++)
    1112             :     {
    1113       49714 :       W = f(E, W);
    1114       49714 :       affii(gel(W, col),gel(b, i));
    1115       49714 :       if (gc_needed(av,1))
    1116             :       {
    1117        2378 :         if (DEBUGMEM>1) pari_warn(warnmem,"Wiedemann: first loop, %ld",i);
    1118        2378 :         W = gerepileupto(av, W);
    1119             :       }
    1120             :     }
    1121         112 :     b = FpX_renormalize(b, m+2);
    1122         112 :     if (lgpol(b)==0) {avma = btop; continue; }
    1123         112 :     M = FpX_halfgcd(b, monomial(gen_1, m, 0), p);
    1124         112 :     Q = FpX_neg(FpX_normalize(gcoeff(M, 2, 1),p),p);
    1125         112 :     W = B; lQ =lg(Q);
    1126         112 :     if (DEBUGLEVEL) err_printf("Wiedemann: deg. minpoly: %ld\n",lQ-3);
    1127         112 :     V = FpC_Fp_mul(W, gel(Q, lQ-2), p);
    1128         112 :     av = avma;
    1129       24374 :     for (i = lQ-3; i > 1; i--)
    1130             :     {
    1131       24262 :       W = f(E, W);
    1132       24262 :       V = ZC_lincomb(gen_1, gel(Q,i), V, W);
    1133       24262 :       if (gc_needed(av,1))
    1134             :       {
    1135        1949 :         if (DEBUGMEM>1) pari_warn(warnmem,"Wiedemann: second loop, %ld",i);
    1136        1949 :         gerepileall(av, 2, &V, &W);
    1137             :       }
    1138             :     }
    1139         112 :     V = FpC_red(V, p);
    1140         112 :     W = FpC_sub(f(E,V), B, p);
    1141         224 :     if (ZV_equal0(W)) return gerepilecopy(ltop, V);
    1142           0 :     av = avma;
    1143           0 :     for (i = 1; i <= n; ++i)
    1144             :     {
    1145           0 :       V = W;
    1146           0 :       W = f(E, V);
    1147           0 :       if (ZV_equal0(W))
    1148           0 :         return gerepilecopy(ltop, shallowtrans(V));
    1149           0 :       gerepileall(av, 2, &V, &W);
    1150             :     }
    1151           0 :     avma = btop;
    1152             :   }
    1153           0 :   return NULL;
    1154             : }
    1155             : 
    1156             : GEN
    1157           0 : zMs_ZC_mul(GEN M, GEN B)
    1158             : {
    1159             :   long i, j;
    1160           0 :   long n = lg(B)-1;
    1161           0 :   GEN V = zerocol(n);
    1162           0 :   for (i = 1; i <= n; ++i)
    1163           0 :     if (signe(gel(B, i)))
    1164             :     {
    1165           0 :       GEN R = gel(M, i), C = gel(R, 1), E = gel(R, 2);
    1166           0 :       long l = lg(C);
    1167           0 :       for (j = 1; j < l; ++j)
    1168             :       {
    1169           0 :         long k = C[j];
    1170           0 :         switch(E[j])
    1171             :         {
    1172             :         case 1:
    1173           0 :           gel(V, k) = gel(V,k)==gen_0 ? gel(B,i) : addii(gel(V, k), gel(B,i));
    1174           0 :           break;
    1175             :         case -1:
    1176           0 :           gel(V, k) = gel(V,k)==gen_0 ? negi(gel(B,i)) : subii(gel(V, k), gel(B,i));
    1177           0 :           break;
    1178             :         default:
    1179           0 :           gel(V, k) = gel(V,k)==gen_0 ? mulis(gel(B, i), E[j]) : addii(gel(V, k), mulis(gel(B, i), E[j]));
    1180           0 :           break;
    1181             :         }
    1182             :       }
    1183             :     }
    1184           0 :   return V;
    1185             : }
    1186             : 
    1187             : GEN
    1188           0 : FpMs_FpC_mul(GEN M, GEN B, GEN p) { return FpC_red(zMs_ZC_mul(M, B), p); }
    1189             : 
    1190             : GEN
    1191       74200 : ZV_zMs_mul(GEN B, GEN M)
    1192             : {
    1193             :   long i, j;
    1194       74200 :   long m = lg(M)-1;
    1195       74200 :   GEN V = cgetg(m+1,t_VEC);
    1196    40371114 :   for (i = 1; i <= m; ++i)
    1197             :   {
    1198    40296914 :     GEN R = gel(M, i), C = gel(R, 1), E = gel(R, 2);
    1199    40296914 :     long l = lg(C);
    1200    40296914 :     GEN z = mulis(gel(B, C[1]), E[1]);
    1201   338856266 :     for (j = 2; j < l; ++j)
    1202             :     {
    1203   298559352 :       long k = C[j];
    1204   298559352 :       switch(E[j])
    1205             :       {
    1206             :       case 1:
    1207   203977662 :         z = addii(z, gel(B,k));
    1208   203977662 :         break;
    1209             :       case -1:
    1210    81072089 :         z = subii(z, gel(B,k));
    1211    81072089 :         break;
    1212             :       default:
    1213    13509601 :         z = addii(z, mulis(gel(B,k), E[j]));
    1214    13509601 :         break;
    1215             :       }
    1216             :     }
    1217    40296914 :     gel(V,i) = z;
    1218             :   }
    1219       74200 :   return V;
    1220             : }
    1221             : 
    1222             : GEN
    1223         112 : FpV_FpMs_mul(GEN B, GEN M, GEN p) { return FpV_red(ZV_zMs_mul(B, M), p); }
    1224             : 
    1225             : GEN
    1226        2097 : ZlM_gauss(GEN a, GEN b, ulong p, long e, GEN C)
    1227             : {
    1228        2097 :   pari_sp av = avma, av2;
    1229        2097 :   GEN xi, xb, pi = gen_1, P;
    1230             :   long i;
    1231        2097 :   if (!C) {
    1232           0 :     C = Flm_inv(ZM_to_Flm(a, p), p);
    1233           0 :     if (!C) pari_err_INV("ZlM_gauss", a);
    1234             :   }
    1235        2097 :   P = utoipos(p);
    1236        2097 :   av2 = avma;
    1237        2097 :   xi = Flm_mul(C, ZM_to_Flm(b, p), p);
    1238        2097 :   xb = Flm_to_ZM(xi);
    1239        5353 :   for (i = 2; i <= e; i++)
    1240             :   {
    1241        3256 :     pi = muliu(pi, p); /* = p^(i-1) */
    1242        3256 :     b = ZM_Z_divexact(ZM_sub(b, ZM_nm_mul(a, xi)), P);
    1243        3256 :     if (gc_needed(av,2))
    1244             :     {
    1245          31 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZlM_gauss. i=%ld",i);
    1246          31 :       gerepileall(av2,3, &pi,&b,&xb);
    1247             :     }
    1248        3256 :     xi = Flm_mul(C, ZM_to_Flm(b, p), p);
    1249        3256 :     xb = ZM_add(xb, nm_Z_mul(xi, pi));
    1250             :   }
    1251        2097 :   return gerepileupto(av, xb);
    1252             : }
    1253             : 
    1254             : struct wrapper_modp_s {
    1255             :   GEN (*f)(void*, GEN);
    1256             :   void *E;
    1257             :   GEN p;
    1258             : };
    1259             : 
    1260             : /* compute f(x) mod p */
    1261             : static GEN
    1262       74088 : wrap_relcomb_modp(void *E, GEN x)
    1263             : {
    1264       74088 :   struct wrapper_modp_s *W = (struct wrapper_modp_s*)E;
    1265       74088 :   return FpC_red(W->f(W->E, x), W->p);
    1266             : }
    1267             : static GEN
    1268           0 : wrap_relcomb(void*E, GEN x) { return zMs_ZC_mul((GEN)E, x); }
    1269             : 
    1270             : static GEN
    1271       74088 : wrap_relker(void*E, GEN x) { return ZV_zMs_mul(x, (GEN)E); }
    1272             : 
    1273             : /* Solve f(X) = B (mod p^e); blackbox version of ZlM_gauss */
    1274             : GEN
    1275           0 : gen_ZpM_Dixon(void *E, GEN (*f)(void*, GEN), GEN B, GEN p, long e)
    1276             : {
    1277             :   struct wrapper_modp_s W;
    1278           0 :   pari_sp av = avma;
    1279           0 :   GEN xb, xi, pi = gen_1;
    1280             :   long i;
    1281           0 :   W.E = E;
    1282           0 :   W.f = f;
    1283           0 :   W.p = p;
    1284           0 :   xi = gen_FpM_Wiedemann((void*)&W, wrap_relcomb_modp, FpC_red(B, p), p); /* f^(-1) B */
    1285           0 :   if (!xi || e == 1 || typ(xi) == t_VEC) return xi;
    1286           0 :   xb = xi;
    1287           0 :   for (i = 2; i <= e; i++)
    1288             :   {
    1289           0 :     pi = mulii(pi, p); /* = p^(i-1) */
    1290           0 :     B = ZC_Z_divexact(ZC_sub(B, f(E, xi)), p);
    1291           0 :     if (gc_needed(av,2))
    1292             :     {
    1293           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"gen_ZpM_Dixon. i=%ld",i);
    1294           0 :       gerepileall(av,3, &pi,&B,&xb);
    1295             :     }
    1296           0 :     xi = gen_FpM_Wiedemann((void*)&W, wrap_relcomb_modp, FpC_red(B, p), p);
    1297           0 :     if (!xi) return NULL;
    1298           0 :     if (typ(xi) == t_VEC) return gerepileupto(av, xi);
    1299           0 :     xb = ZC_add(xb, ZC_Z_mul(xi, pi));
    1300             :   }
    1301           0 :   return gerepileupto(av, xb);
    1302             : }
    1303             : 
    1304             : static GEN
    1305       24857 : vecprow(GEN A, GEN prow)
    1306             : {
    1307       24857 :   return mkvec2(vecsmallpermute(prow,gel(A,1)), gel(A,2));
    1308             : }
    1309             : 
    1310             : /* Solve the equation MX = A. Return either a solution as a t_COL,
    1311             :  * or the index of a column which is linearly dependent from the others as a
    1312             :  * t_VECSMALL with a single component. */
    1313             : GEN
    1314           0 : ZpMs_ZpCs_solve(GEN M, GEN A, long nbrow, GEN p, long e)
    1315             : {
    1316           0 :   pari_sp av = avma;
    1317             :   GEN pcol, prow;
    1318           0 :   long nbi=lg(M)-1, lR;
    1319             :   long i, n;
    1320             :   GEN Mp, Ap, Rp;
    1321             :   pari_timer ti;
    1322           0 :   if (DEBUGLEVEL) timer_start(&ti);
    1323           0 :   RgMs_structelim(M, nbrow, gel(A, 1), &pcol, &prow);
    1324           0 :   if (!pcol) { avma = av; return NULL; }
    1325           0 :   if (DEBUGLEVEL)
    1326           0 :     timer_printf(&ti,"structured elimination (%ld -> %ld)",nbi,lg(pcol)-1);
    1327           0 :   n = lg(pcol)-1;
    1328           0 :   Mp = cgetg(n+1, t_MAT);
    1329           0 :   for(i=1; i<=n; i++)
    1330           0 :     gel(Mp, i) = vecprow(gel(M,pcol[i]), prow);
    1331           0 :   Ap = zCs_to_ZC(vecprow(A, prow), n);
    1332           0 :   if (DEBUGLEVEL) timer_start(&ti);
    1333           0 :   Rp = gen_ZpM_Dixon((void*)Mp,wrap_relcomb, Ap, p, e);
    1334           0 :   if (DEBUGLEVEL) timer_printf(&ti,"linear algebra");
    1335           0 :   if (!Rp) { avma = av; return NULL; }
    1336           0 :   lR = lg(Rp)-1;
    1337           0 :   if (typ(Rp) == t_COL)
    1338             :   {
    1339           0 :     GEN R = zerocol(nbi+1);
    1340           0 :     for (i=1; i<=lR; i++)
    1341           0 :        gel(R,pcol[i]) = gel(Rp,i);
    1342           0 :     return gerepilecopy(av, R);
    1343             :   }
    1344           0 :   for (i = 1; i <= lR; ++i)
    1345           0 :     if (signe(gel(Rp, i)))
    1346           0 :       return gerepileuptoleaf(av, mkvecsmall(pcol[i]));
    1347           0 :   return NULL; /* NOT REACHED */
    1348             : }
    1349             : 
    1350             : GEN
    1351           0 : FpMs_FpCs_solve(GEN M, GEN A, long nbrow, GEN p)
    1352             : {
    1353           0 :   return ZpMs_ZpCs_solve(M, A, nbrow, p, 1);
    1354             : }
    1355             : 
    1356             : GEN
    1357           0 : FpMs_FpCs_solve_safe(GEN M, GEN A, long nbrow, GEN p)
    1358             : {
    1359             :   GEN res;
    1360           0 :   pari_CATCH(e_INV)
    1361             :   {
    1362           0 :     GEN E = pari_err_last();
    1363           0 :     GEN x = err_get_compo(E,2);
    1364           0 :     if (typ(x) != t_INTMOD) pari_err(0,E);
    1365           0 :     if (DEBUGLEVEL)
    1366           0 :       pari_warn(warner,"FpMs_FpCs_solve_safe, impossible inverse %Ps", x);
    1367           0 :     res = NULL;
    1368             :   } pari_TRY {
    1369           0 :     res = ZpMs_ZpCs_solve(M, A, nbrow, p, 1);
    1370           0 :   } pari_ENDCATCH
    1371           0 :   return res;
    1372             : }
    1373             : 
    1374             : static GEN
    1375         112 : FpMs_structelim_back(GEN M, GEN V, GEN prow, GEN p)
    1376             : {
    1377         112 :   long i, j, oldf = 0, f = 0;
    1378         112 :   long lrow = lg(prow), lM = lg(M);
    1379         112 :   GEN W = const_vecsmall(lM-1,1);
    1380         112 :   GEN R = cgetg(lrow, t_VEC);
    1381      192031 :   for (i=1; i<lrow; i++)
    1382      191919 :     gel(R,i) = prow[i]==0 ? NULL: gel(V,prow[i]);
    1383             :   do
    1384             :   {
    1385         357 :     oldf = f;
    1386      402493 :     for(i=1; i<lM; i++)
    1387      402136 :       if (W[i])
    1388             :       {
    1389      135905 :         GEN C = gel(M,i), F = gel(C,1), E = gel(C,2);
    1390      135905 :         long c=0, cj=0, lF = lg(F);
    1391     1191477 :         for(j=1; j<lF; j++)
    1392     1055572 :           if (!gel(R,F[j]))
    1393       75866 :           { c++; cj=j; }
    1394      135905 :         if (c>=2) continue;
    1395      118545 :         if (c==1)
    1396             :         {
    1397       37471 :           pari_sp av = avma;
    1398       37471 :           GEN s = gen_0;
    1399      327208 :           for(j=1; j<lF; j++)
    1400      289737 :             if (j!=cj) s = Fp_add(s, mulis(gel(R,F[j]), E[j]), p);
    1401       37471 :           gel(R,F[cj]) = gerepileupto(av, Fp_div(Fp_neg(s, p), stoi(E[cj]), p));
    1402       37471 :           f++;
    1403             :         }
    1404      118545 :         W[i]=0;
    1405             :       }
    1406         357 :   } while(oldf!=f);
    1407      192031 :   for (i=1; i<lrow; i++)
    1408      191919 :     if (!gel(R,i)) gel(R,i) = gen_0;
    1409         112 :   return R;
    1410             : }
    1411             : 
    1412             : /* Return a linear form Y such that YM is essentially 0 */
    1413             : GEN
    1414         112 : FpMs_leftkernel_elt_col(GEN M, long nbcol, long nbrow, GEN p)
    1415             : {
    1416         112 :   pari_sp av = avma, av2;
    1417             :   GEN pcol, prow;
    1418             :   long i, n;
    1419             :   GEN Mp, B, MB, R, Rp;
    1420             :   pari_timer ti;
    1421             :   struct wrapper_modp_s W;
    1422         112 :   if (DEBUGLEVEL) timer_start(&ti);
    1423         112 :   RgMs_structelim_col(M, nbcol, nbrow, cgetg(1,t_VECSMALL), &pcol, &prow);
    1424         112 :   if (!pcol) { avma = av; return NULL; }
    1425         112 :   if (DEBUGLEVEL)
    1426           0 :     timer_printf(&ti,"structured elimination (%ld -> %ld)",nbcol,lg(pcol)-1);
    1427         112 :   n = lg(pcol)-1;
    1428         112 :   Mp = cgetg(n+1, t_MAT);
    1429       24969 :   for(i=1; i<=n; i++)
    1430       24857 :     gel(Mp, i) = vecprow(gel(M,pcol[i]), prow);
    1431         112 :   W.E = (void*) Mp;
    1432         112 :   W.f = wrap_relker;
    1433         112 :   W.p = p;
    1434         112 :   av2 = avma;
    1435             :   for(;;)
    1436             :   {
    1437         112 :     avma = av2;
    1438         112 :     B = cgetg(n+1,t_VEC);
    1439       24969 :     for(i=1; i<=n; i++)
    1440       24857 :       gel(B,i) = randomi(p);
    1441         112 :     MB = FpV_FpMs_mul(B, Mp, p);
    1442         112 :     if (DEBUGLEVEL) timer_start(&ti);
    1443         112 :     pari_CATCH(e_INV)
    1444             :     {
    1445           0 :       GEN E = pari_err_last();
    1446           0 :       GEN x = err_get_compo(E,2);
    1447           0 :       if (typ(x) != t_INTMOD) pari_err(0,E);
    1448           0 :       if (DEBUGLEVEL)
    1449           0 :         pari_warn(warner,"FpMs_leftkernel_elt, impossible inverse %Ps", x);
    1450           0 :       Rp = NULL;
    1451             :     } pari_TRY {
    1452         112 :       Rp = gen_FpM_Wiedemann((void*)&W, wrap_relcomb_modp, MB, p);
    1453         112 :     } pari_ENDCATCH
    1454         112 :     if (!Rp) continue;
    1455         112 :     if (typ(Rp)==t_COL)
    1456             :     {
    1457         112 :       Rp = FpC_sub(Rp,B,p);
    1458         112 :       if (ZV_equal0(Rp)) continue;
    1459             :     }
    1460         112 :     R = FpMs_structelim_back(M, Rp, prow, p);
    1461         112 :     if (DEBUGLEVEL) timer_printf(&ti,"Wiedemann left kernel");
    1462         112 :     return gerepilecopy(av, R);
    1463           0 :   }
    1464             : }
    1465             : 
    1466             : GEN
    1467          28 : FpMs_leftkernel_elt(GEN M, long nbrow, GEN p)
    1468             : {
    1469          28 :   return FpMs_leftkernel_elt_col(M, lg(M)-1, nbrow, p);
    1470             : }

Generated by: LCOV version 1.11