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 to exceed 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 - F2v.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 25819-e703fe1174) Lines: 517 637 81.2 %
Date: 2020-09-18 06:10:04 Functions: 47 66 71.2 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2012-2019 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             : /**                             F2v                                   **/
      20             : /**                                                                   **/
      21             : /***********************************************************************/
      22             : /* F2v objects are defined as follows:
      23             :    An F2v is a t_VECSMALL:
      24             :    v[0] = codeword
      25             :    v[1] = number of components
      26             :    x[2] = a_0...a_31 x[3] = a_32..a_63, etc.  on 32bit
      27             :    x[2] = a_0...a_63 x[3] = a_64..a_127, etc. on 64bit
      28             : 
      29             :    where the a_i are bits.
      30             : */
      31             : 
      32             : int
      33           0 : F2v_equal0(GEN V)
      34             : {
      35           0 :   long l = lg(V);
      36           0 :   while (--l > 1)
      37           0 :     if (V[l]) return 0;
      38           0 :   return 1;
      39             : }
      40             : 
      41             : GEN
      42      921759 : F2c_to_ZC(GEN x)
      43             : {
      44      921759 :   long l = x[1]+1, lx = lg(x);
      45      921759 :   GEN  z = cgetg(l, t_COL);
      46             :   long i, j, k;
      47     1846635 :   for (i=2, k=1; i<lx; i++)
      48    10664953 :     for (j=0; j<BITS_IN_LONG && k<l; j++,k++)
      49     9740077 :       gel(z,k) = (x[i]&(1UL<<j))? gen_1: gen_0;
      50      921759 :   return z;
      51             : }
      52             : GEN
      53        1610 : F2c_to_mod(GEN x)
      54             : {
      55        1610 :   long l = x[1]+1, lx = lg(x);
      56        1610 :   GEN  z = cgetg(l, t_COL);
      57        1610 :   GEN _0 = mkintmod(gen_0,gen_2);
      58        1610 :   GEN _1 = mkintmod(gen_1,gen_2);
      59             :   long i, j, k;
      60        8020 :   for (i=2, k=1; i<lx; i++)
      61      287047 :     for (j=0; j<BITS_IN_LONG && k<l; j++,k++)
      62      280637 :       gel(z,k) = (x[i]&(1UL<<j))? _1: _0;
      63        1610 :   return z;
      64             : }
      65             : 
      66             : /* x[a..b], a <= b */
      67             : GEN
      68          28 : F2v_slice(GEN x, long a, long b)
      69             : {
      70          28 :   long i,j,k, l = b-a+1;
      71          28 :   GEN z = cgetg(nbits2lg(l), t_VECSMALL);
      72          28 :   z[1] = l;
      73          98 :   for(i=a,k=1,j=BITS_IN_LONG; i<=b; i++,j++)
      74             :   {
      75          70 :     if (j==BITS_IN_LONG) { j=0; z[++k]=0; }
      76          70 :     if (F2v_coeff(x,i)) z[k] |= 1UL<<j;
      77             :   }
      78          28 :   return z;
      79             : }
      80             : /* x[a..b,], a <= b */
      81             : GEN
      82          14 : F2m_rowslice(GEN x, long a, long b)
      83             : {
      84             :   long i, l;
      85          14 :   GEN y = cgetg_copy(x, &l);
      86          42 :   for (i = 1; i < l; i++) gel(y,i) = F2v_slice(gel(x,i),a,b);
      87          14 :   return y;
      88             : }
      89             : 
      90             : GEN
      91      203392 : F2m_to_ZM(GEN z)
      92             : {
      93      203392 :   long i, l = lg(z);
      94      203392 :   GEN x = cgetg(l,t_MAT);
      95     1070130 :   for (i=1; i<l; i++) gel(x,i) = F2c_to_ZC(gel(z,i));
      96      203392 :   return x;
      97             : }
      98             : GEN
      99          84 : F2m_to_mod(GEN z)
     100             : {
     101          84 :   long i, l = lg(z);
     102          84 :   GEN x = cgetg(l,t_MAT);
     103        1680 :   for (i=1; i<l; i++) gel(x,i) = F2c_to_mod(gel(z,i));
     104          84 :   return x;
     105             : }
     106             : 
     107             : GEN
     108           0 : F2v_to_Flv(GEN x)
     109             : {
     110           0 :   long l = x[1]+1, lx = lg(x);
     111           0 :   GEN  z = cgetg(l, t_VECSMALL);
     112             :   long i,j,k;
     113           0 :   for (i=2, k=1; i<lx; i++)
     114           0 :     for (j=0; j<BITS_IN_LONG && k<l; j++,k++)
     115           0 :       z[k] = (x[i]>>j)&1UL;
     116           0 :   return z;
     117             : }
     118             : 
     119             : GEN
     120           0 : F2m_to_Flm(GEN z)
     121             : {
     122           0 :   long i, l = lg(z);
     123           0 :   GEN x = cgetg(l,t_MAT);
     124           0 :   for (i=1; i<l; i++) gel(x,i) = F2v_to_Flv(gel(z,i));
     125           0 :   return x;
     126             : }
     127             : 
     128             : GEN
     129     2796384 : ZV_to_F2v(GEN x)
     130             : {
     131     2796384 :   long l = lg(x)-1;
     132     2796384 :   GEN z = cgetg(nbits2lg(l), t_VECSMALL);
     133             :   long i,j,k;
     134     2796384 :   z[1] = l;
     135    36776141 :   for(i=1,k=1,j=BITS_IN_LONG; i<=l; i++,j++)
     136             :   {
     137    33979757 :     if (j==BITS_IN_LONG) { j=0; z[++k]=0; }
     138    33979757 :     if (mpodd(gel(x,i))) z[k] |= 1UL<<j;
     139             :   }
     140     2796384 :   return z;
     141             : }
     142             : 
     143             : GEN
     144        4788 : RgV_to_F2v(GEN x)
     145             : {
     146        4788 :   long l = lg(x)-1;
     147        4788 :   GEN z = cgetg(nbits2lg(l), t_VECSMALL);
     148             :   long i,j,k;
     149        4788 :   z[1] = l;
     150      846650 :   for(i=1,k=1,j=BITS_IN_LONG; i<=l; i++,j++)
     151             :   {
     152      841862 :     if (j==BITS_IN_LONG) { j=0; z[++k]=0; }
     153      841862 :     if (Rg_to_F2(gel(x,i))) z[k] |= 1UL<<j;
     154             :   }
     155        4788 :   return z;
     156             : }
     157             : 
     158             : GEN
     159         665 : Flv_to_F2v(GEN x)
     160             : {
     161         665 :   long l = lg(x)-1;
     162         665 :   GEN z = cgetg(nbits2lg(l), t_VECSMALL);
     163             :   long i,j,k;
     164         665 :   z[1] = l;
     165        5369 :   for(i=1,k=1,j=BITS_IN_LONG; i<=l; i++,j++)
     166             :   {
     167        4704 :     if (j==BITS_IN_LONG) { j=0; z[++k]=0; }
     168        4704 :     if (x[i]&1L) z[k] |= 1UL<<j;
     169             :   }
     170         665 :   return z;
     171             : }
     172             : 
     173             : GEN
     174     3198994 : ZM_to_F2m(GEN x) { pari_APPLY_same(ZV_to_F2v(gel(x,i))) }
     175             : GEN
     176        5005 : RgM_to_F2m(GEN x) { pari_APPLY_same(RgV_to_F2v(gel(x,i))) }
     177             : GEN
     178           0 : Flm_to_F2m(GEN x) { pari_APPLY_same(Flv_to_F2v(gel(x,i))) }
     179             : 
     180             : GEN
     181      395907 : const_F2v(long m)
     182             : {
     183      395907 :   long i, l = nbits2lg(m);
     184      395907 :   GEN c = cgetg(l, t_VECSMALL);
     185      395907 :   c[1] = m;
     186      798053 :   for (i = 2; i < l; i++) uel(c,i) = -1UL;
     187      395907 :   if (remsBIL(m)) uel(c,l-1) = (1UL<<remsBIL(m))-1UL;
     188      395907 :   return c;
     189             : }
     190             : 
     191             : /* Allow lg(y)<lg(x) */
     192             : void
     193    15198623 : F2v_add_inplace(GEN x, GEN y)
     194             : {
     195    15198623 :   long n = lg(y);
     196    15198623 :   long r = (n-2)&7L, q = n-r, i;
     197    24502718 :   for (i = 2; i < q; i += 8)
     198             :   {
     199     9304095 :     x[  i] ^= y[  i]; x[1+i] ^= y[1+i]; x[2+i] ^= y[2+i]; x[3+i] ^= y[3+i];
     200     9304095 :     x[4+i] ^= y[4+i]; x[5+i] ^= y[5+i]; x[6+i] ^= y[6+i]; x[7+i] ^= y[7+i];
     201             :   }
     202    15198623 :   switch (r)
     203             :   {
     204      907264 :     case 7: x[i] ^= y[i]; i++; case 6: x[i] ^= y[i]; i++;
     205     2337515 :     case 5: x[i] ^= y[i]; i++; case 4: x[i] ^= y[i]; i++;
     206     4871147 :     case 3: x[i] ^= y[i]; i++; case 2: x[i] ^= y[i]; i++;
     207    12418686 :     case 1: x[i] ^= y[i]; i++;
     208             :   }
     209    15198623 : }
     210             : 
     211             : /* Allow lg(y)<lg(x) */
     212             : void
     213           0 : F2v_and_inplace(GEN x, GEN y)
     214             : {
     215           0 :   long n = lg(y);
     216           0 :   long r = (n-2)&7L, q = n-r, i;
     217           0 :   for (i = 2; i < q; i += 8)
     218             :   {
     219           0 :     x[  i] &= y[  i]; x[1+i] &= y[1+i]; x[2+i] &= y[2+i]; x[3+i] &= y[3+i];
     220           0 :     x[4+i] &= y[4+i]; x[5+i] &= y[5+i]; x[6+i] &= y[6+i]; x[7+i] &= y[7+i];
     221             :   }
     222           0 :   switch (r)
     223             :   {
     224           0 :     case 7: x[i] &= y[i]; i++; case 6: x[i] &= y[i]; i++;
     225           0 :     case 5: x[i] &= y[i]; i++; case 4: x[i] &= y[i]; i++;
     226           0 :     case 3: x[i] &= y[i]; i++; case 2: x[i] &= y[i]; i++;
     227           0 :     case 1: x[i] &= y[i]; i++;
     228             :   }
     229           0 : }
     230             : 
     231             : /* Allow lg(y)<lg(x) */
     232             : void
     233           0 : F2v_or_inplace(GEN x, GEN y)
     234             : {
     235           0 :   long n = lg(y);
     236           0 :   long r = (n-2)&7L, q = n-r, i;
     237           0 :   for (i = 2; i < q; i += 8)
     238             :   {
     239           0 :     x[  i] |= y[  i]; x[1+i] |= y[1+i]; x[2+i] |= y[2+i]; x[3+i] |= y[3+i];
     240           0 :     x[4+i] |= y[4+i]; x[5+i] |= y[5+i]; x[6+i] |= y[6+i]; x[7+i] |= y[7+i];
     241             :   }
     242           0 :   switch (r)
     243             :   {
     244           0 :     case 7: x[i] |= y[i]; i++; case 6: x[i] |= y[i]; i++;
     245           0 :     case 5: x[i] |= y[i]; i++; case 4: x[i] |= y[i]; i++;
     246           0 :     case 3: x[i] |= y[i]; i++; case 2: x[i] |= y[i]; i++;
     247           0 :     case 1: x[i] |= y[i]; i++;
     248             :   }
     249           0 : }
     250             : 
     251             : /* Allow lg(y)<lg(x) */
     252             : void
     253           0 : F2v_negimply_inplace(GEN x, GEN y)
     254             : {
     255           0 :   long n = lg(y);
     256           0 :   long r = (n-2)&7L, q = n-r, i;
     257           0 :   for (i = 2; i < q; i += 8)
     258             :   {
     259           0 :     x[  i] &= ~y[  i]; x[1+i] &= ~y[1+i]; x[2+i] &= ~y[2+i]; x[3+i] &= ~y[3+i];
     260           0 :     x[4+i] &= ~y[4+i]; x[5+i] &= ~y[5+i]; x[6+i] &= ~y[6+i]; x[7+i] &= ~y[7+i];
     261             :   }
     262           0 :   switch (r)
     263             :   {
     264           0 :     case 7: x[i] &= ~y[i]; i++; case 6: x[i] &= ~y[i]; i++;
     265           0 :     case 5: x[i] &= ~y[i]; i++; case 4: x[i] &= ~y[i]; i++;
     266           0 :     case 3: x[i] &= ~y[i]; i++; case 2: x[i] &= ~y[i]; i++;
     267           0 :     case 1: x[i] &= ~y[i]; i++;
     268             :   }
     269           0 : }
     270             : 
     271             : ulong
     272           0 : F2v_dotproduct(GEN x, GEN y)
     273             : {
     274           0 :   long i, lx = lg(x);
     275             :   ulong c;
     276           0 :   if (lx <= 2) return 0;
     277           0 :   c = uel(x,2) & uel(y,2);
     278           0 :   for (i=3; i<lx; i++) c ^= uel(x,i) & uel(y,i);
     279             : #ifdef LONG_IS_64BIT
     280           0 :   c ^= c >> 32;
     281             : #endif
     282           0 :   c ^= c >> 16;
     283           0 :   c ^= c >>  8;
     284           0 :   c ^= c >>  4;
     285           0 :   c ^= c >>  2;
     286           0 :   c ^= c >>  1;
     287           0 :   return c & 1;
     288             : }
     289             : 
     290             : ulong
     291           0 : F2v_hamming(GEN H)
     292             : {
     293           0 :   long i, n=0, l=lg(H);
     294           0 :   for (i=2; i<l; i++) n += hammingl(uel(H,i));
     295           0 :   return n;
     296             : }
     297             : 
     298             : GEN
     299       24656 : matid_F2m(long n)
     300             : {
     301       24656 :   GEN y = cgetg(n+1,t_MAT);
     302             :   long i;
     303       24656 :   if (n < 0) pari_err_DOMAIN("matid_F2m", "dimension","<",gen_0,stoi(n));
     304      187912 :   for (i=1; i<=n; i++) { gel(y,i) = zero_F2v(n); F2v_set(gel(y,i),i); }
     305       24656 :   return y;
     306             : }
     307             : 
     308             : GEN
     309           0 : F2m_row(GEN x, long j)
     310             : {
     311           0 :   long i, l = lg(x);
     312           0 :   GEN V = zero_F2v(l-1);
     313           0 :   for(i = 1; i < l; i++)
     314           0 :     if (F2m_coeff(x,j,i))
     315           0 :       F2v_set(V,i);
     316           0 :   return V;
     317             : }
     318             : 
     319             : GEN
     320           0 : F2m_transpose(GEN x)
     321             : {
     322           0 :   long i, dx, lx = lg(x);
     323             :   GEN y;
     324           0 :   if (lx == 1) return cgetg(1,t_MAT);
     325           0 :   dx = coeff(x,1,1); y = cgetg(dx+1,t_MAT);
     326           0 :   for (i=1; i<=dx; i++) gel(y,i) = F2m_row(x,i);
     327           0 :   return y;
     328             : }
     329             : 
     330             : INLINE GEN
     331      507812 : F2m_F2c_mul_i(GEN x, GEN y, long lx, long l)
     332             : {
     333             :   long j;
     334      507812 :   GEN z = NULL;
     335             : 
     336     5896815 :   for (j=1; j<lx; j++)
     337             :   {
     338     5389003 :     if (!F2v_coeff(y,j)) continue;
     339     1234161 :     if (!z) z = vecsmall_copy(gel(x,j));
     340      755686 :     else F2v_add_inplace(z,gel(x,j));
     341             :   }
     342      507812 :   if (!z) z = zero_F2v(l);
     343      507812 :   return z;
     344             : }
     345             : 
     346             : GEN
     347      114196 : F2m_mul(GEN x, GEN y)
     348             : {
     349      114196 :   long i,j,l,lx=lg(x), ly=lg(y);
     350             :   GEN z;
     351      114196 :   if (ly==1) return cgetg(1,t_MAT);
     352      114196 :   z = cgetg(ly,t_MAT);
     353      114196 :   if (lx==1)
     354             :   {
     355           0 :     for (i=1; i<ly; i++) gel(z,i) = mkvecsmall(0);
     356           0 :     return z;
     357             :   }
     358      114196 :   l = coeff(x,1,1);
     359      622008 :   for (j=1; j<ly; j++) gel(z,j) = F2m_F2c_mul_i(x, gel(y,j), lx, l);
     360      114196 :   return z;
     361             : }
     362             : 
     363             : GEN
     364           0 : F2m_F2c_mul(GEN x, GEN y)
     365             : {
     366           0 :   long l, lx = lg(x);
     367           0 :   if (lx==1) return cgetg(1,t_VECSMALL);
     368           0 :   l = coeff(x,1,1);
     369           0 :   return F2m_F2c_mul_i(x, y, lx, l);
     370             : }
     371             : 
     372             : static GEN
     373           0 : _F2m_mul(void *data, GEN x, GEN y) { (void) data; return F2m_mul(x,y); }
     374             : static GEN
     375           0 : _F2m_sqr(void *data, GEN x) { (void) data; return F2m_mul(x,x); }
     376             : GEN
     377           0 : F2m_powu(GEN x, ulong n)
     378             : {
     379           0 :   if (!n) return matid(lg(x)-1);
     380           0 :   return gen_powu(x, n,NULL, &_F2m_sqr, &_F2m_mul);
     381             : }
     382             : 
     383             : static long
     384     2796965 : F2v_find_nonzero(GEN x0, GEN mask0, long m)
     385             : {
     386     2796965 :   ulong *x = (ulong *)x0+2, *mask = (ulong *)mask0+2, e;
     387     2796965 :   long i, l = lg(x0)-2;
     388     3947258 :   for (i = 0; i < l; i++)
     389             :   {
     390     3253288 :     e = *x++ & *mask++;
     391     3253288 :     if (e) return i*BITS_IN_LONG+vals(e)+1;
     392             :   }
     393      693970 :   return m+1;
     394             : }
     395             : 
     396             : /* in place, destroy x */
     397             : GEN
     398      345833 : F2m_ker_sp(GEN x, long deplin)
     399             : {
     400             :   GEN y, c, d;
     401             :   long i, j, k, r, m, n;
     402             : 
     403      345833 :   n = lg(x)-1;
     404      345833 :   m = mael(x,1,1); r=0;
     405             : 
     406      345833 :   d = cgetg(n+1, t_VECSMALL);
     407      345833 :   c = const_F2v(m);
     408     2730022 :   for (k=1; k<=n; k++)
     409             :   {
     410     2417899 :     GEN xk = gel(x,k);
     411     2417899 :     j = F2v_find_nonzero(xk, c, m);
     412     2417898 :     if (j>m)
     413             :     {
     414      549681 :       if (deplin) {
     415       33710 :         GEN c = zero_F2v(n);
     416      150159 :         for (i=1; i<k; i++)
     417      116449 :           if (F2v_coeff(xk, d[i]))
     418       61583 :             F2v_set(c, i);
     419       33710 :         F2v_set(c, k);
     420       33710 :         return c;
     421             :       }
     422      515971 :       r++; d[k] = 0;
     423             :     }
     424             :     else
     425             :     {
     426     1868217 :       F2v_clear(c,j); d[k] = j;
     427     1868218 :       F2v_clear(xk, j);
     428    32567205 :       for (i=k+1; i<=n; i++)
     429             :       {
     430    30698980 :         GEN xi = gel(x,i);
     431    30698980 :         if (F2v_coeff(xi,j)) F2v_add_inplace(xi, xk);
     432             :       }
     433     1868225 :       F2v_set(xk, j);
     434             :     }
     435             :   }
     436      312123 :   if (deplin) return NULL;
     437             : 
     438      310688 :   y = zero_F2m_copy(n,r);
     439      826659 :   for (j=k=1; j<=r; j++,k++)
     440             :   {
     441     1969898 :     GEN C = gel(y,j); while (d[k]) k++;
     442     5906865 :     for (i=1; i<k; i++)
     443     5390894 :       if (d[i] && F2m_coeff(x,d[i],k))
     444     1683760 :         F2v_set(C,i);
     445      515971 :     F2v_set(C, k);
     446             :   }
     447      310689 :   return y;
     448             : }
     449             : 
     450             : /* not memory clean */
     451             : GEN
     452       22022 : F2m_ker(GEN x) { return F2m_ker_sp(F2m_copy(x), 0); }
     453             : GEN
     454           0 : F2m_deplin(GEN x) { return F2m_ker_sp(F2m_copy(x), 1); }
     455             : 
     456             : ulong
     457        1631 : F2m_det_sp(GEN x) { return !F2m_ker_sp(x, 1); }
     458             : 
     459             : ulong
     460           0 : F2m_det(GEN x)
     461             : {
     462           0 :   pari_sp av = avma;
     463           0 :   ulong d = F2m_det_sp(F2m_copy(x));
     464           0 :   return gc_ulong(av, d);
     465             : }
     466             : 
     467             : /* Destroy x */
     468             : GEN
     469       50074 : F2m_gauss_pivot(GEN x, long *rr)
     470             : {
     471             :   GEN c, d;
     472             :   long i, j, k, r, m, n;
     473             : 
     474       50074 :   n = lg(x)-1; if (!n) { *rr=0; return NULL; }
     475       50074 :   m = mael(x,1,1); r=0;
     476             : 
     477       50074 :   d = cgetg(n+1, t_VECSMALL);
     478       50074 :   c = const_F2v(m);
     479      429142 :   for (k=1; k<=n; k++)
     480             :   {
     481      379068 :     GEN xk = gel(x,k);
     482      379068 :     j = F2v_find_nonzero(xk, c, m);
     483      379068 :     if (j>m) { r++; d[k] = 0; }
     484             :     else
     485             :     {
     486      234776 :       F2v_clear(c,j); d[k] = j;
     487     3289167 :       for (i=k+1; i<=n; i++)
     488             :       {
     489     3054391 :         GEN xi = gel(x,i);
     490     3054391 :         if (F2v_coeff(xi,j)) F2v_add_inplace(xi, xk);
     491             :       }
     492             :     }
     493             :   }
     494             : 
     495       50074 :   *rr = r; return gc_const((pari_sp)d, d);
     496             : }
     497             : 
     498             : long
     499          63 : F2m_rank(GEN x)
     500             : {
     501          63 :   pari_sp av = avma;
     502             :   long r;
     503          63 :   (void)F2m_gauss_pivot(F2m_copy(x),&r);
     504          63 :   return gc_long(av, lg(x)-1 - r);
     505             : }
     506             : 
     507             : static GEN
     508          14 : F2m_inv_upper_1_ind(GEN A, long index)
     509             : {
     510          14 :   pari_sp av = avma;
     511          14 :   long n = lg(A)-1, i = index, j;
     512          14 :   GEN u = const_vecsmall(n, 0);
     513          14 :   u[i] = 1;
     514          21 :   for (i--; i>0; i--)
     515             :   {
     516           7 :     ulong m = F2m_coeff(A,i,i+1) & uel(u,i+1); /* j = i+1 */
     517           7 :     for (j=i+2; j<=n; j++) m ^= F2m_coeff(A,i,j) & uel(u,j);
     518           7 :     u[i] = m & 1;
     519             :   }
     520          14 :   return gerepileuptoleaf(av, Flv_to_F2v(u));
     521             : }
     522             : static GEN
     523           7 : F2m_inv_upper_1(GEN A)
     524             : {
     525             :   long i, l;
     526           7 :   GEN B = cgetg_copy(A, &l);
     527          21 :   for (i = 1; i < l; i++) gel(B,i) = F2m_inv_upper_1_ind(A, i);
     528           7 :   return B;
     529             : }
     530             : 
     531             : static GEN
     532      163235 : F2_get_col(GEN b, GEN d, long li, long aco)
     533             : {
     534      163235 :   long i, l = nbits2lg(aco);
     535      163235 :   GEN u = cgetg(l, t_VECSMALL);
     536      163235 :   u[1] = aco;
     537     2192742 :   for (i = 1; i <= li; i++)
     538     2029507 :     if (d[i]) /* d[i] can still be 0 if li > aco */
     539             :     {
     540     1960669 :       if (F2v_coeff(b, i))
     541      610936 :         F2v_set(u, d[i]);
     542             :       else
     543     1349733 :         F2v_clear(u, d[i]);
     544             :     }
     545      163235 :   return u;
     546             : }
     547             : 
     548             : /* destroy a, b */
     549             : GEN
     550       24691 : F2m_gauss_sp(GEN a, GEN b)
     551             : {
     552       24691 :   long i, j, k, l, li, bco, aco = lg(a)-1;
     553             :   GEN u, d;
     554             : 
     555       24691 :   if (!aco) return cgetg(1,t_MAT);
     556       24691 :   li = gel(a,1)[1];
     557       24691 :   d = zero_Flv(li);
     558       24691 :   bco = lg(b)-1;
     559      179274 :   for (i=1; i<=aco; i++)
     560             :   {
     561      154604 :     GEN ai = vecsmall_copy(gel(a,i));
     562      154604 :     if (!d[i] && F2v_coeff(ai, i))
     563       79934 :       k = i;
     564             :     else
     565      596568 :       for (k = 1; k <= li; k++) if (!d[k] && F2v_coeff(ai,k)) break;
     566             :     /* found a pivot on row k */
     567      154604 :     if (k > li) return NULL;
     568      154583 :     d[k] = i;
     569             : 
     570             :     /* Clear k-th row but column-wise instead of line-wise */
     571             :     /*  a_ij -= a_i1*a1j/a_11
     572             :        line-wise grouping:  L_j -= a_1j/a_11*L_1
     573             :        column-wise:         C_i -= a_i1/a_11*C_1
     574             :     */
     575      154583 :     F2v_clear(ai,k);
     576     2098137 :     for (l=1; l<=aco; l++)
     577             :     {
     578     1943554 :       GEN al = gel(a,l);
     579     1943554 :       if (F2v_coeff(al,k)) F2v_add_inplace(al,ai);
     580             :     }
     581     2115511 :     for (l=1; l<=bco; l++)
     582             :     {
     583     1960928 :       GEN bl = gel(b,l);
     584     1960928 :       if (F2v_coeff(bl,k)) F2v_add_inplace(bl,ai);
     585             :     }
     586             :   }
     587       24670 :   u = cgetg(bco+1,t_MAT);
     588      187905 :   for (j = 1; j <= bco; j++) gel(u,j) = F2_get_col(gel(b,j), d, li, aco);
     589       24670 :   return u;
     590             : }
     591             : 
     592             : GEN
     593          35 : F2m_gauss(GEN a, GEN b)
     594             : {
     595          35 :   pari_sp av = avma;
     596          35 :   if (lg(a) == 1) return cgetg(1,t_MAT);
     597          35 :   return gerepileupto(av, F2m_gauss_sp(F2m_copy(a), F2m_copy(b)));
     598             : }
     599             : GEN
     600          14 : F2m_F2c_gauss(GEN a, GEN b)
     601             : {
     602          14 :   pari_sp av = avma;
     603          14 :   GEN z = F2m_gauss(a, mkmat(b));
     604          14 :   if (!z) return gc_NULL(av);
     605           7 :   if (lg(z) == 1) { set_avma(av); return cgetg(1,t_VECSMALL); }
     606           7 :   return gerepileuptoleaf(av, gel(z,1));
     607             : }
     608             : 
     609             : GEN
     610          35 : F2m_inv(GEN a)
     611             : {
     612          35 :   pari_sp av = avma;
     613          35 :   if (lg(a) == 1) return cgetg(1,t_MAT);
     614          35 :   return gerepileupto(av, F2m_gauss_sp(F2m_copy(a), matid_F2m(gel(a,1)[1])));
     615             : }
     616             : 
     617             : GEN
     618           7 : F2m_invimage_i(GEN A, GEN B)
     619             : {
     620             :   GEN d, x, X, Y;
     621           7 :   long i, j, nY, nA = lg(A)-1, nB = lg(B)-1;
     622           7 :   x = F2m_ker_sp(shallowconcat(A, B), 0);
     623             :   /* AX = BY, Y in strict upper echelon form with pivots = 1.
     624             :    * We must find T such that Y T = Id_nB then X T = Z. This exists iff
     625             :    * Y has at least nB columns and full rank */
     626           7 :   nY = lg(x)-1;
     627           7 :   if (nY < nB) return NULL;
     628             : 
     629             :   /* implicitly: Y = rowslice(x, nA+1, nA+nB), nB rows */
     630           7 :   d = cgetg(nB+1, t_VECSMALL);
     631          21 :   for (i = nB, j = nY; i >= 1; i--, j--)
     632             :   {
     633          14 :     for (; j>=1; j--)
     634          14 :       if (F2m_coeff(x,nA+i,j)) { d[i] = j; break; } /* Y[i,j] */
     635          14 :     if (!j) return NULL;
     636             :   }
     637           7 :   x = vecpermute(x, d);
     638             : 
     639           7 :   X = F2m_rowslice(x, 1, nA);
     640           7 :   Y = F2m_rowslice(x, nA+1, nA+nB);
     641           7 :   return F2m_mul(X, F2m_inv_upper_1(Y));
     642             : }
     643             : GEN
     644           0 : F2m_invimage(GEN A, GEN B)
     645             : {
     646           0 :   pari_sp av = avma;
     647           0 :   GEN X = F2m_invimage_i(A,B);
     648           0 :   if (!X) return gc_NULL(av);
     649           0 :   return gerepileupto(av, X);
     650             : }
     651             : 
     652             : GEN
     653       21549 : F2m_F2c_invimage(GEN A, GEN y)
     654             : {
     655       21549 :   pari_sp av = avma;
     656       21549 :   long i, l = lg(A);
     657             :   GEN M, x;
     658             : 
     659       21549 :   if (l==1) return NULL;
     660       21549 :   if (lg(y) != lgcols(A)) pari_err_DIM("F2m_F2c_invimage");
     661       21549 :   M = cgetg(l+1,t_MAT);
     662      191165 :   for (i=1; i<l; i++) gel(M,i) = gel(A,i);
     663       21549 :   gel(M,l) = y; M = F2m_ker(M);
     664       21549 :   i = lg(M)-1; if (!i) return gc_NULL(av);
     665             : 
     666       21549 :   x = gel(M,i);
     667       21549 :   if (!F2v_coeff(x,l)) return gc_NULL(av);
     668       21549 :   F2v_clear(x, x[1]); x[1]--; /* remove last coord */
     669       21549 :   return gerepileuptoleaf(av, x);
     670             : }
     671             : 
     672             : /*  Block Lanczos algorithm for kernel of sparse matrix (F2Ms)
     673             :     Based on lanczos.cpp by Jason Papadopoulos
     674             :     <https://github.com/sagemath/FlintQS/blob/master/src/lanczos.cpp>
     675             :     Copyright Jason Papadopoulos 2006
     676             :     Released under the GNU General Public License v2 or later version.
     677             : */
     678             : 
     679             : /* F2Ms are vector of vecsmall which represents non-zero entries of columns
     680             :  * F2w are vecsmall whoses entries are columns of a n x BIL F2m
     681             :  * F2wB are F2w in the special case where dim = BIL.
     682             :  */
     683             : 
     684             : #define BIL BITS_IN_LONG
     685             : 
     686             : static GEN
     687         156 : F2w_transpose_F2m(GEN x)
     688             : {
     689         156 :   long i, j, l = lg(x)-1;
     690         156 :   GEN z = cgetg(BIL+1, t_MAT);
     691        9756 :   for (j = 1; j <= BIL; j++)
     692        9600 :     gel(z,j) = zero_F2v(l);
     693      264532 :   for (i = 1; i <= l; i++)
     694             :   {
     695      264376 :     ulong xi = uel(x,i);
     696    16784568 :     for(j = 1; j <= BIL; j++)
     697    16520192 :       if (xi&(1UL<<(j-1)))
     698     4753976 :         F2v_set(gel(z, j), i);
     699             :   }
     700         156 :   return z;
     701             : }
     702             : 
     703             : static GEN
     704        6096 : F2wB_mul(GEN a, GEN b)
     705             : {
     706             :   long i, j;
     707        6096 :   GEN c = cgetg(BIL+1, t_VECSMALL);
     708      378960 :   for (i = 1; i <= BIL; i++)
     709             :   {
     710      372864 :     ulong s = 0, ai = a[i];
     711    19800840 :     for (j = 1; ai; j++, ai>>=1)
     712    19427976 :       if (ai & 1)
     713     9777182 :         s ^= b[j];
     714      372864 :     c[i] = s;
     715             :   }
     716        6096 :   return c;
     717             : }
     718             : 
     719             : static void
     720        4064 : precompute_F2w_F2wB(GEN x, GEN c)
     721             : {
     722             :   ulong z, xk;
     723             :   ulong i, j, k, index;
     724        4064 :   x++; c++;
     725       35136 :   for (j = 0; j < (BIL>>3); j++)
     726             :   {
     727     7985504 :     for (i = 0; i < 256; i++)
     728             :     {
     729     7954432 :       k = 0;
     730     7954432 :       index = i;
     731     7954432 :       z = 0;
     732    63666528 :       while (index) {
     733    55712096 :         xk = x[k];
     734    55712096 :         if (index & 1)
     735    31817728 :           z ^= xk;
     736    55712096 :         index >>= 1;
     737    55712096 :         k++;
     738             :       }
     739     7954432 :       c[i] = z;
     740             :     }
     741       31072 :     x += 8; c += 256;
     742             :   }
     743        4064 : }
     744             : 
     745             : static void
     746        4064 : F2w_F2wB_mul_add_inplace(GEN v, GEN x, GEN y)
     747             : {
     748        4064 :   long i, n = lg(y)-1;
     749             :   ulong word;
     750        4064 :   GEN c = cgetg(1+(BIL<<5), t_VECSMALL);
     751        4064 :   precompute_F2w_F2wB(x, c);
     752        4064 :   c++;
     753     7403488 :   for (i = 1; i <= n; i++)
     754             :   {
     755     7399424 :     word = v[i];
     756    14798848 :     y[i] ^=  c[ 0*256 + ((word>> 0) & 0xff) ]
     757     7399424 :            ^ c[ 1*256 + ((word>> 8) & 0xff) ]
     758     7399424 :            ^ c[ 2*256 + ((word>>16) & 0xff) ]
     759     7399424 :            ^ c[ 3*256 + ((word>>24) & 0xff) ]
     760             : #ifdef LONG_IS_64BIT
     761     7030112 :            ^ c[ 4*256 + ((word>>32) & 0xff) ]
     762     7030112 :            ^ c[ 5*256 + ((word>>40) & 0xff) ]
     763     7030112 :            ^ c[ 6*256 + ((word>>48) & 0xff) ]
     764     7030112 :            ^ c[ 7*256 + ((word>>56)       ) ]
     765             : #endif
     766             :            ;
     767             :   }
     768        4064 : }
     769             : 
     770             : /* Return x*y~, which is a F2wB */
     771             : static GEN
     772        3091 : F2w_transmul(GEN x, GEN y)
     773             : {
     774        3091 :   long i, j, n = lg(x)-1;
     775        3091 :   GEN z = zero_zv(BIL);
     776        3091 :   pari_sp av = avma;
     777        3091 :   GEN c = zero_zv(BIL<<5) + 1;
     778        3091 :   GEN xy = z + 1;
     779             : 
     780     5618384 :   for (i = 1; i <= n; i++)
     781             :   {
     782     5615293 :     ulong xi = x[i];
     783     5615293 :     ulong yi = y[i];
     784     5615293 :     c[ 0*256 + ( xi        & 0xff) ] ^= yi;
     785     5615293 :     c[ 1*256 + ((xi >>  8) & 0xff) ] ^= yi;
     786     5615293 :     c[ 2*256 + ((xi >> 16) & 0xff) ] ^= yi;
     787     5615293 :     c[ 3*256 + ((xi >> 24) & 0xff) ] ^= yi;
     788             : #ifdef LONG_IS_64BIT
     789     5335366 :     c[ 4*256 + ((xi >> 32) & 0xff) ] ^= yi;
     790     5335366 :     c[ 5*256 + ((xi >> 40) & 0xff) ] ^= yi;
     791     5335366 :     c[ 6*256 + ((xi >> 48) & 0xff) ] ^= yi;
     792     5335366 :     c[ 7*256 + ((xi >> 56)       ) ] ^= yi;
     793             : #endif
     794             :   }
     795       27819 :   for(i = 0; i < 8; i++)
     796             :   {
     797       24728 :     ulong a0 = 0, a1 = 0, a2 = 0, a3 = 0;
     798             : #ifdef LONG_IS_64BIT
     799       22544 :     ulong a4 = 0, a5 = 0, a6 = 0, a7 = 0;
     800             : #endif
     801     6355096 :     for (j = 0; j < 256; j++) {
     802     6330368 :       if ((j >> i) & 1) {
     803     3165184 :         a0 ^= c[0*256 + j];
     804     3165184 :         a1 ^= c[1*256 + j];
     805     3165184 :         a2 ^= c[2*256 + j];
     806     3165184 :         a3 ^= c[3*256 + j];
     807             : #ifdef LONG_IS_64BIT
     808     2885632 :         a4 ^= c[4*256 + j];
     809     2885632 :         a5 ^= c[5*256 + j];
     810     2885632 :         a6 ^= c[6*256 + j];
     811     2885632 :         a7 ^= c[7*256 + j];
     812             : #endif
     813             :       }
     814             :     }
     815       24728 :     xy[ 0] = a0; xy[ 8] = a1; xy[16] = a2; xy[24] = a3;
     816             : #ifdef LONG_IS_64BIT
     817       22544 :     xy[32] = a4; xy[40] = a5; xy[48] = a6; xy[56] = a7;
     818             : #endif
     819       24728 :     xy++;
     820             :   }
     821        3091 :   return gc_const(av, z);
     822             : }
     823             : 
     824             : static GEN
     825        1018 : identity_F2wB(void)
     826             : {
     827             :   long i;
     828        1018 :   GEN M = cgetg(BIL+1, t_VECSMALL);
     829       63290 :   for (i = 1; i <= BIL; i++)
     830       62272 :     uel(M,i) = 1UL<<(i-1);
     831        1018 :   return M;
     832             : }
     833             : 
     834             : static GEN
     835        1018 : find_nonsingular_sub(GEN t, GEN last_s, GEN *pt_s)
     836             : {
     837        1018 :   long i, j, dim = 0;
     838             :   ulong mask, row_i, row_j;
     839        1018 :   long last_dim = lg(last_s)-1;
     840        1018 :   GEN s = cgetg(BIL+1, t_VECSMALL);
     841        1018 :   GEN M1 = identity_F2wB();
     842        1018 :   pari_sp av = avma;
     843        1018 :   GEN cols = cgetg(BIL+1, t_VECSMALL);
     844        1018 :   GEN M0 = zv_copy(t);
     845             : 
     846        1018 :   mask = 0;
     847       62454 :   for (i = 1; i <= last_dim; i++)
     848             :   {
     849       61436 :     cols[BIL + 1 - i] = last_s[i];
     850       61436 :     mask |= 1UL<<(last_s[i]-1);
     851             :   }
     852       63290 :   for (i = j = 1; i <= BIL; i++)
     853       62272 :     if (!(mask & (1UL<<(i-1))))
     854         836 :       cols[j++] = i;
     855             : 
     856             :   /* compute the inverse of t[][] */
     857             : 
     858       63290 :   for (i = 1; i <= BIL; i++)
     859             :   {
     860       62272 :     mask = 1UL<<(cols[i]-1);
     861       62272 :     row_i = cols[i];
     862      124666 :     for (j = i; j <= BIL; j++)
     863             :     {
     864      122938 :       row_j = cols[j];
     865      122938 :       if (uel(M0,row_j) & mask)
     866             :       {
     867       60544 :         swap(gel(M0, row_j), gel(M0, row_i));
     868       60544 :         swap(gel(M1, row_j), gel(M1, row_i));
     869       60544 :         break;
     870             :       }
     871             :     }
     872       62272 :     if (j <= BIL)
     873             :     {
     874     3846464 :       for (j = 1; j <= BIL; j++)
     875             :       {
     876     3785920 :         row_j = cols[j];
     877     3785920 :         if (row_i != row_j && (M0[row_j] & mask))
     878             :         {
     879     1833193 :           uel(M0,row_j) ^= uel(M0,row_i);
     880     1833193 :           uel(M1,row_j) ^= uel(M1,row_i);
     881             :         }
     882             :       }
     883       60544 :       s[++dim] = cols[i];
     884       60544 :       continue;
     885             :     }
     886        1728 :     for (j = i; j <= BIL; j++)
     887             :     {
     888        1728 :       row_j = cols[j];
     889        1728 :       if (uel(M1,row_j) & mask)
     890             :       {
     891        1728 :         swap(gel(M0, row_j), gel(M0, row_i));
     892        1728 :         swap(gel(M1, row_j), gel(M1, row_i));
     893        1728 :         break;
     894             :       }
     895             :     }
     896        1728 :     if (j > BIL) return NULL;
     897      109056 :     for (j = 1; j <= BIL; j++)
     898             :     {
     899      107328 :       row_j = cols[j];
     900      107328 :       if (row_i != row_j && (M1[row_j] & mask))
     901             :       {
     902           0 :         uel(M0,row_j) ^= uel(M0,row_i);
     903           0 :         uel(M1,row_j) ^= uel(M1,row_i);
     904             :       }
     905             :     }
     906        1728 :     M0[row_i] = M1[row_i] = 0;
     907             :   }
     908        1018 :   mask = 0;
     909       61562 :   for (i = 1; i <= dim; i++)
     910       60544 :     mask |= 1UL<<(s[i]-1);
     911       62454 :   for (i = 1; i <= last_dim; i++)
     912       61436 :     mask |= 1UL<<(last_s[i]-1);
     913        1018 :   if (mask != (ulong)(-1))
     914           2 :     return NULL; /* Failure */
     915        1016 :   setlg(s, dim+1);
     916        1016 :   set_avma(av);
     917        1016 :   *pt_s = s;
     918        1016 :   return M1;
     919             : }
     920             : 
     921             : /* Compute x * A~ */
     922             : static GEN
     923        1176 : F2w_F2Ms_transmul(GEN x, GEN A, long nbrow)
     924             : {
     925        1176 :   long i, j, l = lg(A);
     926        1176 :   GEN b = zero_zv(nbrow);
     927     2103408 :   for (i = 1; i < l; i++)
     928             :   {
     929     2102232 :     GEN c = gel(A,i);
     930     2102232 :     long lc = lg(c);
     931     2102232 :     ulong xi = x[i];
     932    39643268 :     for (j = 1; j < lc; j++)
     933    37541036 :       b[c[j]] ^= xi;
     934             :   }
     935        1176 :   return b;
     936             : }
     937             : 
     938             : /* Compute x * A */
     939             : static GEN
     940        1098 : F2w_F2Ms_mul(GEN x, GEN A)
     941             : {
     942        1098 :   long i, j, l = lg(A);
     943        1098 :   GEN b = cgetg(l, t_VECSMALL);
     944     1978896 :   for (i = 1; i < l; i++)
     945             :   {
     946     1977798 :     GEN c = gel(A,i);
     947     1977798 :     long lc = lg(c);
     948     1977798 :     ulong s = 0;
     949    37322584 :     for (j = 1; j < lc; j++)
     950    35344786 :       s ^= x[c[j]];
     951     1977798 :     b[i] = s;
     952             :   }
     953        1098 :   return b;
     954             : }
     955             : 
     956             : static void
     957        2032 : F2wB_addid_inplace(GEN f)
     958             : {
     959             :   long i;
     960      126320 :   for (i = 1; i <= BIL; i++)
     961      124288 :     uel(f,i) ^= 1UL<<(i-1);
     962        2032 : }
     963             : 
     964             : static void
     965        2032 : F2w_mask_inplace(GEN f, ulong m)
     966             : {
     967        2032 :   long i, l = lg(f);
     968     1914032 :   for (i = 1; i < l; i++)
     969     1912000 :     uel(f,i) &= m;
     970        2032 : }
     971             : 
     972             : static GEN
     973          41 : block_lanczos(GEN B, ulong nbrow)
     974             : {
     975          41 :   pari_sp av = avma, av2;
     976             :   GEN v0, v1, v2, vnext, x, w;
     977             :   GEN winv0, winv1, winv2;
     978             :   GEN vt_a_v0, vt_a_v1, vt_a2_v0, vt_a2_v1;
     979             :   GEN d, e, f, f2, s0;
     980             :   long i, iter;
     981          41 :   long n = lg(B)-1;
     982             :   long dim0;
     983             :   ulong mask0, mask1;
     984          41 :   v1 = zero_zv(n);
     985          41 :   v2 = zero_zv(n);
     986          41 :   vt_a_v1 = zero_zv(BIL);
     987          41 :   vt_a2_v1 = zero_zv(BIL);
     988          41 :   winv1 = zero_zv(BIL);
     989          41 :   winv2 = zero_zv(BIL);
     990          41 :   s0 = identity_zv(BIL);
     991          41 :   mask1 = (ulong)(-1);
     992             : 
     993          41 :   x = random_zv(n);
     994          41 :   w = F2w_F2Ms_mul(F2w_F2Ms_transmul(x, B, nbrow), B);
     995          41 :   v0 = w;
     996          41 :   av2 = avma;
     997          41 :   for (iter=1;;iter++)
     998             :   {
     999        1057 :     vnext = F2w_F2Ms_mul(F2w_F2Ms_transmul(v0, B, nbrow), B);
    1000        1057 :     vt_a_v0  = F2w_transmul(v0, vnext);
    1001        1057 :     if (zv_equal0(vt_a_v0)) break; /* success */
    1002        1018 :     vt_a2_v0 = F2w_transmul(vnext, vnext);
    1003        1018 :     winv0 = find_nonsingular_sub(vt_a_v0, s0, &s0);
    1004        1018 :     if (!winv0) return gc_NULL(av); /* failure */
    1005        1016 :     dim0 = lg(s0)-1;
    1006        1016 :     mask0 = 0;
    1007       61558 :     for (i = 1; i <= dim0; i++)
    1008       60542 :       mask0 |= 1UL<<(s0[i]-1);
    1009        1016 :     d = cgetg(BIL+1, t_VECSMALL);
    1010       63160 :     for (i = 1; i <= BIL; i++)
    1011       62144 :       d[i] = (vt_a2_v0[i] & mask0) ^ vt_a_v0[i];
    1012             : 
    1013        1016 :     d = F2wB_mul(winv0, d);
    1014        1016 :     F2wB_addid_inplace(d);
    1015        1016 :     e = F2wB_mul(winv1, vt_a_v0);
    1016        1016 :     F2w_mask_inplace(e, mask0);
    1017        1016 :     f = F2wB_mul(vt_a_v1, winv1);
    1018        1016 :     F2wB_addid_inplace(f);
    1019        1016 :     f = F2wB_mul(winv2, f);
    1020        1016 :     f2 = cgetg(BIL+1, t_VECSMALL);
    1021       63160 :     for (i = 1; i <= BIL; i++)
    1022       62144 :       f2[i] = ((vt_a2_v1[i] & mask1) ^ vt_a_v1[i]) & mask0;
    1023             : 
    1024        1016 :     f = F2wB_mul(f, f2);
    1025        1016 :     F2w_mask_inplace(vnext, mask0);
    1026        1016 :     F2w_F2wB_mul_add_inplace(v0, d, vnext);
    1027        1016 :     F2w_F2wB_mul_add_inplace(v1, e, vnext);
    1028        1016 :     F2w_F2wB_mul_add_inplace(v2, f, vnext);
    1029        1016 :     d = F2wB_mul(winv0, F2w_transmul(v0, w));
    1030        1016 :     F2w_F2wB_mul_add_inplace(v0, d, x);
    1031        1016 :     v2 = v1; v1 = v0; v0 = vnext;
    1032        1016 :     winv2 = winv1; winv1 = winv0;
    1033        1016 :     vt_a_v1 = vt_a_v0;
    1034        1016 :     vt_a2_v1 = vt_a2_v0;
    1035        1016 :     mask1 = mask0;
    1036        1016 :     gerepileall(av2, 9, &x, &s0, &v0, &v1, &v2,
    1037             :                         &winv1, &winv2, &vt_a_v1, &vt_a2_v1);
    1038             :   }
    1039          39 :   if (DEBUGLEVEL >= 5)
    1040           0 :     err_printf("Lanczos halted after %ld iterations\n", iter);
    1041          39 :   v1 = F2w_F2Ms_transmul(x, B, nbrow);
    1042          39 :   v2 = F2w_F2Ms_transmul(v0, B, nbrow);
    1043          39 :   x  = shallowconcat(F2w_transpose_F2m(x), F2w_transpose_F2m(v0));
    1044          39 :   v1 = shallowconcat(F2w_transpose_F2m(v1), F2w_transpose_F2m(v2));
    1045          39 :   s0 = gel(F2m_indexrank(x), 2);
    1046          39 :   x = shallowextract(x, s0);
    1047          39 :   v1 = shallowextract(v1, s0);
    1048          39 :   return F2m_mul(x, F2m_ker(v1));
    1049             : }
    1050             : 
    1051             : static GEN
    1052        2192 : F2v_inflate(GEN v, GEN p, long n)
    1053             : {
    1054        2192 :   long i, l = lg(p)-1;
    1055        2192 :   GEN w = zero_F2v(n);
    1056     3648838 :   for (i=1; i<=l; i++)
    1057     3646646 :     if (F2v_coeff(v,i))
    1058     1820765 :       F2v_set(w, p[i]);
    1059        2192 :   return w;
    1060             : }
    1061             : 
    1062             : static GEN
    1063          39 : F2m_inflate(GEN x, GEN p, long n)
    1064        2231 : { pari_APPLY_same(F2v_inflate(gel(x,i), p, n)) }
    1065             : 
    1066             : GEN
    1067         382 : F2Ms_ker(GEN M, long nbrow)
    1068             : {
    1069         382 :   pari_sp av = avma;
    1070         382 :   long nbcol = lg(M)-1;
    1071             :   GEN Mp, R, Rp, p;
    1072         382 :   if (nbrow <= 640)
    1073         343 :     return gerepileupto(av, F2m_ker(F2Ms_to_F2m(M, nbrow)));
    1074          39 :   p = F2Ms_colelim(M, nbrow);
    1075          39 :   Mp = vecpermute(M, p);
    1076             :   do
    1077             :   {
    1078          41 :     R = block_lanczos(Mp, nbrow);
    1079          41 :   } while(!R);
    1080          39 :   Rp = F2m_inflate(R, p, nbcol);
    1081          39 :   return gerepilecopy(av, Rp);
    1082             : }
    1083             : 
    1084             : GEN
    1085           0 : F2m_to_F2Ms(GEN M)
    1086             : {
    1087           0 :   long ncol = lg(M)-1;
    1088           0 :   GEN B = cgetg(ncol+1, t_MAT);
    1089             :   long i, j, k;
    1090           0 :   for(i = 1; i <= ncol; i++)
    1091             :   {
    1092           0 :     GEN D, V = gel(M,i);
    1093           0 :     long n = F2v_hamming(V), l = V[1];
    1094           0 :     D = cgetg(n+1, t_VECSMALL);
    1095           0 :     for (j=1, k=1; j<=l; j++)
    1096           0 :       if( F2v_coeff(V,j))
    1097           0 :         D[k++] = j;
    1098           0 :     gel(B, i) = D;
    1099             :   }
    1100           0 :   return B;
    1101             : }
    1102             : 
    1103             : GEN
    1104         343 : F2Ms_to_F2m(GEN M, long nrow)
    1105             : {
    1106         343 :   long i, j, l = lg(M);
    1107         343 :   GEN B = cgetg(l, t_MAT);
    1108       68446 :   for(i = 1; i < l; i++)
    1109             :   {
    1110       68103 :     GEN Bi = zero_F2v(nrow), Mi = gel(M,i);
    1111       68103 :     long l = lg(Mi);
    1112      953176 :     for (j = 1; j < l; j++)
    1113      885073 :       F2v_set(Bi, Mi[j]);
    1114       68103 :     gel(B, i) = Bi;
    1115             :   }
    1116         343 :   return B;
    1117             : }

Generated by: LCOV version 1.13