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 - ZV.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.10.0 lcov report (development 21741-70cf009) Lines: 756 830 91.1 %
Date: 2018-01-21 06:18:30 Functions: 114 120 95.0 %
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             : static int
      18     1215335 : check_ZV(GEN x, long l)
      19             : {
      20             :   long i;
      21     8092101 :   for (i=1; i<l; i++)
      22     6876850 :     if (typ(gel(x,i)) != t_INT) return 0;
      23     1215251 :   return 1;
      24             : }
      25             : void
      26      504456 : RgV_check_ZV(GEN A, const char *s)
      27             : {
      28      504456 :   if (!RgV_is_ZV(A)) pari_err_TYPE(stack_strcat(s," [integer vector]"), A);
      29      504449 : }
      30             : void
      31      317711 : RgM_check_ZM(GEN A, const char *s)
      32             : {
      33      317711 :   long n = lg(A);
      34      317711 :   if (n != 1)
      35             :   {
      36      317676 :     long j, m = lgcols(A);
      37     1532927 :     for (j=1; j<n; j++)
      38     1215335 :       if (!check_ZV(gel(A,j), m))
      39          84 :         pari_err_TYPE(stack_strcat(s," [integer matrix]"), A);
      40             :   }
      41      317627 : }
      42             : 
      43             : static long
      44    35908062 : ZV_max_lg_i(GEN x, long m)
      45             : {
      46    35908062 :   long i, prec = 2;
      47    35908062 :   for (i=1; i<m; i++) { long l = lgefint(gel(x,i)); if (l > prec) prec = l; }
      48    35908062 :   return prec;
      49             : }
      50             : 
      51             : long
      52       10612 : ZV_max_lg(GEN x)
      53       10612 : { return ZV_max_lg_i(x, lg(x)); }
      54             : 
      55             : static long
      56    10201069 : ZM_max_lg_i(GEN x, long n, long m)
      57             : {
      58    10201069 :   long prec = 2;
      59    10201069 :   if (n != 1)
      60             :   {
      61             :     long j;
      62    46098519 :     for (j=1; j<n; j++)
      63             :     {
      64    35897450 :       long l = ZV_max_lg_i(gel(x,j), m);
      65    35897450 :       if (l > prec) prec = l;
      66             :     }
      67             :   }
      68    10201069 :   return prec;
      69             : }
      70             : 
      71             : long
      72        3222 : ZM_max_lg(GEN x)
      73             : {
      74        3222 :   long n = lg(x);
      75        3222 :   if (n==1) return 2;
      76        3222 :   return ZM_max_lg_i(x, n, lgcols(x));
      77             : }
      78             : 
      79             : GEN
      80        3080 : ZM_supnorm(GEN x)
      81             : {
      82        3080 :   long i, j, h, lx = lg(x);
      83        3080 :   GEN s = gen_0;
      84        3080 :   if (lx == 1) return gen_1;
      85        3080 :   h = lgcols(x);
      86       19138 :   for (j=1; j<lx; j++)
      87             :   {
      88       16058 :     GEN xj = gel(x,j);
      89      238854 :     for (i=1; i<h; i++)
      90             :     {
      91      222796 :       GEN c = gel(xj,i);
      92      222796 :       if (abscmpii(c, s) > 0) s = c;
      93             :     }
      94             :   }
      95        3080 :   return absi(s);
      96             : }
      97             : 
      98             : /********************************************************************/
      99             : /**                                                                **/
     100             : /**                           MULTIPLICATION                       **/
     101             : /**                                                                **/
     102             : /********************************************************************/
     103             : /* x non-empty ZM, y a compatible nc (dimension > 0). */
     104             : static GEN
     105      745882 : ZM_nc_mul_i(GEN x, GEN y, long c, long l)
     106             : {
     107             :   long i, j;
     108             :   pari_sp av;
     109      745882 :   GEN z = cgetg(l,t_COL), s;
     110             : 
     111     6541732 :   for (i=1; i<l; i++)
     112             :   {
     113     5795850 :     av = avma; s = muliu(gcoeff(x,i,1),y[1]);
     114   182640026 :     for (j=2; j<c; j++)
     115   176844176 :       if (y[j]) s = addii(s, muliu(gcoeff(x,i,j),y[j]));
     116     5795850 :     gel(z,i) = gerepileuptoint(av,s);
     117             :   }
     118      745882 :   return z;
     119             : }
     120             : 
     121             : /* x ZV, y a compatible zc. */
     122             : GEN
     123        1981 : ZV_zc_mul(GEN x, GEN y)
     124             : {
     125        1981 :   long j, l = lg(x);
     126        1981 :   pari_sp av = avma;
     127        1981 :   GEN s = mulis(gel(x,1),y[1]);
     128      166355 :   for (j=2; j<l; j++)
     129      164374 :     if (y[j]) s = addii(s, mulis(gel(x,j),y[j]));
     130        1981 :   return gerepileuptoint(av,s);
     131             : }
     132             : 
     133             : /* x non-empty ZM, y a compatible zc (dimension > 0). */
     134             : static GEN
     135     2981262 : ZM_zc_mul_i(GEN x, GEN y, long c, long l)
     136             : {
     137             :   long i, j;
     138     2981262 :   GEN z = cgetg(l,t_COL);
     139             : 
     140    17662418 :   for (i=1; i<l; i++)
     141             :   {
     142    14681156 :     pari_sp av = avma;
     143    14681156 :     GEN s = mulis(gcoeff(x,i,1),y[1]);
     144   254198570 :     for (j=2; j<c; j++)
     145   239517414 :       if (y[j]) s = addii(s, mulis(gcoeff(x,i,j),y[j]));
     146    14681156 :     gel(z,i) = gerepileuptoint(av,s);
     147             :   }
     148     2981262 :   return z;
     149             : }
     150             : GEN
     151     1598618 : ZM_zc_mul(GEN x, GEN y) {
     152     1598618 :   long l = lg(x);
     153     1598618 :   if (l == 1) return cgetg(1, t_COL);
     154     1598618 :   return ZM_zc_mul_i(x,y, l, lgcols(x));
     155             : }
     156             : 
     157             : /* y non-empty ZM, x a compatible zv (dimension > 0). */
     158             : GEN
     159        1309 : zv_ZM_mul(GEN x, GEN y) {
     160        1309 :   long i,j, lx = lg(x), ly = lg(y);
     161             :   GEN z;
     162        1309 :   if (lx == 1) return zerovec(ly-1);
     163        1309 :   z = cgetg(ly,t_VEC);
     164        3178 :   for (j=1; j<ly; j++)
     165             :   {
     166        1869 :     pari_sp av = avma;
     167        1869 :     GEN s = mulsi(x[1], gcoeff(y,1,j));
     168        3521 :     for (i=2; i<lx; i++)
     169        1652 :       if (x[i]) s = addii(s, mulsi(x[i], gcoeff(y,i,j)));
     170        1869 :     gel(z,j) = gerepileuptoint(av,s);
     171             :   }
     172        1309 :   return z;
     173             : }
     174             : 
     175             : /* x ZM, y a compatible zm (dimension > 0). */
     176             : GEN
     177      675188 : ZM_zm_mul(GEN x, GEN y)
     178             : {
     179      675188 :   long j, c, l = lg(x), ly = lg(y);
     180      675188 :   GEN z = cgetg(ly, t_MAT);
     181      675188 :   if (l == 1) return z;
     182      675188 :   c = lgcols(x);
     183      675188 :   for (j = 1; j < ly; j++) gel(z,j) = ZM_zc_mul_i(x, gel(y,j), l,c);
     184      675188 :   return z;
     185             : }
     186             : /* x ZM, y a compatible zn (dimension > 0). */
     187             : GEN
     188      693909 : ZM_nm_mul(GEN x, GEN y)
     189             : {
     190      693909 :   long j, c, l = lg(x), ly = lg(y);
     191      693909 :   GEN z = cgetg(ly, t_MAT);
     192      693909 :   if (l == 1) return z;
     193      693909 :   c = lgcols(x);
     194      693909 :   for (j = 1; j < ly; j++) gel(z,j) = ZM_nc_mul_i(x, gel(y,j), l,c);
     195      693909 :   return z;
     196             : }
     197             : 
     198             : /* Strassen-Winograd algorithm */
     199             : 
     200             : /*
     201             :   Return A[ma+1..ma+da, na+1..na+ea] - B[mb+1..mb+db, nb+1..nb+eb]
     202             :   as an (m x n)-matrix, padding the input with zeroes as necessary.
     203             : */
     204             : static GEN
     205       96376 : add_slices(long m, long n,
     206             :            GEN A, long ma, long da, long na, long ea,
     207             :            GEN B, long mb, long db, long nb, long eb)
     208             : {
     209       96376 :   long min_d = minss(da, db), min_e = minss(ea, eb), i, j;
     210       96376 :   GEN M = cgetg(n + 1, t_MAT), C;
     211             : 
     212     1921894 :   for (j = 1; j <= min_e; j++) {
     213     1825518 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     214    47557127 :     for (i = 1; i <= min_d; i++)
     215    91463218 :       gel(C, i) = addii(gcoeff(A, ma + i, na + j),
     216    45731609 :                         gcoeff(B, mb + i, nb + j));
     217     1955480 :     for (; i <= da; i++)
     218      129962 :       gel(C, i) = gcoeff(A, ma + i, na + j);
     219     1825518 :     for (; i <= db; i++)
     220           0 :       gel(C, i) = gcoeff(B, mb + i, nb + j);
     221     1825518 :     for (; i <= m; i++)
     222           0 :       gel(C, i) = gen_0;
     223             :   }
     224      110223 :   for (; j <= ea; j++) {
     225       13847 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     226      277174 :     for (i = 1; i <= da; i++)
     227      263327 :       gel(C, i) = gcoeff(A, ma + i, na + j);
     228       13847 :     for (; i <= m; i++)
     229           0 :       gel(C, i) = gen_0;
     230             :   }
     231       96376 :   for (; j <= eb; j++) {
     232           0 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     233           0 :     for (i = 1; i <= db; i++)
     234           0 :       gel(C, i) = gcoeff(B, mb + i, nb + j);
     235           0 :     for (; i <= m; i++)
     236           0 :       gel(C, i) = gen_0;
     237             :   }
     238       96376 :   for (; j <= n; j++)
     239           0 :     gel(M, j) = zerocol(m);
     240       96376 :   return M;
     241             : }
     242             : 
     243             : /*
     244             :   Return A[ma+1..ma+da, na+1..na+ea] - B[mb+1..mb+db, nb+1..nb+eb]
     245             :   as an (m x n)-matrix, padding the input with zeroes as necessary.
     246             : */
     247             : static GEN
     248       84329 : subtract_slices(long m, long n,
     249             :                 GEN A, long ma, long da, long na, long ea,
     250             :                 GEN B, long mb, long db, long nb, long eb)
     251             : {
     252       84329 :   long min_d = minss(da, db), min_e = minss(ea, eb), i, j;
     253       84329 :   GEN M = cgetg(n + 1, t_MAT), C;
     254             : 
     255     1847635 :   for (j = 1; j <= min_e; j++) {
     256     1763306 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     257    51054988 :     for (i = 1; i <= min_d; i++)
     258    98583364 :       gel(C, i) = subii(gcoeff(A, ma + i, na + j),
     259    49291682 :                         gcoeff(B, mb + i, nb + j));
     260     1889468 :     for (; i <= da; i++)
     261      126162 :       gel(C, i) = gcoeff(A, ma + i, na + j);
     262     1952212 :     for (; i <= db; i++)
     263      188906 :       gel(C, i) = negi(gcoeff(B, mb + i, nb + j));
     264     1763306 :     for (; i <= m; i++)
     265           0 :       gel(C, i) = gen_0;
     266             :   }
     267       84329 :   for (; j <= ea; j++) {
     268           0 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     269           0 :     for (i = 1; i <= da; i++)
     270           0 :       gel(C, i) = gcoeff(A, ma + i, na + j);
     271           0 :     for (; i <= m; i++)
     272           0 :       gel(C, i) = gen_0;
     273             :   }
     274       91553 :   for (; j <= eb; j++) {
     275        7224 :     gel(M, j) = C = cgetg(m + 1, t_COL);
     276      197438 :     for (i = 1; i <= db; i++)
     277      190214 :       gel(C, i) = negi(gcoeff(B, mb + i, nb + j));
     278        7224 :     for (; i <= m; i++)
     279           0 :       gel(C, i) = gen_0;
     280             :   }
     281       91553 :   for (; j <= n; j++)
     282        7224 :     gel(M, j) = zerocol(m);
     283       84329 :   return M;
     284             : }
     285             : 
     286             : static GEN ZM_mul_i(GEN x, GEN y, long l, long lx, long ly);
     287             : 
     288             : /* Strassen-Winograd matrix product A (m x n) * B (n x p) */
     289             : static GEN
     290       12047 : ZM_mul_sw(GEN A, GEN B, long m, long n, long p)
     291             : {
     292       12047 :   pari_sp av = avma;
     293       12047 :   long m1 = (m + 1)/2, m2 = m/2,
     294       12047 :     n1 = (n + 1)/2, n2 = n/2,
     295       12047 :     p1 = (p + 1)/2, p2 = p/2;
     296             :   GEN A11, A12, A22, B11, B21, B22,
     297             :     S1, S2, S3, S4, T1, T2, T3, T4,
     298             :     M1, M2, M3, M4, M5, M6, M7,
     299             :     V1, V2, V3, C11, C12, C21, C22, C;
     300             : 
     301       12047 :   T2 = subtract_slices(n1, p2, B, 0, n1, p1, p2, B, n1, n2, p1, p2);
     302       12047 :   S1 = subtract_slices(m2, n1, A, m1, m2, 0, n1, A, 0, m2, 0, n1);
     303       12047 :   M2 = ZM_mul_i(S1, T2, m2 + 1, n1 + 1, p2 + 1);
     304       12047 :   if (gc_needed(av, 1))
     305           0 :     gerepileall(av, 2, &T2, &M2);  /* destroy S1 */
     306       12047 :   T3 = subtract_slices(n1, p1, T2, 0, n1, 0, p2, B, 0, n1, 0, p1);
     307       12047 :   if (gc_needed(av, 1))
     308           0 :     gerepileall(av, 2, &M2, &T3);  /* destroy T2 */
     309       12047 :   S2 = add_slices(m2, n1, A, m1, m2, 0, n1, A, m1, m2, n1, n2);
     310       12047 :   T1 = subtract_slices(n1, p1, B, 0, n1, p1, p2, B, 0, n1, 0, p2);
     311       12047 :   M3 = ZM_mul_i(S2, T1, m2 + 1, n1 + 1, p2 + 1);
     312       12047 :   if (gc_needed(av, 1))
     313           0 :     gerepileall(av, 4, &M2, &T3, &S2, &M3);  /* destroy T1 */
     314       12047 :   S3 = subtract_slices(m1, n1, S2, 0, m2, 0, n1, A, 0, m1, 0, n1);
     315       12047 :   if (gc_needed(av, 1))
     316           0 :     gerepileall(av, 4, &M2, &T3, &M3, &S3);  /* destroy S2 */
     317       12047 :   A11 = matslice(A, 1, m1, 1, n1);
     318       12047 :   B11 = matslice(B, 1, n1, 1, p1);
     319       12047 :   M1 = ZM_mul_i(A11, B11, m1 + 1, n1 + 1, p1 + 1);
     320       12047 :   if (gc_needed(av, 1))
     321           0 :     gerepileall(av, 5, &M2, &T3, &M3, &S3, &M1);  /* destroy A11, B11 */
     322       12047 :   A12 = matslice(A, 1, m1, n1 + 1, n);
     323       12047 :   B21 = matslice(B, n1 + 1, n, 1, p1);
     324       12047 :   M4 = ZM_mul_i(A12, B21, m1 + 1, n2 + 1, p1 + 1);
     325       12047 :   if (gc_needed(av, 1))
     326           0 :     gerepileall(av, 6, &M2, &T3, &M3, &S3, &M1, &M4);  /* destroy A12, B21 */
     327       12047 :   C11 = add_slices(m1, p1, M1, 0, m1, 0, p1, M4, 0, m1, 0, p1);
     328       12047 :   if (gc_needed(av, 1))
     329           0 :     gerepileall(av, 6, &M2, &T3, &M3, &S3, &M1, &C11);  /* destroy M4 */
     330       12047 :   M5 = ZM_mul_i(S3, T3, m1 + 1, n1 + 1, p1 + 1);
     331       12047 :   S4 = subtract_slices(m1, n2, A, 0, m1, n1, n2, S3, 0, m1, 0, n2);
     332       12047 :   if (gc_needed(av, 1))
     333           6 :     gerepileall(av, 7, &M2, &T3, &M3, &M1, &C11, &M5, &S4);  /* destroy S3 */
     334       12047 :   T4 = add_slices(n2, p1, B, n1, n2, 0, p1, T3, 0, n2, 0, p1);
     335       12047 :   if (gc_needed(av, 1))
     336           0 :     gerepileall(av, 7, &M2, &M3, &M1, &C11, &M5, &S4, &T4);  /* destroy T3 */
     337       12047 :   V1 = subtract_slices(m1, p1, M1, 0, m1, 0, p1, M5, 0, m1, 0, p1);
     338       12047 :   if (gc_needed(av, 1))
     339           0 :     gerepileall(av, 6, &M2, &M3, &S4, &T4, &C11, &V1);  /* destroy M1, M5 */
     340       12047 :   B22 = matslice(B, n1 + 1, n, p1 + 1, p);
     341       12047 :   M6 = ZM_mul_i(S4, B22, m1 + 1, n2 + 1, p2 + 1);
     342       12047 :   if (gc_needed(av, 1))
     343           0 :     gerepileall(av, 6, &M2, &M3, &T4, &C11, &V1, &M6);  /* destroy S4, B22 */
     344       12047 :   A22 = matslice(A, m1 + 1, m, n1 + 1, n);
     345       12047 :   M7 = ZM_mul_i(A22, T4, m2 + 1, n2 + 1, p1 + 1);
     346       12047 :   if (gc_needed(av, 1))
     347           0 :     gerepileall(av, 6, &M2, &M3, &C11, &V1, &M6, &M7);  /* destroy A22, T4 */
     348       12047 :   V3 = add_slices(m1, p2, V1, 0, m1, 0, p2, M3, 0, m2, 0, p2);
     349       12047 :   C12 = add_slices(m1, p2, V3, 0, m1, 0, p2, M6, 0, m1, 0, p2);
     350       12047 :   if (gc_needed(av, 1))
     351           6 :     gerepileall(av, 6, &M2, &M3, &C11, &V1, &M7, &C12);  /* destroy V3, M6 */
     352       12047 :   V2 = add_slices(m2, p1, V1, 0, m2, 0, p1, M2, 0, m2, 0, p2);
     353       12047 :   if (gc_needed(av, 1))
     354           0 :     gerepileall(av, 5, &M3, &C11, &M7, &C12, &V2);  /* destroy V1, M2 */
     355       12047 :   C21 = add_slices(m2, p1, V2, 0, m2, 0, p1, M7, 0, m2, 0, p1);
     356       12047 :   if (gc_needed(av, 1))
     357           0 :     gerepileall(av, 5, &M3, &C11, &C12, &V2, &C21);  /* destroy M7 */
     358       12047 :   C22 = add_slices(m2, p2, V2, 0, m2, 0, p2, M3, 0, m2, 0, p2);
     359       12047 :   if (gc_needed(av, 1))
     360           0 :     gerepileall(av, 4, &C11, &C12, &C21, &C22);  /* destroy V2, M3 */
     361       12047 :   C = mkmat2(mkcol2(C11, C21), mkcol2(C12, C22));
     362       12047 :   return gerepilecopy(av, shallowmatconcat(C));
     363             : }
     364             : 
     365             : /* x[i,]*y. Assume lg(x) > 1 and 0 < i < lgcols(x) */
     366             : static GEN
     367   145461449 : ZMrow_ZC_mul_i(GEN x, GEN y, long i, long lx)
     368             : {
     369   145461449 :   pari_sp av = avma;
     370   145461449 :   GEN c = mulii(gcoeff(x,i,1), gel(y,1)), ZERO = gen_0;
     371             :   long k;
     372  1971740603 :   for (k = 2; k < lx; k++)
     373             :   {
     374  1826279154 :     GEN t = mulii(gcoeff(x,i,k), gel(y,k));
     375  1826279154 :     if (t != ZERO) c = addii(c, t);
     376             :   }
     377   145461449 :   return gerepileuptoint(av, c);
     378             : }
     379             : GEN
     380    25447200 : ZMrow_ZC_mul(GEN x, GEN y, long i)
     381    25447200 : { return ZMrow_ZC_mul_i(x, y, i, lg(x)); }
     382             : 
     383             : /* return x * y, 1 < lx = lg(x), l = lgcols(x) */
     384             : static GEN
     385    22419852 : ZM_ZC_mul_i(GEN x, GEN y, long lx, long l)
     386             : {
     387    22419852 :   GEN z = cgetg(l,t_COL);
     388             :   long i;
     389    22419852 :   for (i=1; i<l; i++) gel(z,i) = ZMrow_ZC_mul_i(x,y,i,lx);
     390    22419852 :   return z;
     391             : }
     392             : 
     393             : static GEN
     394     5088672 : ZM_mul_classical(GEN x, GEN y, long l, long lx, long ly)
     395             : {
     396             :   long j;
     397     5088672 :   GEN z = cgetg(ly, t_MAT);
     398    22095654 :   for (j = 1; j < ly; j++)
     399    17006982 :     gel(z, j) = ZM_ZC_mul_i(x, gel(y, j), lx, l);
     400     5088672 :   return z;
     401             : }
     402             : 
     403             : /* Strassen-Winograd used for dim >= ZM_sw_bound */
     404             : static GEN
     405     5097128 : ZM_mul_i(GEN x, GEN y, long l, long lx, long ly)
     406             : {
     407     5097128 :   long s = maxss(ZM_max_lg_i(x,lx,l), ZM_max_lg_i(y,ly,lx));
     408     5097128 :   long ZM_sw_bound = s > 60 ? 2: s > 25 ? 4: s>15 ? 8 : s > 8 ? 16 : 32;
     409     5097128 :   if (l <= ZM_sw_bound || lx <= ZM_sw_bound || ly <= ZM_sw_bound)
     410     5085088 :     return ZM_mul_classical(x, y, l, lx, ly);
     411             :   else
     412       12040 :     return ZM_mul_sw(x, y, l - 1, lx - 1, ly - 1);
     413             : }
     414             : 
     415             : GEN
     416     5072432 : ZM_mul(GEN x, GEN y)
     417             : {
     418     5072432 :   long lx=lg(x), ly=lg(y);
     419     5072432 :   if (ly==1) return cgetg(1,t_MAT);
     420     5013723 :   if (lx==1) return zeromat(0, ly-1);
     421     5012799 :   return ZM_mul_i(x, y, lgcols(x), lx, ly);
     422             : }
     423             : 
     424             : GEN
     425      246347 : QM_mul(GEN x, GEN y)
     426             : {
     427      246347 :   GEN dx, nx = Q_primitive_part(x, &dx);
     428      246347 :   GEN dy, ny = Q_primitive_part(y, &dy);
     429      246347 :   GEN z = ZM_mul(nx, ny);
     430      246347 :   if (dx || dy)
     431             :   {
     432      189937 :     GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
     433      189937 :     if (!gequal1(d)) z = ZM_Q_mul(z, d);
     434             :   }
     435      246347 :   return z;
     436             : }
     437             : 
     438             : /* assume result is symmetric */
     439             : GEN
     440           0 : ZM_multosym(GEN x, GEN y)
     441             : {
     442           0 :   long j, lx, ly = lg(y);
     443             :   GEN M;
     444           0 :   if (ly == 1) return cgetg(1,t_MAT);
     445           0 :   lx = lg(x); /* = lgcols(y) */
     446           0 :   if (lx == 1) return cgetg(1,t_MAT);
     447             :   /* ly = lgcols(x) */
     448           0 :   M = cgetg(ly, t_MAT);
     449           0 :   for (j=1; j<ly; j++)
     450             :   {
     451           0 :     GEN z = cgetg(ly,t_COL), yj = gel(y,j);
     452             :     long i;
     453           0 :     for (i=1; i<j; i++) gel(z,i) = gcoeff(M,j,i);
     454           0 :     for (i=j; i<ly; i++)gel(z,i) = ZMrow_ZC_mul_i(x,yj,i,lx);
     455           0 :     gel(M,j) = z;
     456             :   }
     457           0 :   return M;
     458             : }
     459             : 
     460             : /* compute m*diagonal(d), assume lg(d) = lg(m). Shallow */
     461             : GEN
     462          14 : ZM_mul_diag(GEN m, GEN d)
     463             : {
     464             :   long j, l;
     465          14 :   GEN y = cgetg_copy(m, &l);
     466          56 :   for (j=1; j<l; j++)
     467             :   {
     468          42 :     GEN c = gel(d,j);
     469          42 :     gel(y,j) = equali1(c)? gel(m,j): ZC_Z_mul(gel(m,j), c);
     470             :   }
     471          14 :   return y;
     472             : }
     473             : /* compute diagonal(d)*m, assume lg(d) = lg(m~). Shallow */
     474             : GEN
     475       86153 : ZM_diag_mul(GEN d, GEN m)
     476             : {
     477       86153 :   long i, j, l = lg(d), lm = lg(m);
     478       86153 :   GEN y = cgetg(lm, t_MAT);
     479       86153 :   for (j=1; j<lm; j++) gel(y,j) = cgetg(l, t_COL);
     480      278154 :   for (i=1; i<l; i++)
     481             :   {
     482      192001 :     GEN c = gel(d,i);
     483      192001 :     if (equali1(c))
     484        9058 :       for (j=1; j<lm; j++) gcoeff(y,i,j) = gcoeff(m,i,j);
     485             :     else
     486      182943 :       for (j=1; j<lm; j++) gcoeff(y,i,j) = mulii(gcoeff(m,i,j), c);
     487             :   }
     488       86153 :   return y;
     489             : }
     490             : 
     491             : /* assume lx > 1 is lg(x) = lg(y) */
     492             : static GEN
     493    12047146 : ZV_dotproduct_i(GEN x, GEN y, long lx)
     494             : {
     495    12047146 :   pari_sp av = avma;
     496    12047146 :   GEN c = mulii(gel(x,1), gel(y,1));
     497             :   long i;
     498    96831840 :   for (i = 2; i < lx; i++)
     499             :   {
     500    84784694 :     GEN t = mulii(gel(x,i), gel(y,i));
     501    84784694 :     if (t != gen_0) c = addii(c, t);
     502             :   }
     503    12047146 :   return gerepileuptoint(av, c);
     504             : }
     505             : 
     506             : /* x~ * y, assuming result is symmetric */
     507             : GEN
     508      160293 : ZM_transmultosym(GEN x, GEN y)
     509             : {
     510      160293 :   long i, j, l, ly = lg(y);
     511             :   GEN M;
     512      160293 :   if (ly == 1) return cgetg(1,t_MAT);
     513             :   /* lg(x) = ly */
     514      160293 :   l = lgcols(y); /* = lgcols(x) */
     515      160293 :   M = cgetg(ly, t_MAT);
     516     1251162 :   for (i=1; i<ly; i++)
     517             :   {
     518     1090869 :     GEN xi = gel(x,i), c = cgetg(ly,t_COL);
     519     1090869 :     gel(M,i) = c;
     520     4913427 :     for (j=1; j<i; j++)
     521     3822558 :       gcoeff(M,i,j) = gel(c,j) = ZV_dotproduct_i(xi,gel(y,j),l);
     522     1090869 :     gel(c,i) = ZV_dotproduct_i(xi,gel(y,i),l);
     523             :   }
     524      160293 :   return M;
     525             : }
     526             : /* x~ * y */
     527             : GEN
     528         497 : ZM_transmul(GEN x, GEN y)
     529             : {
     530         497 :   long i, j, l, lx, ly = lg(y);
     531             :   GEN M;
     532         497 :   if (ly == 1) return cgetg(1,t_MAT);
     533         497 :   lx = lg(x);
     534         497 :   l = lgcols(y);
     535         497 :   if (lgcols(x) != l) pari_err_OP("operation 'ZM_transmul'", x,y);
     536         497 :   M = cgetg(ly, t_MAT);
     537        1617 :   for (i=1; i<ly; i++)
     538             :   {
     539        1120 :     GEN yi = gel(y,i), c = cgetg(lx,t_COL);
     540        1120 :     gel(M,i) = c;
     541        1120 :     for (j=1; j<lx; j++) gel(c,j) = ZV_dotproduct_i(yi,gel(x,j),l);
     542             :   }
     543         497 :   return M;
     544             : }
     545             : 
     546             : static GEN
     547        3591 : ZM_sqr_i(GEN x, long l)
     548             : {
     549        3591 :   long s = ZM_max_lg_i(x,l,l);
     550        3591 :   long ZM_sw_bound = s > 60 ? 2: s > 25 ? 4: s>15 ? 8 : s > 8 ? 16 : 32;
     551        3591 :   if (l <= ZM_sw_bound)
     552        3584 :     return ZM_mul_classical(x, x, l, l, l);
     553             :   else
     554           7 :     return ZM_mul_sw(x, x, l - 1, l - 1, l - 1);
     555             : }
     556             : 
     557             : GEN
     558        3591 : ZM_sqr(GEN x)
     559             : {
     560        3591 :   long lx=lg(x);
     561        3591 :   if (lx==1) return cgetg(1,t_MAT);
     562        3591 :   return ZM_sqr_i(x, lx);
     563             : }
     564             : GEN
     565     5504074 : ZM_ZC_mul(GEN x, GEN y)
     566             : {
     567     5504074 :   long lx = lg(x);
     568     5504074 :   return lx==1? cgetg(1,t_COL): ZM_ZC_mul_i(x, y, lx, lgcols(x));
     569             : }
     570             : 
     571             : GEN
     572     1101968 : ZC_Z_div(GEN x, GEN c)
     573     1101968 : { pari_APPLY_type(t_COL, Qdivii(gel(x,i), c)) }
     574             : 
     575             : GEN
     576       12194 : ZM_Z_div(GEN x, GEN c)
     577       12194 : { pari_APPLY_same(ZC_Z_div(gel(x, i), c)) }
     578             : 
     579             : GEN
     580     1181294 : ZC_Q_mul(GEN A, GEN z)
     581             : {
     582     1181294 :   pari_sp av = avma;
     583     1181294 :   long i, l = lg(A);
     584             :   GEN d, n, Ad, B, u;
     585     1181294 :   if (typ(z)==t_INT) return ZC_Z_mul(A,z);
     586     1181035 :   n = gel(z, 1); d = gel(z, 2);
     587     1181035 :   Ad = FpC_red(A, d);
     588     1181035 :   u = gcdii(d, FpV_factorback(Ad, NULL, d));
     589     1181035 :   B = cgetg(l, t_COL);
     590     1181035 :   if (equali1(u))
     591             :   {
     592      183759 :     for(i=1; i<l; i++)
     593      148510 :       gel(B, i) = mkfrac(mulii(n, gel(A,i)), d);
     594             :   } else
     595             :   {
     596    11535598 :     for(i=1; i<l; i++)
     597             :     {
     598    10389812 :       GEN di = gcdii(gel(Ad, i), u), ni = mulii(n, diviiexact(gel(A,i), di));
     599    10389812 :       if (equalii(d, di))
     600     7796715 :         gel(B, i) = ni;
     601             :       else
     602     2593097 :         gel(B, i) = mkfrac(ni, diviiexact(d, di));
     603             :     }
     604             :   }
     605     1181035 :   return gerepilecopy(av, B);
     606             : }
     607             : 
     608             : GEN
     609      299043 : ZM_Q_mul(GEN x, GEN z)
     610             : {
     611      299043 :   if (typ(z)==t_INT) return ZM_Z_mul(x,z);
     612      197661 :   pari_APPLY_same(ZC_Q_mul(gel(x, i), z));
     613             : }
     614             : 
     615             : long
     616   151102651 : zv_dotproduct(GEN x, GEN y)
     617             : {
     618   151102651 :   long i, lx = lg(x);
     619             :   ulong c;
     620   151102651 :   if (lx == 1) return 0;
     621   151102651 :   c = uel(x,1)*uel(y,1);
     622  2329284041 :   for (i = 2; i < lx; i++)
     623  2178181390 :     c += uel(x,i)*uel(y,i);
     624   151102651 :   return c;
     625             : }
     626             : 
     627             : GEN
     628      144858 : ZV_ZM_mul(GEN x, GEN y)
     629             : {
     630      144858 :   long i, lx = lg(x), ly = lg(y);
     631             :   GEN z;
     632      144858 :   if (lx == 1) return zerovec(ly-1);
     633      144844 :   z = cgetg(ly, t_VEC);
     634      144844 :   for (i = 1; i < ly; i++) gel(z,i) = ZV_dotproduct_i(x, gel(y,i), lx);
     635      144844 :   return z;
     636             : }
     637             : 
     638             : GEN
     639           0 : ZC_ZV_mul(GEN x, GEN y)
     640             : {
     641           0 :   long i,j, lx=lg(x), ly=lg(y);
     642             :   GEN z;
     643           0 :   if (ly==1) return cgetg(1,t_MAT);
     644           0 :   z = cgetg(ly,t_MAT);
     645           0 :   for (j=1; j < ly; j++)
     646             :   {
     647           0 :     gel(z,j) = cgetg(lx,t_COL);
     648           0 :     for (i=1; i<lx; i++) gcoeff(z,i,j) = mulii(gel(x,i),gel(y,j));
     649             :   }
     650           0 :   return z;
     651             : }
     652             : 
     653             : GEN
     654     3800000 : ZV_dotsquare(GEN x)
     655             : {
     656             :   long i, lx;
     657             :   pari_sp av;
     658             :   GEN z;
     659     3800000 :   lx = lg(x);
     660     3800000 :   if (lx == 1) return gen_0;
     661     3800000 :   av = avma; z = sqri(gel(x,1));
     662     3800000 :   for (i=2; i<lx; i++) z = addii(z, sqri(gel(x,i)));
     663     3800000 :   return gerepileuptoint(av,z);
     664             : }
     665             : 
     666             : GEN
     667     9364840 : ZV_dotproduct(GEN x,GEN y)
     668             : {
     669             :   long lx;
     670     9364840 :   if (x == y) return ZV_dotsquare(x);
     671     6721412 :   lx = lg(x);
     672     6721412 :   if (lx == 1) return gen_0;
     673     6721412 :   return ZV_dotproduct_i(x, y, lx);
     674             : }
     675             : 
     676             : static GEN
     677         280 : _ZM_mul(void *data /*ignored*/, GEN x, GEN y)
     678         280 : { (void)data; return ZM_mul(x,y); }
     679             : static GEN
     680        2604 : _ZM_sqr(void *data /*ignored*/, GEN x)
     681        2604 : { (void)data; return ZM_sqr(x); }
     682             : GEN
     683           0 : ZM_pow(GEN x, GEN n)
     684             : {
     685           0 :   pari_sp av = avma;
     686           0 :   if (!signe(n)) return matid(lg(x)-1);
     687           0 :   return gerepileupto(av, gen_pow(x, n, NULL, &_ZM_sqr, &_ZM_mul));
     688             : }
     689             : GEN
     690        2170 : ZM_powu(GEN x, ulong n)
     691             : {
     692        2170 :   pari_sp av = avma;
     693        2170 :   if (!n) return matid(lg(x)-1);
     694        2170 :   return gerepileupto(av, gen_powu(x, n, NULL, &_ZM_sqr, &_ZM_mul));
     695             : }
     696             : /********************************************************************/
     697             : /**                                                                **/
     698             : /**                           ADD, SUB                             **/
     699             : /**                                                                **/
     700             : /********************************************************************/
     701             : static GEN
     702    12967405 : ZC_add_i(GEN x, GEN y, long lx)
     703             : {
     704    12967405 :   GEN A = cgetg(lx, t_COL);
     705             :   long i;
     706    12967413 :   for (i=1; i<lx; i++) gel(A,i) = addii(gel(x,i), gel(y,i));
     707    12967411 :   return A;
     708             : }
     709             : GEN
     710     6886076 : ZC_add(GEN x, GEN y) { return ZC_add_i(x, y, lg(x)); }
     711             : GEN
     712      103017 : ZC_Z_add(GEN x, GEN y)
     713             : {
     714      103017 :   long k, lx = lg(x);
     715      103017 :   GEN z = cgetg(lx, t_COL);
     716      103017 :   if (lx == 1) pari_err_TYPE2("+",x,y);
     717      103017 :   gel(z,1) = addii(y,gel(x,1));
     718      103017 :   for (k = 2; k < lx; k++) gel(z,k) = icopy(gel(x,k));
     719      103017 :   return z;
     720             : }
     721             : 
     722             : static GEN
     723     2010755 : ZC_sub_i(GEN x, GEN y, long lx)
     724             : {
     725             :   long i;
     726     2010755 :   GEN A = cgetg(lx, t_COL);
     727     2010755 :   for (i=1; i<lx; i++) gel(A,i) = subii(gel(x,i), gel(y,i));
     728     2010755 :   return A;
     729             : }
     730             : GEN
     731     1194561 : ZC_sub(GEN x, GEN y) { return ZC_sub_i(x, y, lg(x)); }
     732             : GEN
     733           0 : ZC_Z_sub(GEN x, GEN y)
     734             : {
     735           0 :   long k, lx = lg(x);
     736           0 :   GEN z = cgetg(lx, t_COL);
     737           0 :   if (lx == 1) pari_err_TYPE2("+",x,y);
     738           0 :   gel(z,1) = subii(gel(x,1), y);
     739           0 :   for (k = 2; k < lx; k++) gel(z,k) = icopy(gel(x,k));
     740           0 :   return z;
     741             : }
     742             : GEN
     743       94663 : Z_ZC_sub(GEN a, GEN x)
     744             : {
     745       94663 :   long k, lx = lg(x);
     746       94663 :   GEN z = cgetg(lx, t_COL);
     747       94663 :   if (lx == 1) pari_err_TYPE2("-",a,x);
     748       94663 :   gel(z,1) = subii(a, gel(x,1));
     749       94663 :   for (k = 2; k < lx; k++) gel(z,k) = negi(gel(x,k));
     750       94663 :   return z;
     751             : }
     752             : 
     753             : GEN
     754     1156242 : ZM_add(GEN x, GEN y)
     755             : {
     756     1156242 :   long lx = lg(x), l, j;
     757             :   GEN z;
     758     1156242 :   if (lx == 1) return cgetg(1, t_MAT);
     759     1149928 :   z = cgetg(lx, t_MAT); l = lgcols(x);
     760     1149928 :   for (j = 1; j < lx; j++) gel(z,j) = ZC_add_i(gel(x,j), gel(y,j), l);
     761     1149928 :   return z;
     762             : }
     763             : GEN
     764      711475 : ZM_sub(GEN x, GEN y)
     765             : {
     766      711475 :   long lx = lg(x), l, j;
     767             :   GEN z;
     768      711475 :   if (lx == 1) return cgetg(1, t_MAT);
     769      711475 :   z = cgetg(lx, t_MAT); l = lgcols(x);
     770      711475 :   for (j = 1; j < lx; j++) gel(z,j) = ZC_sub_i(gel(x,j), gel(y,j), l);
     771      711475 :   return z;
     772             : }
     773             : /********************************************************************/
     774             : /**                                                                **/
     775             : /**                         LINEAR COMBINATION                     **/
     776             : /**                                                                **/
     777             : /********************************************************************/
     778             : /* return X/c assuming division is exact */
     779             : GEN
     780     2225946 : ZC_Z_divexact(GEN x, GEN c)
     781     2225946 : { pari_APPLY_type(t_COL, diviiexact(gel(x,i), c)) }
     782             : 
     783             : GEN
     784      846978 : ZM_Z_divexact(GEN x, GEN c)
     785      846978 : { pari_APPLY_same(ZC_Z_divexact(gel(x,i), c)) }
     786             : 
     787             : GEN
     788    10486127 : ZC_Z_mul(GEN x, GEN c)
     789             : {
     790    10486127 :   if (!signe(c)) return zerocol(lg(x)-1);
     791    10186843 :   if (is_pm1(c)) return (signe(c) > 0)? ZC_copy(x): ZC_neg(x);
     792     8245537 :   pari_APPLY_type(t_COL, mulii(gel(x,i), c))
     793             : }
     794             : 
     795             : GEN
     796       15064 : ZC_z_mul(GEN x, long c)
     797             : {
     798       15064 :   if (!c) return zerocol(lg(x)-1);
     799        8449 :   if (c == 1) return ZC_copy(x);
     800        3577 :   if (c ==-1) return ZC_neg(x);
     801        2121 :   pari_APPLY_type(t_COL, mulsi(c, gel(x,i)))
     802             : }
     803             : 
     804             : GEN
     805       10746 : zv_z_mul(GEN x, long n)
     806       10746 : { pari_APPLY_long(x[i]*n) }
     807             : 
     808             : /* return a ZM */
     809             : GEN
     810      693909 : nm_Z_mul(GEN X, GEN c)
     811             : {
     812      693909 :   long i, j, h, l = lg(X), s = signe(c);
     813             :   GEN A;
     814      693909 :   if (l == 1) return cgetg(1, t_MAT);
     815      693909 :   h = lgcols(X);
     816      693909 :   if (!s) return zeromat(h-1, l-1);
     817      693909 :   if (is_pm1(c)) {
     818           0 :     if (s > 0) return Flm_to_ZM(X);
     819           0 :     X = Flm_to_ZM(X); ZM_togglesign(X); return X;
     820             :   }
     821      693909 :   A = cgetg(l, t_MAT);
     822     1439791 :   for (j = 1; j < l; j++)
     823             :   {
     824      745882 :     GEN a = cgetg(h, t_COL), x = gel(X, j);
     825      745882 :     for (i = 1; i < h; i++) gel(a,i) = muliu(c, x[i]);
     826      745882 :     gel(A,j) = a;
     827             :   }
     828      693909 :   return A;
     829             : }
     830             : GEN
     831      746501 : ZM_Z_mul(GEN X, GEN c)
     832             : {
     833      746501 :   long i, j, h, l = lg(X);
     834             :   GEN A;
     835      746501 :   if (l == 1) return cgetg(1, t_MAT);
     836      731493 :   h = lgcols(X);
     837      731493 :   if (!signe(c)) return zeromat(h-1, l-1);
     838      730065 :   if (is_pm1(c)) return (signe(c) > 0)? ZM_copy(X): ZM_neg(X);
     839      557084 :   A = cgetg(l, t_MAT);
     840     6138385 :   for (j = 1; j < l; j++)
     841             :   {
     842     5581301 :     GEN a = cgetg(h, t_COL), x = gel(X, j);
     843     5581301 :     for (i = 1; i < h; i++) gel(a,i) = mulii(c, gel(x,i));
     844     5581301 :     gel(A,j) = a;
     845             :   }
     846      557084 :   return A;
     847             : }
     848             : 
     849             : /* X <- X + v Y (elementary col operation) */
     850             : void
     851    38775284 : ZC_lincomb1_inplace(GEN X, GEN Y, GEN v)
     852             : {
     853    38775284 :   long i, m = lgefint(v);
     854    77550568 :   if (m == 2) return; /* v = 0 */
     855    38775284 :   for (i = lg(X)-1; i; i--) gel(X,i) = addmulii_inplace(gel(X,i), gel(Y,i), v);
     856             : }
     857             : void
     858     8549483 : Flc_lincomb1_inplace(GEN X, GEN Y, ulong v, ulong q)
     859             : {
     860             :   long i;
     861    17098966 :   if (!v) return; /* v = 0 */
     862     8549483 :   for (i = lg(X)-1; i; i--) X[i] = Fl_add(X[i], Fl_mul(Y[i], v, q), q);
     863             : }
     864             : 
     865             : /* X + v Y, wasteful if (v = 0) */
     866             : static GEN
     867     2102499 : ZC_lincomb1(GEN v, GEN x, GEN y)
     868     2102499 : { pari_APPLY_type(t_COL, addmulii(gel(x,i), gel(y,i), v)) }
     869             : 
     870             : /* -X + vY */
     871             : static GEN
     872      382243 : ZC_lincomb_1(GEN v, GEN x, GEN y)
     873      382243 : { pari_APPLY_type(t_COL, mulsubii(gel(y,i), v, gel(x,i))) }
     874             : 
     875             : /* X,Y compatible ZV; u,v in Z. Returns A = u*X + v*Y */
     876             : GEN
     877     9042596 : ZC_lincomb(GEN u, GEN v, GEN X, GEN Y)
     878             : {
     879             :   long su, sv;
     880             :   GEN A;
     881             : 
     882     9042596 :   su = signe(u); if (!su) return ZC_Z_mul(Y, v);
     883     9042400 :   sv = signe(v); if (!sv) return ZC_Z_mul(X, u);
     884     9042176 :   if (is_pm1(v))
     885             :   {
     886     1622248 :     if (is_pm1(u))
     887             :     {
     888     1008729 :       if (su != sv) A = ZC_sub(X, Y);
     889      336981 :       else          A = ZC_add(X, Y);
     890     1008729 :       if (su < 0) ZV_togglesign(A); /* in place but was created above */
     891             :     }
     892             :     else
     893             :     {
     894      613519 :       if (sv > 0) A = ZC_lincomb1 (u, Y, X);
     895      243647 :       else        A = ZC_lincomb_1(u, Y, X);
     896             :     }
     897             :   }
     898     7419928 :   else if (is_pm1(u))
     899             :   {
     900     1871223 :     if (su > 0) A = ZC_lincomb1 (v, X, Y);
     901      138596 :     else        A = ZC_lincomb_1(v, X, Y);
     902             :   }
     903             :   else
     904             :   { /* not cgetg_copy: x may be a t_VEC */
     905     5548705 :     long i, lx = lg(X);
     906     5548705 :     A = cgetg(lx,t_COL);
     907     5548705 :     for (i=1; i<lx; i++) gel(A,i) = lincombii(u,v,gel(X,i),gel(Y,i));
     908             :   }
     909     9042176 :   return A;
     910             : }
     911             : 
     912             : /********************************************************************/
     913             : /**                                                                **/
     914             : /**                           CONVERSIONS                          **/
     915             : /**                                                                **/
     916             : /********************************************************************/
     917             : GEN
     918       47383 : ZV_to_nv(GEN x)
     919       47383 : { pari_APPLY_ulong(itou(gel(x,i))) }
     920             : 
     921             : GEN
     922       37068 : zm_to_ZM(GEN x)
     923       37068 : { pari_APPLY_type(t_MAT, zc_to_ZC(gel(x,i))) }
     924             : 
     925             : GEN
     926          98 : zmV_to_ZMV(GEN x)
     927          98 : { pari_APPLY_type(t_VEC, zm_to_ZM(gel(x,i))) }
     928             : 
     929             : /* same as Flm_to_ZM but do not assume positivity */
     930             : GEN
     931         868 : ZM_to_zm(GEN x)
     932         868 : { pari_APPLY_same(ZV_to_zv(gel(x,i))) }
     933             : 
     934             : GEN
     935      223314 : zv_to_Flv(GEN x, ulong p)
     936      223314 : { pari_APPLY_ulong(umodsu(x[i], p)) }
     937             : 
     938             : GEN
     939       15792 : zm_to_Flm(GEN x, ulong p)
     940       15792 : { pari_APPLY_same(zv_to_Flv(gel(x,i),p)) }
     941             : 
     942             : GEN
     943          21 : ZMV_to_zmV(GEN x)
     944          21 : { pari_APPLY_type(t_VEC, ZM_to_zm(gel(x,i))) }
     945             : 
     946             : /********************************************************************/
     947             : /**                                                                **/
     948             : /**                         COPY, NEGATION                         **/
     949             : /**                                                                **/
     950             : /********************************************************************/
     951             : GEN
     952     3763560 : ZC_copy(GEN x)
     953             : {
     954     3763560 :   long i, lx = lg(x);
     955     3763560 :   GEN y = cgetg(lx, t_COL);
     956    31785728 :   for (i=1; i<lx; i++)
     957             :   {
     958    28022168 :     GEN c = gel(x,i);
     959    28022168 :     gel(y,i) = lgefint(c) == 2? gen_0: icopy(c);
     960             :   }
     961     3763560 :   return y;
     962             : }
     963             : 
     964             : GEN
     965      211834 : ZM_copy(GEN x)
     966      211834 : { pari_APPLY_same(ZC_copy(gel(x,i))) }
     967             : 
     968             : void
     969      126769 : ZV_neg_inplace(GEN M)
     970             : {
     971      126769 :   long l = lg(M);
     972      126769 :   while (--l > 0) gel(M,l) = negi(gel(M,l));
     973      126769 : }
     974             : GEN
     975     1868555 : ZC_neg(GEN x)
     976     1868555 : { pari_APPLY_type(t_COL, negi(gel(x,i))) }
     977             : GEN
     978       23490 : zv_neg(GEN x)
     979       23490 : { pari_APPLY_long(-x[i]) }
     980             : GEN
     981         268 : zv_neg_inplace(GEN M)
     982             : {
     983         268 :   long l = lg(M);
     984         268 :   while (--l > 0) M[l] = -M[l];
     985         268 :   return M;
     986             : }
     987             : GEN
     988      273469 : ZM_neg(GEN x)
     989      273469 : { pari_APPLY_same(ZC_neg(gel(x,i))) }
     990             : 
     991             : void
     992     1398559 : ZV_togglesign(GEN M)
     993             : {
     994     1398559 :   long l = lg(M);
     995     1398559 :   while (--l > 0) togglesign_safe(&gel(M,l));
     996     1398559 : }
     997             : void
     998           0 : ZM_togglesign(GEN M)
     999             : {
    1000           0 :   long l = lg(M);
    1001           0 :   while (--l > 0) ZV_togglesign(gel(M,l));
    1002           0 : }
    1003             : 
    1004             : /********************************************************************/
    1005             : /**                                                                **/
    1006             : /**                        "DIVISION" mod HNF                      **/
    1007             : /**                                                                **/
    1008             : /********************************************************************/
    1009             : /* Reduce ZC x modulo ZM y in HNF, may return x itself (not a copy) */
    1010             : GEN
    1011      851990 : ZC_hnfremdiv(GEN x, GEN y, GEN *Q)
    1012             : {
    1013      851990 :   long i, l = lg(x);
    1014             :   GEN q;
    1015             : 
    1016      851990 :   if (Q) *Q = cgetg(l,t_COL);
    1017      851990 :   if (l == 1) return cgetg(1,t_COL);
    1018     5054202 :   for (i = l-1; i>0; i--)
    1019             :   {
    1020     4202212 :     q = diviiround(gel(x,i), gcoeff(y,i,i));
    1021     4202212 :     if (signe(q)) {
    1022     2272065 :       togglesign(q);
    1023     2272065 :       x = ZC_lincomb(gen_1, q, x, gel(y,i));
    1024             :     }
    1025     4202212 :     if (Q) gel(*Q, i) = q;
    1026             :   }
    1027      851990 :   return x;
    1028             : }
    1029             : 
    1030             : /* x = y Q + R, may return some columns of x (not copies) */
    1031             : GEN
    1032       51574 : ZM_hnfdivrem(GEN x, GEN y, GEN *Q)
    1033             : {
    1034       51574 :   long lx = lg(x), i;
    1035       51574 :   GEN R = cgetg(lx, t_MAT);
    1036       51574 :   if (Q)
    1037             :   {
    1038       15964 :     GEN q = cgetg(lx, t_MAT); *Q = q;
    1039       15964 :     for (i=1; i<lx; i++) gel(R,i) = ZC_hnfremdiv(gel(x,i),y,(GEN*)(q+i));
    1040             :   }
    1041             :   else
    1042      121143 :     for (i=1; i<lx; i++)
    1043             :     {
    1044       85533 :       pari_sp av = avma;
    1045       85533 :       GEN z = ZC_hnfrem(gel(x,i),y);
    1046       85533 :       gel(R,i) = (avma == av)? ZC_copy(z): gerepileupto(av, z);
    1047             :     }
    1048       51574 :   return R;
    1049             : }
    1050             : 
    1051             : 
    1052             : /********************************************************************/
    1053             : /**                                                                **/
    1054             : /**                               TESTS                            **/
    1055             : /**                                                                **/
    1056             : /********************************************************************/
    1057             : int
    1058    17351224 : zv_equal0(GEN V)
    1059             : {
    1060    17351224 :   long l = lg(V);
    1061    45112312 :   while (--l > 0)
    1062    22458154 :     if (V[l]) return 0;
    1063     5302934 :   return 1;
    1064             : }
    1065             : 
    1066             : int
    1067      844280 : ZV_equal0(GEN V)
    1068             : {
    1069      844280 :   long l = lg(V);
    1070     3401709 :   while (--l > 0)
    1071     2370029 :     if (signe(gel(V,l))) return 0;
    1072      187400 :   return 1;
    1073             : }
    1074             : 
    1075             : static int
    1076     9191853 : ZV_equal_lg(GEN V, GEN W, long l)
    1077             : {
    1078    26309388 :   while (--l > 0)
    1079    15108269 :     if (!equalii(gel(V,l), gel(W,l))) return 0;
    1080     2009266 :   return 1;
    1081             : }
    1082             : int
    1083     7962080 : ZV_equal(GEN V, GEN W)
    1084             : {
    1085     7962080 :   long l = lg(V);
    1086     7962080 :   if (lg(W) != l) return 0;
    1087     7962080 :   return ZV_equal_lg(V, W, l);
    1088             : }
    1089             : int
    1090      399212 : ZM_equal(GEN A, GEN B)
    1091             : {
    1092      399212 :   long i, m, l = lg(A);
    1093      399212 :   if (lg(B) != l) return 0;
    1094      398813 :   if (l == 1) return 1;
    1095      398813 :   m = lgcols(A);
    1096      398813 :   if (lgcols(B) != m) return 0;
    1097     1556598 :   for (i = 1; i < l; i++)
    1098     1229773 :     if (!ZV_equal_lg(gel(A,i), gel(B,i), m)) return 0;
    1099      326825 :   return 1;
    1100             : }
    1101             : int
    1102          56 : ZM_equal0(GEN A)
    1103             : {
    1104          56 :   long i, j, m, l = lg(A);
    1105          56 :   if (l == 1) return 1;
    1106          56 :   m = lgcols(A);
    1107          63 :   for (j = 1; j < l; j++)
    1108         217 :     for (i = 1; i < m; i++)
    1109         210 :       if (signe(gcoeff(A,i,j))) return 0;
    1110           0 :   return 1;
    1111             : }
    1112             : int
    1113    20974855 : zv_equal(GEN V, GEN W)
    1114             : {
    1115    20974855 :   long l = lg(V);
    1116    20974855 :   if (lg(W) != l) return 0;
    1117   199274316 :   while (--l > 0)
    1118   158090919 :     if (V[l] != W[l]) return 0;
    1119    20570407 :   return 1;
    1120             : }
    1121             : 
    1122             : int
    1123         385 : zvV_equal(GEN V, GEN W)
    1124             : {
    1125         385 :   long l = lg(V);
    1126         385 :   if (lg(W) != l) return 0;
    1127       58982 :   while (--l > 0)
    1128       58373 :     if (!zv_equal(gel(V,l),gel(W,l))) return 0;
    1129         231 :   return 1;
    1130             : }
    1131             : 
    1132             : int
    1133      114730 : ZM_ishnf(GEN x)
    1134             : {
    1135      114730 :   long i,j, lx = lg(x);
    1136      364207 :   for (i=1; i<lx; i++)
    1137             :   {
    1138      276231 :     GEN xii = gcoeff(x,i,i);
    1139      276231 :     if (signe(xii) <= 0) return 0;
    1140      675384 :     for (j=1; j<i; j++)
    1141      412558 :       if (signe(gcoeff(x,i,j))) return 0;
    1142      717440 :     for (j=i+1; j<lx; j++)
    1143             :     {
    1144      467963 :       GEN xij = gcoeff(x,i,j);
    1145      467963 :       if (signe(xij)<0 || cmpii(xij,xii)>=0) return 0;
    1146             :     }
    1147             :   }
    1148       87976 :   return 1;
    1149             : }
    1150             : int
    1151      108109 : ZM_isidentity(GEN x)
    1152             : {
    1153      108109 :   long i,j, lx = lg(x);
    1154             : 
    1155      108109 :   if (lx == 1) return 1;
    1156      108109 :   if (lx != lgcols(x)) return 0;
    1157      798321 :   for (j=1; j<lx; j++)
    1158             :   {
    1159      690226 :     GEN c = gel(x,j);
    1160     3897702 :     for (i=1; i<j; )
    1161     2517250 :       if (signe(gel(c,i++))) return 0;
    1162             :     /* i = j */
    1163      690226 :       if (!equali1(gel(c,i++))) return 0;
    1164     3897688 :     for (   ; i<lx; )
    1165     2517257 :       if (signe(gel(c,i++))) return 0;
    1166             :   }
    1167      108095 :   return 1;
    1168             : }
    1169             : int
    1170       60254 : ZM_isdiagonal(GEN x)
    1171             : {
    1172       60254 :   long i,j, lx = lg(x);
    1173       60254 :   if (lx == 1) return 1;
    1174       60254 :   if (lx != lgcols(x)) return 0;
    1175             : 
    1176      186729 :   for (j=1; j<lx; j++)
    1177             :   {
    1178      162085 :     GEN c = gel(x,j);
    1179      356768 :     for (i=1; i<j; i++)
    1180      230286 :       if (signe(gel(c,i))) return 0;
    1181      520277 :     for (i++; i<lx; i++)
    1182      393802 :       if (signe(gel(c,i))) return 0;
    1183             :   }
    1184       24644 :   return 1;
    1185             : }
    1186             : int
    1187      100384 : ZM_isscalar(GEN x, GEN s)
    1188             : {
    1189      100384 :   long i, j, lx = lg(x);
    1190             : 
    1191      100384 :   if (lx == 1) return 1;
    1192      100384 :   if (!s) s = gcoeff(x,1,1);
    1193      100384 :   if (equali1(s)) return ZM_isidentity(x);
    1194       17174 :   if (lx != lgcols(x)) return 0;
    1195      144903 :   for (j=1; j<lx; j++)
    1196             :   {
    1197      129447 :     GEN c = gel(x,j);
    1198     1231443 :     for (i=1; i<j; )
    1199      973456 :       if (signe(gel(c,i++))) return 0;
    1200             :     /* i = j */
    1201      128540 :       if (!equalii(gel(c,i++), s)) return 0;
    1202     1236827 :     for (   ; i<lx; )
    1203      981311 :       if (signe(gel(c,i++))) return 0;
    1204             :   }
    1205       15456 :   return 1;
    1206             : }
    1207             : 
    1208             : long
    1209       59269 : ZC_is_ei(GEN x)
    1210             : {
    1211       59269 :   long i, j = 0, l = lg(x);
    1212      555093 :   for (i = 1; i < l; i++)
    1213             :   {
    1214      495824 :     GEN c = gel(x,i);
    1215      495824 :     long s = signe(c);
    1216      495824 :     if (!s) continue;
    1217       59255 :     if (s < 0 || !is_pm1(c) || j) return 0;
    1218       59255 :     j = i;
    1219             :   }
    1220       59269 :   return j;
    1221             : }
    1222             : 
    1223             : /********************************************************************/
    1224             : /**                                                                **/
    1225             : /**                       MISCELLANEOUS                            **/
    1226             : /**                                                                **/
    1227             : /********************************************************************/
    1228             : /* assume lg(x) = lg(y), x,y in Z^n */
    1229             : int
    1230      660112 : ZV_cmp(GEN x, GEN y)
    1231             : {
    1232      660112 :   long fl,i, lx = lg(x);
    1233     1319243 :   for (i=1; i<lx; i++)
    1234     1061396 :     if (( fl = cmpii(gel(x,i), gel(y,i)) )) return fl;
    1235      257847 :   return 0;
    1236             : }
    1237             : /* assume lg(x) = lg(y), x,y in Z^n */
    1238             : int
    1239        3557 : ZV_abscmp(GEN x, GEN y)
    1240             : {
    1241        3557 :   long fl,i, lx = lg(x);
    1242       19608 :   for (i=1; i<lx; i++)
    1243       19528 :     if (( fl = abscmpii(gel(x,i), gel(y,i)) )) return fl;
    1244          80 :   return 0;
    1245             : }
    1246             : 
    1247             : long
    1248     1015342 : zv_content(GEN x)
    1249             : {
    1250     1015342 :   long i, s, l = lg(x);
    1251     1015342 :   if (l == 1) return 0;
    1252     1015335 :   s = labs(x[1]);
    1253     1015335 :   for (i=2; i<l && s!=1; i++) s = cgcd(x[i],s);
    1254     1015335 :   return s;
    1255             : }
    1256             : GEN
    1257      157129 : ZV_content(GEN x)
    1258             : {
    1259      157129 :   long i, l = lg(x);
    1260      157129 :   pari_sp av = avma;
    1261             :   GEN c;
    1262      157129 :   if (l == 1) return gen_0;
    1263      157129 :   if (l == 2) return absi(gel(x,1));
    1264      119980 :   c = gel(x,1);
    1265      296597 :   for (i = 2; i < l; i++)
    1266             :   {
    1267      206976 :     c = gcdii(c, gel(x,i));
    1268      206976 :     if (is_pm1(c)) { avma = av; return gen_1; }
    1269             :   }
    1270       89621 :   return gerepileuptoint(av, c);
    1271             : }
    1272             : 
    1273             : GEN
    1274      978571 : ZM_det_triangular(GEN mat)
    1275             : {
    1276             :   pari_sp av;
    1277      978571 :   long i,l = lg(mat);
    1278             :   GEN s;
    1279             : 
    1280      978571 :   if (l<3) return l<2? gen_1: icopy(gcoeff(mat,1,1));
    1281      898490 :   av = avma; s = gcoeff(mat,1,1);
    1282      898490 :   for (i=2; i<l; i++) s = mulii(s,gcoeff(mat,i,i));
    1283      898490 :   return gerepileuptoint(av,s);
    1284             : }
    1285             : 
    1286             : /* assumes no overflow */
    1287             : long
    1288      363314 : zv_prod(GEN v)
    1289             : {
    1290      363314 :   long n, i, l = lg(v);
    1291      363314 :   if (l == 1) return 1;
    1292      270305 :   n = v[1]; for (i = 2; i < l; i++) n *= v[i];
    1293      270305 :   return n;
    1294             : }
    1295             : 
    1296             : static GEN
    1297   272771977 : _mulii(void *E, GEN a, GEN b)
    1298   272771977 : { (void) E; return mulii(a, b); }
    1299             : 
    1300             : /* product of ulongs */
    1301             : GEN
    1302      364852 : zv_prod_Z(GEN v)
    1303             : {
    1304      364852 :   pari_sp av = avma;
    1305      364852 :   long k, m, n = lg(v)-1;
    1306             :   GEN V;
    1307      364852 :   switch(n) {
    1308       19509 :     case 0: return gen_1;
    1309       56553 :     case 1: return utoi(v[1]);
    1310      125426 :     case 2: return muluu(v[1], v[2]);
    1311             :   }
    1312      163364 :   m = n >> 1;
    1313      163364 :   V = cgetg(m + (odd(n)? 2: 1), t_VEC);
    1314      163364 :   for (k = 1; k <= m; k++) gel(V,k) = muluu(v[k<<1], v[(k<<1)-1]);
    1315      163364 :   if (odd(n)) gel(V,k) = utoipos(v[n]);
    1316      163364 :   return gerepileuptoint(av, gen_product(V, NULL, &_mulii));
    1317             : }
    1318             : GEN
    1319    14694393 : vecsmall_prod(GEN v)
    1320             : {
    1321    14694393 :   pari_sp av = avma;
    1322    14694393 :   long k, m, n = lg(v)-1;
    1323             :   GEN V;
    1324    14694393 :   switch (n) {
    1325           0 :     case 0: return gen_1;
    1326           0 :     case 1: return stoi(v[1]);
    1327          21 :     case 2: return mulss(v[1], v[2]);
    1328             :   }
    1329    14694372 :   m = n >> 1;
    1330    14694372 :   V = cgetg(m + (odd(n)? 2: 1), t_VEC);
    1331    14694372 :   for (k = 1; k <= m; k++) gel(V,k) = mulss(v[k<<1], v[(k<<1)-1]);
    1332    14694372 :   if (odd(n)) gel(V,k) = stoi(v[n]);
    1333    14694372 :   return gerepileuptoint(av, gen_product(V, NULL, &_mulii));
    1334             : }
    1335             : 
    1336             : GEN
    1337     7471457 : ZV_prod(GEN v)
    1338             : {
    1339     7471457 :   pari_sp av = avma;
    1340     7471457 :   long i, l = lg(v);
    1341             :   GEN n;
    1342     7471457 :   if (l == 1) return gen_1;
    1343     7462259 :   if (l > 7) return gerepileuptoint(av, gen_product(v, NULL, _mulii));
    1344      114097 :   n = gel(v,1);
    1345      114097 :   if (l == 2) return icopy(n);
    1346       62101 :   for (i = 2; i < l; i++) n = mulii(n, gel(v,i));
    1347       62101 :   return gerepileuptoint(av, n);
    1348             : }
    1349             : /* assumes no overflow */
    1350             : long
    1351        1596 : zv_sum(GEN v)
    1352             : {
    1353        1596 :   long n, i, l = lg(v);
    1354        1596 :   if (l == 1) return 0;
    1355        1575 :   n = v[1]; for (i = 2; i < l; i++) n += v[i];
    1356        1575 :   return n;
    1357             : }
    1358             : GEN
    1359          56 : ZV_sum(GEN v)
    1360             : {
    1361          56 :   pari_sp av = avma;
    1362          56 :   long i, l = lg(v);
    1363             :   GEN n;
    1364          56 :   if (l == 1) return gen_0;
    1365          56 :   n = gel(v,1);
    1366          56 :   if (l == 2) return icopy(n);
    1367          56 :   for (i = 2; i < l; i++) n = addii(n, gel(v,i));
    1368          56 :   return gerepileuptoint(av, n);
    1369             : }
    1370             : 
    1371             : /********************************************************************/
    1372             : /**                                                                **/
    1373             : /**         GRAM SCHMIDT REDUCTION (integer matrices)              **/
    1374             : /**                                                                **/
    1375             : /********************************************************************/
    1376             : 
    1377             : /* L[k,] += q * L[l,], l < k. Inefficient if q = 0 */
    1378             : static void
    1379       52048 : Zupdate_row(long k, long l, GEN q, GEN L, GEN B)
    1380             : {
    1381       52048 :   long i, qq = itos_or_0(q);
    1382       52048 :   if (!qq)
    1383             :   {
    1384        2662 :     for(i=1;i<l;i++)  gcoeff(L,k,i) = addii(gcoeff(L,k,i),mulii(q,gcoeff(L,l,i)));
    1385        2662 :     gcoeff(L,k,l) = addii(gcoeff(L,k,l), mulii(q,B));
    1386       54710 :     return;
    1387             :   }
    1388       49386 :   if (qq == 1) {
    1389       13664 :     for (i=1;i<l; i++) gcoeff(L,k,i) = addii(gcoeff(L,k,i),gcoeff(L,l,i));
    1390       13664 :     gcoeff(L,k,l) = addii(gcoeff(L,k,l), B);
    1391       35722 :   } else if (qq == -1) {
    1392       10897 :     for (i=1;i<l; i++) gcoeff(L,k,i) = subii(gcoeff(L,k,i),gcoeff(L,l,i));
    1393       10897 :     gcoeff(L,k,l) = addii(gcoeff(L,k,l), negi(B));
    1394             :   } else {
    1395       24825 :     for(i=1;i<l;i++) gcoeff(L,k,i) = addii(gcoeff(L,k,i),mulsi(qq,gcoeff(L,l,i)));
    1396       24825 :     gcoeff(L,k,l) = addii(gcoeff(L,k,l), mulsi(qq,B));
    1397             :   }
    1398             : }
    1399             : 
    1400             : /* update L[k,] */
    1401             : static void
    1402      245816 : ZRED(long k, long l, GEN x, GEN L, GEN B)
    1403             : {
    1404      245816 :   GEN q = truedivii(addii(B,shifti(gcoeff(L,k,l),1)), shifti(B,1));
    1405      491632 :   if (!signe(q)) return;
    1406       52048 :   q = negi(q);
    1407       52048 :   Zupdate_row(k,l,q,L,B);
    1408       52048 :   gel(x,k) = ZC_lincomb(gen_1, q, gel(x,k), gel(x,l));
    1409             : }
    1410             : 
    1411             : /* Gram-Schmidt reduction, x a ZM */
    1412             : static void
    1413      311919 : ZincrementalGS(GEN x, GEN L, GEN B, long k)
    1414             : {
    1415             :   long i, j;
    1416     1109028 :   for (j=1; j<=k; j++)
    1417             :   {
    1418      797109 :     pari_sp av = avma;
    1419      797109 :     GEN u = ZV_dotproduct(gel(x,k), gel(x,j));
    1420     1776781 :     for (i=1; i<j; i++)
    1421             :     {
    1422      979672 :       u = subii(mulii(gel(B,i+1), u), mulii(gcoeff(L,k,i), gcoeff(L,j,i)));
    1423      979672 :       u = diviiexact(u, gel(B,i));
    1424             :     }
    1425      797109 :     gcoeff(L,k,j) = gerepileuptoint(av, u);
    1426             :   }
    1427      311919 :   gel(B,k+1) = gcoeff(L,k,k); gcoeff(L,k,k) = gen_1;
    1428      311919 : }
    1429             : 
    1430             : /* Variant reducemodinvertible(ZC v, ZM y), when y singular.
    1431             :  * Very inefficient if y is not LLL-reduced of maximal rank */
    1432             : static GEN
    1433       62512 : ZC_reducemodmatrix_i(GEN v, GEN y)
    1434             : {
    1435       62512 :   GEN B, L, x = shallowconcat(y, v);
    1436       62512 :   long k, lx = lg(x), nx = lx-1;
    1437             : 
    1438       62512 :   B = scalarcol_shallow(gen_1, lx);
    1439       62512 :   L = zeromatcopy(nx, nx);
    1440       62512 :   for (k=1; k <= nx; k++) ZincrementalGS(x, L, B, k);
    1441       62512 :   for (k = nx-1; k >= 1; k--) ZRED(nx,k, x,L,gel(B,k+1));
    1442       62512 :   return gel(x,nx);
    1443             : }
    1444             : GEN
    1445       62512 : ZC_reducemodmatrix(GEN v, GEN y) {
    1446       62512 :   pari_sp av = avma;
    1447       62512 :   return gerepilecopy(av, ZC_reducemodmatrix_i(v,y));
    1448             : }
    1449             : 
    1450             : /* Variant reducemodinvertible(ZM v, ZM y), when y singular.
    1451             :  * Very inefficient if y is not LLL-reduced of maximal rank */
    1452             : static GEN
    1453       21280 : ZM_reducemodmatrix_i(GEN v, GEN y)
    1454             : {
    1455             :   GEN B, L, V;
    1456       21280 :   long j, k, lv = lg(v), nx = lg(y), lx = nx+1;
    1457             : 
    1458       21280 :   V = cgetg(lv, t_MAT);
    1459       21280 :   B = scalarcol_shallow(gen_1, lx);
    1460       21280 :   L = zeromatcopy(nx, nx);
    1461       21280 :   for (k=1; k < nx; k++) ZincrementalGS(y, L, B, k);
    1462       61551 :   for (j = 1; j < lg(v); j++)
    1463             :   {
    1464       40271 :     GEN x = shallowconcat(y, gel(v,j));
    1465       40271 :     ZincrementalGS(x, L, B, nx); /* overwrite last */
    1466       40271 :     for (k = nx-1; k >= 1; k--) ZRED(nx,k, x,L,gel(B,k+1));
    1467       40271 :     gel(V,j) = gel(x,nx);
    1468             :   }
    1469       21280 :   return V;
    1470             : }
    1471             : GEN
    1472       21280 : ZM_reducemodmatrix(GEN v, GEN y) {
    1473       21280 :   pari_sp av = avma;
    1474       21280 :   return gerepilecopy(av, ZM_reducemodmatrix_i(v,y));
    1475             : }
    1476             : 
    1477             : GEN
    1478        9907 : ZC_reducemodlll(GEN x,GEN y)
    1479             : {
    1480        9907 :   pari_sp av = avma;
    1481        9907 :   GEN z = ZC_reducemodmatrix(x, ZM_lll(y, 0.75, LLL_INPLACE));
    1482        9907 :   return gerepilecopy(av, z);
    1483             : }
    1484             : GEN
    1485           0 : ZM_reducemodlll(GEN x,GEN y)
    1486             : {
    1487           0 :   pari_sp av = avma;
    1488           0 :   GEN z = ZM_reducemodmatrix(x, ZM_lll(y, 0.75, LLL_INPLACE));
    1489           0 :   return gerepilecopy(av, z);
    1490             : }

Generated by: LCOV version 1.11