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 - alglin3.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 25819-e703fe1174) Lines: 529 571 92.6 %
Date: 2020-09-18 06:10:04 Functions: 47 53 88.7 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2012  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             : /********************************************************************/
      15             : /**                                                                **/
      16             : /**                         LINEAR ALGEBRA                         **/
      17             : /**                          (third part)                          **/
      18             : /**                                                                **/
      19             : /********************************************************************/
      20             : #include "pari.h"
      21             : #include "paripriv.h"
      22             : 
      23             : /*******************************************************************/
      24             : /*                                                                 */
      25             : /*                               SUM                               */
      26             : /*                                                                 */
      27             : /*******************************************************************/
      28             : 
      29             : GEN
      30       94160 : vecsum(GEN v)
      31             : {
      32       94160 :   pari_sp av = avma;
      33             :   long i, l;
      34             :   GEN p;
      35       94160 :   if (!is_vec_t(typ(v)))
      36           7 :     pari_err_TYPE("vecsum", v);
      37       94153 :   l = lg(v);
      38       94153 :   if (l == 1) return gen_0;
      39       94146 :   p = gel(v,1);
      40       94146 :   if (l == 2) return gcopy(p);
      41      109229 :   for (i=2; i<l; i++)
      42             :   {
      43       69686 :     p = gadd(p, gel(v,i));
      44       69686 :     if (gc_needed(av, 2))
      45             :     {
      46           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"sum");
      47           0 :       p = gerepileupto(av, p);
      48             :     }
      49             :   }
      50       39543 :   return gerepileupto(av, p);
      51             : }
      52             : 
      53             : /*******************************************************************/
      54             : /*                                                                 */
      55             : /*                         TRANSPOSE                               */
      56             : /*                                                                 */
      57             : /*******************************************************************/
      58             : /* A[x0,]~ */
      59             : static GEN
      60    11021422 : row_transpose(GEN A, long x0)
      61             : {
      62    11021422 :   long i, lB = lg(A);
      63    11021422 :   GEN B  = cgetg(lB, t_COL);
      64   115471405 :   for (i=1; i<lB; i++) gel(B, i) = gcoeff(A, x0, i);
      65    11021422 :   return B;
      66             : }
      67             : static GEN
      68       18473 : row_transposecopy(GEN A, long x0)
      69             : {
      70       18473 :   long i, lB = lg(A);
      71       18473 :   GEN B  = cgetg(lB, t_COL);
      72      158830 :   for (i=1; i<lB; i++) gel(B, i) = gcopy(gcoeff(A, x0, i));
      73       18473 :   return B;
      74             : }
      75             : 
      76             : /* No copy*/
      77             : GEN
      78     3036235 : shallowtrans(GEN x)
      79             : {
      80             :   long i, dx, lx;
      81             :   GEN y;
      82     3036235 :   switch(typ(x))
      83             :   {
      84        1624 :     case t_VEC: y = leafcopy(x); settyp(y,t_COL); break;
      85       19351 :     case t_COL: y = leafcopy(x); settyp(y,t_VEC); break;
      86     3015260 :     case t_MAT:
      87     3015260 :       lx = lg(x); if (lx==1) return cgetg(1,t_MAT);
      88     3015260 :       dx = lgcols(x); y = cgetg(dx,t_MAT);
      89    14036682 :       for (i = 1; i < dx; i++) gel(y,i) = row_transpose(x,i);
      90     3015260 :       break;
      91           0 :     default: pari_err_TYPE("shallowtrans",x);
      92             :       return NULL;/*LCOV_EXCL_LINE*/
      93             :   }
      94     3036235 :   return y;
      95             : }
      96             : 
      97             : GEN
      98       41657 : gtrans(GEN x)
      99             : {
     100             :   long i, dx, lx;
     101             :   GEN y;
     102       41657 :   switch(typ(x))
     103             :   {
     104       34937 :     case t_VEC: y = gcopy(x); settyp(y,t_COL); break;
     105        4766 :     case t_COL: y = gcopy(x); settyp(y,t_VEC); break;
     106        1947 :     case t_MAT:
     107        1947 :       lx = lg(x); if (lx==1) return cgetg(1,t_MAT);
     108        1940 :       dx = lgcols(x); y = cgetg(dx,t_MAT);
     109       20413 :       for (i = 1; i < dx; i++) gel(y,i) = row_transposecopy(x,i);
     110        1940 :       break;
     111           7 :     default: pari_err_TYPE("gtrans",x);
     112             :       return NULL;/*LCOV_EXCL_LINE*/
     113             :   }
     114       41643 :   return y;
     115             : }
     116             : 
     117             : /*******************************************************************/
     118             : /*                                                                 */
     119             : /*                           EXTRACTION                            */
     120             : /*                                                                 */
     121             : /*******************************************************************/
     122             : 
     123             : static long
     124         182 : str_to_long(char *s, char **pt)
     125             : {
     126         182 :   long a = atol(s);
     127         182 :   while (isspace((int)*s)) s++;
     128         182 :   if (*s == '-' || *s == '+') s++;
     129         385 :   while (isdigit((int)*s) || isspace((int)*s)) s++;
     130         182 :   *pt = s; return a;
     131             : }
     132             : 
     133             : static int
     134         112 : get_range(char *s, long *a, long *b, long *cmpl, long lx)
     135             : {
     136         112 :   long max = lx - 1;
     137             : 
     138         112 :   *a = 1; *b = max;
     139         112 :   if (*s == '^') { *cmpl = 1; s++; } else *cmpl = 0;
     140         112 :   if (!*s) return 0;
     141         112 :   if (*s != '.')
     142             :   {
     143         105 :     *a = str_to_long(s, &s);
     144         105 :     if (*a < 0) *a += lx;
     145         105 :     if (*a<1 || *a>max) return 0;
     146             :   }
     147         112 :   if (*s == '.')
     148             :   {
     149         105 :     s++; if (*s != '.') return 0;
     150         105 :     do s++; while (isspace((int)*s));
     151         105 :     if (*s)
     152             :     {
     153          77 :       *b = str_to_long(s, &s);
     154          77 :       if (*b < 0) *b += lx;
     155          77 :       if (*b<1 || *b>max || *s) return 0;
     156             :     }
     157          98 :     return 1;
     158             :   }
     159           7 :   if (*s) return 0;
     160           7 :   *b = *a; return 1;
     161             : }
     162             : 
     163             : static int
     164          35 : extract_selector_ok(long lx, GEN L)
     165             : {
     166             :   long i, l;
     167          35 :   switch (typ(L))
     168             :   {
     169           7 :     case t_INT: {
     170             :       long maxj;
     171           7 :       if (!signe(L)) return 1;
     172           7 :       l = lgefint(L)-1;
     173           7 :       maxj = BITS_IN_LONG - bfffo(*int_MSW(L));
     174           7 :       return ((l-2) * BITS_IN_LONG + maxj < lx);
     175             :     }
     176           7 :     case t_STR: {
     177             :       long first, last, cmpl;
     178           7 :       return get_range(GSTR(L), &first, &last, &cmpl, lx);
     179             :     }
     180          14 :     case t_VEC: case t_COL:
     181          14 :       l = lg(L);
     182          28 :       for (i=1; i<l; i++)
     183             :       {
     184          21 :         long j = itos(gel(L,i));
     185          21 :         if (j>=lx || j<=0) return 0;
     186             :       }
     187           7 :       return 1;
     188           7 :     case t_VECSMALL:
     189           7 :       l = lg(L);
     190          21 :       for (i=1; i<l; i++)
     191             :       {
     192          14 :         long j = L[i];
     193          14 :         if (j>=lx || j<=0) return 0;
     194             :       }
     195           7 :       return 1;
     196             :   }
     197           0 :   return 0;
     198             : }
     199             : 
     200             : GEN
     201        9849 : shallowmatextract(GEN x, GEN l1, GEN l2)
     202             : {
     203        9849 :   long i, j, n1 = lg(l1), n2 = lg(l2);
     204        9849 :   GEN M = cgetg(n2, t_MAT);
     205       62972 :   for(i=1; i < n2; i++)
     206             :   {
     207       53123 :     long ii = l2[i];
     208       53123 :     GEN C = cgetg(n1, t_COL);
     209      725354 :     for (j=1; j < n1; j++)
     210             :     {
     211      672231 :       long jj = l1[j];
     212      672231 :       gel(C, j) = gmael(x, ii, jj);
     213             :     }
     214       53123 :     gel(M, i) = C;
     215             :   }
     216        9849 :   return M;
     217             : }
     218             : 
     219             : GEN
     220       36198 : shallowextract(GEN x, GEN L)
     221             : {
     222       36198 :   long i,j, tl = typ(L), tx = typ(x), lx = lg(x);
     223             :   GEN y;
     224             : 
     225       36198 :   switch(tx)
     226             :   {
     227       36191 :     case t_VEC:
     228             :     case t_COL:
     229             :     case t_MAT:
     230       36191 :     case t_VECSMALL: break;
     231           7 :     default: pari_err_TYPE("extract",x);
     232             : 
     233             :   }
     234       36191 :   if (tl==t_INT)
     235             :   { /* extract components of x as per the bits of mask L */
     236             :     long k, l, ix, iy, maxj;
     237             :     GEN Ld;
     238        3276 :     if (!signe(L)) return cgetg(1,tx);
     239        3269 :     y = new_chunk(lx);
     240        3269 :     l = lgefint(L)-1; ix = iy = 1;
     241        3269 :     maxj = BITS_IN_LONG - bfffo(*int_MSW(L));
     242        3269 :     if ((l-2) * BITS_IN_LONG + maxj >= lx)
     243           7 :       pari_err_TYPE("vecextract [mask too large]", L);
     244        3637 :     for (k = 2, Ld = int_LSW(L); k < l; k++, Ld = int_nextW(Ld))
     245             :     {
     246         375 :       ulong B = *Ld;
     247       20439 :       for (j = 0; j < BITS_IN_LONG; j++, B >>= 1, ix++)
     248       20064 :         if (B & 1) y[iy++] = x[ix];
     249             :     }
     250             :     { /* k = l */
     251        3262 :       ulong B = *Ld;
     252       28404 :       for (j = 0; j < maxj; j++, B >>= 1, ix++)
     253       25142 :         if (B & 1) y[iy++] = x[ix];
     254             :     }
     255        3262 :     y[0] = evaltyp(tx) | evallg(iy);
     256        3262 :     return y;
     257             :   }
     258       32915 :   if (tl==t_STR)
     259             :   {
     260         105 :     char *s = GSTR(L);
     261             :     long first, last, cmpl, d;
     262         105 :     if (! get_range(s, &first, &last, &cmpl, lx))
     263           7 :       pari_err_TYPE("vecextract [incorrect range]", L);
     264          98 :     if (lx == 1) return cgetg(1,tx);
     265          98 :     d = last - first;
     266          98 :     if (cmpl)
     267             :     {
     268          21 :       if (d >= 0)
     269             :       {
     270          14 :         y = cgetg(lx - (1+d),tx);
     271         469 :         for (j=1; j<first; j++) gel(y,j) = gel(x,j);
     272         266 :         for (i=last+1; i<lx; i++,j++) gel(y,j) = gel(x,i);
     273             :       }
     274             :       else
     275             :       {
     276           7 :         y = cgetg(lx - (1-d),tx);
     277          14 :         for (j=1,i=lx-1; i>first; i--,j++) gel(y,j) = gel(x,i);
     278          14 :         for (i=last-1; i>0; i--,j++) gel(y,j) = gel(x,i);
     279             :       }
     280             :     }
     281             :     else
     282             :     {
     283          77 :       if (d >= 0)
     284             :       {
     285          35 :         y = cgetg(d+2,tx);
     286         112 :         for (i=first,j=1; i<=last; i++,j++) gel(y,j) = gel(x,i);
     287             :       }
     288             :       else
     289             :       {
     290          42 :         y = cgetg(2-d,tx);
     291         203 :         for (i=first,j=1; i>=last; i--,j++) gel(y,j) = gel(x,i);
     292             :       }
     293             :     }
     294          98 :     return y;
     295             :   }
     296             : 
     297       32810 :   if (is_vec_t(tl))
     298             :   {
     299          77 :     long ll=lg(L); y=cgetg(ll,tx);
     300         196 :     for (i=1; i<ll; i++)
     301             :     {
     302         133 :       j = itos(gel(L,i));
     303         133 :       if (j<=0) pari_err_COMPONENT("vecextract","<=",gen_0,stoi(j));
     304         126 :       if (j>=lx) pari_err_COMPONENT("vecextract",">=",stoi(lx),stoi(j));
     305         119 :       gel(y,i) = gel(x,j);
     306             :     }
     307          63 :     return y;
     308             :   }
     309       32733 :   if (tl == t_VECSMALL)
     310             :   {
     311       32726 :     long ll=lg(L); y=cgetg(ll,tx);
     312      145751 :     for (i=1; i<ll; i++)
     313             :     {
     314      113025 :       j = L[i];
     315      113025 :       if (j<=0) pari_err_COMPONENT("vecextract","<=",gen_0,stoi(j));
     316      113025 :       if (j>=lx) pari_err_COMPONENT("vecextract",">=",stoi(lx),stoi(j));
     317      113025 :       gel(y,i) = gel(x,j);
     318             :     }
     319       32726 :     return y;
     320             :   }
     321           7 :   pari_err_TYPE("vecextract [mask]", L);
     322             :   return NULL; /* LCOV_EXCL_LINE */
     323             : }
     324             : 
     325             : /* does the component selector l select 0 component ? */
     326             : static int
     327          70 : select_0(GEN l)
     328             : {
     329          70 :   switch(typ(l))
     330             :   {
     331          14 :     case t_INT:
     332          14 :       return (!signe(l));
     333          35 :     case t_VEC: case t_COL: case t_VECSMALL:
     334          35 :       return (lg(l) == 1);
     335             :   }
     336          21 :   return 0;
     337             : }
     338             : 
     339             : GEN
     340       28035 : extract0(GEN x, GEN l1, GEN l2)
     341             : {
     342       28035 :   pari_sp av = avma, av2;
     343             :   GEN y;
     344       28035 :   if (! l2)
     345             :   {
     346       27965 :     y = shallowextract(x, l1);
     347       27923 :     if (lg(y) == 1 || typ(y) == t_VECSMALL) return y;
     348       27916 :     av2 = avma;
     349       27916 :     y = gcopy(y);
     350             :   }
     351             :   else
     352             :   {
     353          70 :     if (typ(x) != t_MAT) pari_err_TYPE("extract",x);
     354          70 :     y = shallowextract(x,l2);
     355          70 :     if (select_0(l1)) { set_avma(av); return zeromat(0, lg(y)-1); }
     356          56 :     if (lg(y) == 1 && lg(x) > 1)
     357             :     {
     358          35 :       if (!extract_selector_ok(lgcols(x), l1))
     359           7 :         pari_err_TYPE("vecextract [incorrect mask]", l1);
     360          28 :       set_avma(av); return cgetg(1, t_MAT);
     361             :     }
     362          21 :     y = shallowextract(shallowtrans(y), l1);
     363          21 :     av2 = avma;
     364          21 :     y = gtrans(y);
     365             :   }
     366       27937 :   stackdummy(av, av2);
     367       27937 :   return y;
     368             : }
     369             : 
     370             : static long
     371        5167 : vecslice_parse_arg(long lA, long *y1, long *y2, long *skip)
     372             : {
     373        5167 :   *skip=0;
     374        5167 :   if (*y1==LONG_MAX)
     375             :   {
     376         210 :     if (*y2!=LONG_MAX)
     377             :     {
     378         119 :       if (*y2<0) *y2 += lA;
     379         119 :       if (*y2<0 || *y2==LONG_MAX || *y2>=lA)
     380           0 :         pari_err_DIM("_[..]");
     381         119 :       *skip=*y2;
     382             :     }
     383         210 :     *y1 = 1; *y2 = lA-1;
     384             :   }
     385        4957 :   else if (*y2==LONG_MAX) *y2 = *y1;
     386        5167 :   if (*y1<=0) *y1 += lA;
     387        5167 :   if (*y2<0) *y2 += lA;
     388        5167 :   if (*y1<=0 || *y1>*y2+1 || *y2>=lA) pari_err_DIM("_[..]");
     389        5153 :   return *y2 - *y1 + 2 - !!*skip;
     390             : }
     391             : 
     392             : static GEN
     393        5650 : vecslice_i(GEN A, long t, long lB, long y1, long skip)
     394             : {
     395        5650 :   GEN B = cgetg(lB, t);
     396             :   long i;
     397      144843 :   for (i=1; i<lB; i++, y1++)
     398             :   {
     399      139193 :     if (y1 == skip) { i--; continue; }
     400      139074 :     gel(B,i) = gcopy(gel(A,y1));
     401             :   }
     402        5650 :   return B;
     403             : }
     404             : 
     405             : static GEN
     406          14 : rowslice_i(GEN A, long lB, long x1, long y1, long skip)
     407             : {
     408          14 :   GEN B = cgetg(lB, t_VEC);
     409             :   long i;
     410          77 :   for (i=1; i<lB; i++, y1++)
     411             :   {
     412          63 :     if (y1 == skip) { i--; continue; }
     413          56 :     gel(B,i) = gcopy(gcoeff(A,x1,y1));
     414             :   }
     415          14 :   return B;
     416             : }
     417             : 
     418             : static GEN
     419           0 : rowsmallslice_i(GEN A, long lB, long x1, long y1, long skip)
     420             : {
     421           0 :   GEN B = cgetg(lB, t_VECSMALL);
     422             :   long i;
     423           0 :   for (i=1; i<lB; i++, y1++)
     424             :   {
     425           0 :     if (y1 == skip) { i--; continue; }
     426           0 :     B[i] = coeff(A,x1,y1);
     427             :   }
     428           0 :   return B;
     429             : }
     430             : 
     431             : static GEN
     432          28 : vecsmallslice_i(GEN A, long t, long lB, long y1, long skip)
     433             : {
     434          28 :   GEN B = cgetg(lB, t);
     435             :   long i;
     436         126 :   for (i=1; i<lB; i++, y1++)
     437             :   {
     438          98 :     if (y1 == skip) { i--; continue; }
     439          91 :     B[i] = A[y1];
     440             :   }
     441          28 :   return B;
     442             : }
     443             : GEN
     444        4824 : vecslice0(GEN A, long y1, long y2)
     445             : {
     446        4824 :   long skip, lB, t = typ(A);
     447        4824 :   switch(t)
     448             :   {
     449        4726 :     case t_VEC: case t_COL:
     450        4726 :       lB = vecslice_parse_arg(lg(A), &y1, &y2, &skip);
     451        4712 :       return vecslice_i(A, t,lB,y1,skip);
     452          28 :     case t_VECSMALL:
     453          28 :       lB = vecslice_parse_arg(lg(A), &y1, &y2, &skip);
     454          28 :       return vecsmallslice_i(A, t,lB,y1,skip);
     455          63 :     case t_LIST:
     456          63 :       if (list_typ(A) == t_LIST_RAW)
     457             :       {
     458          63 :         GEN y, z = list_data(A);
     459          63 :         long l = z? lg(z): 1;
     460          63 :         lB = vecslice_parse_arg(l, &y1, &y2, &skip);
     461          63 :         y = mklist(); if (!z) return y;
     462          63 :         list_data(y) = vecslice_i(z, t_VEC,lB,y1,skip);
     463          63 :         return y;
     464             :       }
     465             :     default:
     466           7 :       pari_err_TYPE("_[_.._]",A);
     467             :       return NULL;/*LCOV_EXCL_LINE*/
     468             :   }
     469             : }
     470             : 
     471             : GEN
     472         182 : matslice0(GEN A, long x1, long x2, long y1, long y2)
     473             : {
     474             :   GEN B;
     475         182 :   long i, lB, lA = lg(A), rA, t, skip, rskip, rlB;
     476         182 :   long is_col = y1!=LONG_MAX && y2==LONG_MAX;
     477         182 :   long is_row = x1!=LONG_MAX && x2==LONG_MAX;
     478             :   GEN (*slice)(GEN A, long t, long lB, long y1, long skip);
     479         182 :   if (typ(A)!=t_MAT) pari_err_TYPE("_[_.._,_.._]",A);
     480         182 :   lB = vecslice_parse_arg(lA, &y1, &y2, &skip);
     481         182 :   if (is_col) return vecslice0(gel(A, y1), x1, x2);
     482         168 :   rA = lg(A)==1 ? 1: lgcols(A);
     483         168 :   rlB = vecslice_parse_arg(rA, &x1, &x2, &rskip);
     484         168 :   t = lg(A)==1 ? t_COL: typ(gel(A,1));
     485         168 :   if (is_row) return t == t_COL ? rowslice_i(A, lB, x1, y1, skip):
     486           0 :                                   rowsmallslice_i(A, lB, x1, y1, skip);
     487         154 :   slice = t == t_COL? &vecslice_i: &vecsmallslice_i;
     488             : 
     489         154 :   B = cgetg(lB, t_MAT);
     490        1043 :   for (i=1; i<lB; i++, y1++)
     491             :   {
     492         889 :     if (y1 == skip) { i--; continue; }
     493         875 :     gel(B,i) = slice(gel(A,y1),t,rlB, x1, rskip);
     494             :   }
     495         154 :   return B;
     496             : }
     497             : 
     498             : GEN
     499       10090 : vecrange(GEN a, GEN b)
     500             : {
     501             :   GEN y;
     502             :   long i, l;
     503       10090 :   if (typ(a)!=t_INT) pari_err_TYPE("[_.._]",a);
     504       10083 :   if (typ(b)!=t_INT) pari_err_TYPE("[_.._]",b);
     505       10076 :   if (cmpii(a,b)>0) return cgetg(1,t_VEC);
     506       10069 :   l = itos(subii(b,a))+1;
     507       10069 :   a = setloop(a);
     508       10069 :   y = cgetg(l+1, t_VEC);
     509      371427 :   for (i=1; i<=l; a = incloop(a), i++) gel(y,i) = icopy(a);
     510       10069 :   return y;
     511             : }
     512             : 
     513             : GEN
     514           0 : vecrangess(long a, long b)
     515             : {
     516             :   GEN y;
     517             :   long i, l;
     518           0 :   if (a>b) return cgetg(1,t_VEC);
     519           0 :   l = b-a+1;
     520           0 :   y = cgetg(l+1, t_VEC);
     521           0 :   for (i=1; i<=l; a++, i++) gel(y,i) = stoi(a);
     522           0 :   return y;
     523             : }
     524             : 
     525             : GEN
     526          96 : genindexselect(void *E, long (*f)(void* E, GEN x), GEN A)
     527             : {
     528             :   long l, i, lv;
     529             :   GEN v, z;
     530             :   pari_sp av;
     531          96 :   clone_lock(A);
     532          96 :   switch(typ(A))
     533             :   {
     534           7 :     case t_LIST:
     535           7 :       z = list_data(A);
     536           7 :       l = z? lg(z): 1;
     537           7 :       break;
     538          82 :     case t_VEC: case t_COL: case t_MAT:
     539          82 :       l = lg(A);
     540          82 :       z = A;
     541          82 :       break;
     542           7 :     default:
     543           7 :       pari_err_TYPE("select",A);
     544             :       return NULL;/*LCOV_EXCL_LINE*/
     545             :   }
     546          89 :   v = cgetg(l, t_VECSMALL);
     547          89 :   av = avma;
     548       12801 :   for (i = lv = 1; i < l; i++) {
     549       12712 :     if (f(E, gel(z,i))) v[lv++] = i;
     550       12712 :     set_avma(av);
     551             :   }
     552          89 :   clone_unlock_deep(A); fixlg(v, lv); return v;
     553             : }
     554             : static GEN
     555          88 : extract_copy(GEN A, GEN v)
     556             : {
     557          88 :   long i, l = lg(v);
     558          88 :   GEN B = cgetg(l, typ(A));
     559         387 :   for (i = 1; i < l; i++) gel(B,i) = gcopy(gel(A,v[i]));
     560          88 :   return B;
     561             : }
     562             : /* as genselect, but treat A [ t_VEC,t_COL, or t_MAT] as a t_VEC */
     563             : GEN
     564           0 : vecselect(void *E, long (*f)(void* E, GEN x), GEN A)
     565             : {
     566             :   GEN v;
     567           0 :   clone_lock(A);
     568           0 :   v = genindexselect(E, f, A);
     569           0 :   A = extract_copy(A, v); settyp(A, t_VEC);
     570           0 :   clone_unlock_deep(A); return A;
     571             : }
     572             : GEN
     573          83 : genselect(void *E, long (*f)(void* E, GEN x), GEN A)
     574             : {
     575             :   GEN y, z, v;/* v left on stack for efficiency */
     576          83 :   clone_lock(A);
     577          83 :   switch(typ(A))
     578             :   {
     579          14 :     case t_LIST:
     580          14 :       z = list_data(A);
     581          14 :       if (!z) y = mklist();
     582             :       else
     583             :       {
     584             :         GEN B;
     585          14 :         y = cgetg(3, t_LIST);
     586          14 :         v = genindexselect(E, f, z);
     587          14 :         B = extract_copy(z, v);
     588          14 :         y[1] = lg(B)-1;
     589          14 :         list_data(y) = B;
     590             :       }
     591          14 :       break;
     592          62 :     case t_VEC: case t_COL: case t_MAT:
     593          62 :       v = genindexselect(E, f, A);
     594          62 :       y = extract_copy(A, v);
     595          62 :       break;
     596           7 :     default:
     597           7 :       pari_err_TYPE("select",A);
     598             :       return NULL;/*LCOV_EXCL_LINE*/
     599             :   }
     600          76 :   clone_unlock_deep(A); return y;
     601             : }
     602             : 
     603             : static void
     604        9655 : check_callgen1(GEN f, const char *s)
     605             : {
     606        9655 :   if (typ(f) != t_CLOSURE || closure_is_variadic(f)  || closure_arity(f) < 1)
     607           0 :     pari_err_TYPE(s, f);
     608        9655 : }
     609             : 
     610             : GEN
     611         103 : select0(GEN f, GEN x, long flag)
     612             : {
     613         103 :   check_callgen1(f, "select");
     614         103 :   switch(flag)
     615             :   {
     616          83 :     case 0: return genselect((void *) f, gp_callbool, x);
     617          20 :     case 1: return genindexselect((void *) f, gp_callbool, x);
     618           0 :     default: pari_err_FLAG("select");
     619             :              return NULL;/*LCOV_EXCL_LINE*/
     620             :   }
     621             : }
     622             : 
     623             : GEN
     624       23917 : parselect_worker(GEN d, GEN C)
     625             : {
     626       23917 :   return gequal0(closure_callgen1(C, d))? gen_0: gen_1;
     627             : }
     628             : 
     629             : GEN
     630          24 : parselect(GEN C, GEN D, long flag)
     631             : {
     632             :   pari_sp av;
     633          24 :   long lv, l = lg(D), i;
     634             :   GEN V, W, worker;
     635          24 :   check_callgen1(C, "parselect");
     636          24 :   if (!is_vec_t(typ(D))) pari_err_TYPE("parselect",D);
     637          24 :   W = cgetg(l, t_VECSMALL); av = avma;
     638          24 :   worker = snm_closure(is_entry("_parselect_worker"), mkvec(C));
     639          24 :   V = gen_parapply(worker, D);
     640       24048 :   for (lv=1, i=1; i<l; i++)
     641       24024 :     if (signe(gel(V,i))) W[lv++] = i;
     642          24 :   fixlg(W, lv);
     643          24 :   set_avma(av);
     644          24 :   return flag? W: extract_copy(D, W);
     645             : }
     646             : 
     647             : GEN
     648           0 : veccatapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     649             : {
     650           0 :   pari_sp av = avma;
     651           0 :   GEN v = vecapply(E, f, x);
     652           0 :   return lg(v) == 1? v: gerepilecopy(av, shallowconcat1(v));
     653             : }
     654             : 
     655             : static GEN
     656          14 : vecapply2(void *E, GEN (*f)(void* E, GEN x), GEN x)
     657             : {
     658             :   long i, lx;
     659          14 :   GEN y = cgetg_copy(x, &lx); y[1] = x[1];
     660          56 :   for (i=2; i<lx; i++) gel(y,i) = f(E, gel(x,i));
     661          14 :   return y;
     662             : }
     663             : static GEN
     664      123781 : vecapply1(void *E, GEN (*f)(void* E, GEN x), GEN x)
     665             : {
     666             :   long i, lx;
     667      123781 :   GEN y = cgetg_copy(x, &lx);
     668      720286 :   for (i=1; i<lx; i++) gel(y,i) = f(E, gel(x,i));
     669      123781 :   return y;
     670             : }
     671             : static GEN
     672           7 : mapapply1(void *E, GEN (*f)(void* E, GEN x), GEN x)
     673             : {
     674             :   long i, lx;
     675           7 :   GEN y = cgetg_copy(x, &lx);
     676          28 :   for (i=1; i<lx; i++)
     677             :   {
     678          21 :     GEN xi = gel(x, i);
     679          21 :     gel(y,i) = mkvec2(mkvec2(gcopy(gmael(xi, 1, 1)), f(E,gmael(xi, 1, 2))),
     680          21 :                              gcopy(gel(xi, 2)));
     681             :   }
     682           7 :   return y;
     683             : }
     684             : /* as genapply, but treat A [ t_VEC,t_COL, or t_MAT] as a t_VEC */
     685             : GEN
     686      114863 : vecapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     687             : {
     688             :   GEN y;
     689      114863 :   clone_lock(x); y = vecapply1(E,f,x);
     690      114863 :   clone_unlock_deep(x); settyp(y, t_VEC); return y;
     691             : }
     692             : GEN
     693        8904 : genapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     694             : {
     695        8904 :   long i, lx, tx = typ(x);
     696             :   GEN y, z;
     697        8904 :   if (is_scalar_t(tx)) return f(E, x);
     698        8904 :   clone_lock(x);
     699        8904 :   switch(tx) {
     700           7 :     case t_POL: y = normalizepol(vecapply2(E,f,x)); break;
     701           7 :     case t_SER:
     702           7 :       y = ser_isexactzero(x)? gcopy(x): normalize(vecapply2(E,f,x));
     703           7 :       break;
     704          28 :     case t_LIST:
     705             :       {
     706          28 :         long t = list_typ(x);
     707          28 :         z = list_data(x);
     708          28 :         if (!z)
     709           7 :           y = mklist_typ(t);
     710             :         else
     711             :         {
     712          21 :           y = cgetg(3, t_LIST);
     713          21 :           y[1] = evaltyp(t)|evallg(lg(z)-1);
     714             :           switch(t)
     715             :           {
     716          14 :           case t_LIST_RAW:
     717          14 :             list_data(y) = vecapply1(E,f,z);
     718          14 :             break;
     719           7 :           case t_LIST_MAP:
     720           7 :             list_data(y) = mapapply1(E,f,z);
     721           7 :             break;
     722             :           }
     723          28 :         }
     724             :       }
     725          28 :       break;
     726          42 :     case t_MAT:
     727          42 :       y = cgetg_copy(x, &lx);
     728         126 :       for (i = 1; i < lx; i++) gel(y,i) = vecapply1(E,f,gel(x,i));
     729          42 :       break;
     730             : 
     731        8820 :     case t_VEC: case t_COL: y = vecapply1(E,f,x); break;
     732           0 :     default:
     733           0 :       pari_err_TYPE("apply",x);
     734             :       return NULL;/*LCOV_EXCL_LINE*/
     735             :   }
     736        8904 :   clone_unlock_deep(x); return y;
     737             : }
     738             : 
     739             : GEN
     740        8904 : apply0(GEN f, GEN x)
     741             : {
     742        8904 :   check_callgen1(f, "apply");
     743        8904 :   return genapply((void *) f, gp_call, x);
     744             : }
     745             : 
     746             : GEN
     747         350 : vecselapply(void *Epred, long (*pred)(void* E, GEN x), void *Efun,
     748             :                          GEN (*fun)(void* E, GEN x), GEN A)
     749             : {
     750             :   GEN y;
     751         350 :   long i, l = lg(A), nb=1;
     752         350 :   clone_lock(A); y = cgetg(l, t_VEC);
     753      168252 :   for (i=1; i<l; i++)
     754      167902 :     if (pred(Epred, gel(A,i))) gel(y,nb++) = fun(Efun, gel(A,i));
     755         350 :   fixlg(y,nb); clone_unlock_deep(A); return y;
     756             : }
     757             : 
     758             : GEN
     759           0 : veccatselapply(void *Epred, long (*pred)(void* E, GEN x), void *Efun,
     760             :                             GEN (*fun)(void* E, GEN x), GEN A)
     761             : {
     762           0 :   pari_sp av = avma;
     763           0 :   GEN v = vecselapply(Epred, pred, Efun, fun, A);
     764           0 :   return lg(v) == 1? v: gerepilecopy(av, shallowconcat1(v));
     765             : }
     766             : 
     767             : /* suitable for gerepileupto */
     768             : GEN
     769         849 : parapply_slice_worker(GEN D, GEN worker)
     770             : {
     771             :   long l, i;
     772         849 :   GEN w = cgetg_copy(D, &l);
     773       73493 :   for (i = 1; i < l; i++) gel(w,i) = closure_callgen1(worker, gel(D,i));
     774         774 :   return w;
     775             : }
     776             : 
     777             : /* B <- {A[i] : i = r (mod m)}, 1 <= r <= m */
     778             : static void
     779       71893 : arithprogset(GEN B, GEN A, long r, long m)
     780             : {
     781       71893 :   long i, k, l = lg(A);
     782      163210 :   for (k = 1, i = r; i < l; i += m, k++) gel(B, k) = gel(A, i);
     783       71893 :   setlg(B, k);
     784       71893 : }
     785             : GEN
     786        6622 : gen_parapply_slice(GEN worker, GEN D, long mmin)
     787             : {
     788        6622 :   long l, r, n = lg(D)-1, m = minss(mmin, n), pending = 0;
     789        6622 :   GEN L = cgetg(n / m + 2, t_VEC), va = mkvec(L), V = cgetg_copy(D, &l);
     790             :   struct pari_mt pt;
     791        6622 :   mt_queue_start_lim(&pt, worker, m);
     792      143786 :   for (r = 1; r <= m || pending; r++)
     793             :   {
     794             :     long workid;
     795             :     GEN done;
     796      137164 :     if (r <= m) arithprogset(L, D, r, m);
     797      137164 :     mt_queue_submit(&pt, r, r <= m? va: NULL);
     798      137164 :     done = mt_queue_get(&pt, &workid, &pending);
     799      137164 :     if (done)
     800             :     {
     801       71893 :       long j, k, J = lg(done)-1;
     802      163210 :       for (j = 1, k = workid; j <= J; j++, k +=m) gel(V, k) = gel(done, j);
     803             :     }
     804             :   }
     805        6622 :   mt_queue_end(&pt); return V;
     806             : }
     807             : 
     808             : GEN
     809       14617 : gen_parapply(GEN worker, GEN D)
     810             : {
     811       14617 :   long l = lg(D), i, pending = 0;
     812             :   GEN W, V;
     813             :   struct pari_mt pt;
     814             : 
     815       14617 :   if (l == 1) return cgetg(1, typ(D));
     816       14617 :   W = cgetg(2, t_VEC);
     817       14617 :   V = cgetg(l, typ(D));
     818       14617 :   mt_queue_start_lim(&pt, worker, l-1);
     819       76859 :   for (i = 1; i < l || pending; i++)
     820             :   {
     821             :     long workid;
     822             :     GEN done;
     823       62242 :     if (i < l) gel(W,1) = gel(D,i);
     824       62242 :     mt_queue_submit(&pt, i, i<l? W: NULL);
     825       62242 :     done = mt_queue_get(&pt, &workid, &pending);
     826       62242 :     if (done) gel(V,workid) = done;
     827             :   }
     828       14617 :   mt_queue_end(&pt); return V;
     829             : }
     830             : 
     831             : GEN
     832         624 : parapply(GEN C, GEN D)
     833             : {
     834         624 :   pari_sp av = avma;
     835         624 :   check_callgen1(C, "parapply");
     836         624 :   if (!is_vec_t(typ(D))) pari_err_TYPE("parapply",D);
     837         624 :   return gerepileupto(av, gen_parapply(C, D));
     838             : }
     839             : 
     840             : GEN
     841          28 : genfold(void *E, GEN (*f)(void* E, GEN x, GEN y), GEN x)
     842             : {
     843          28 :   pari_sp av = avma;
     844             :   GEN z;
     845          28 :   long i, l = lg(x);
     846          28 :   if (!is_vec_t(typ(x))|| l==1  ) pari_err_TYPE("fold",x);
     847          28 :   clone_lock(x);
     848          28 :   z = gel(x,1);
     849         119 :   for (i=2; i<l; i++)
     850             :   {
     851          91 :     z = f(E,z,gel(x,i));
     852          91 :     if (gc_needed(av, 2))
     853             :     {
     854           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"fold");
     855           0 :       z = gerepilecopy(av, z);
     856             :     }
     857             :   }
     858          28 :   clone_unlock_deep(x);
     859          28 :   return gerepilecopy(av, z);
     860             : }
     861             : 
     862             : GEN
     863          28 : fold0(GEN f, GEN x)
     864             : {
     865          28 :   if (typ(f) != t_CLOSURE || closure_arity(f) < 2) pari_err_TYPE("apply",f);
     866          28 :   return genfold((void *) f, gp_call2, x);
     867             : }
     868             : /*******************************************************************/
     869             : /*                                                                 */
     870             : /*                     SCALAR-MATRIX OPERATIONS                    */
     871             : /*                                                                 */
     872             : /*******************************************************************/
     873             : GEN
     874       93433 : gtomat(GEN x)
     875             : {
     876             :   long lx, i;
     877             :   GEN y;
     878             : 
     879       93433 :   if (!x) return cgetg(1, t_MAT);
     880       93412 :   switch(typ(x))
     881             :   {
     882          28 :     case t_LIST:
     883          28 :       if (list_typ(x)==t_LIST_MAP)
     884          14 :         return maptomat(x);
     885          14 :       x = list_data(x);
     886          14 :       if (!x) return cgetg(1, t_MAT);
     887             :       /* fall through */
     888             :     case t_VEC: {
     889        2702 :       lx=lg(x); y=cgetg(lx,t_MAT);
     890        2702 :       if (lx == 1) break;
     891        2702 :       if (typ(gel(x,1)) == t_COL) {
     892        2394 :         long h = lgcols(x);
     893        5649 :         for (i=2; i<lx; i++) {
     894        3255 :           if (typ(gel(x,i)) != t_COL || lg(gel(x,i)) != h) break;
     895             :         }
     896        2394 :         if (i == lx) { /* matrix with h-1 rows */
     897        2394 :           y = cgetg(lx, t_MAT);
     898        8043 :           for (i=1 ; i<lx; i++) gel(y,i) = gcopy(gel(x,i));
     899        2394 :           return y;
     900             :         }
     901             :       }
     902         938 :       for (i=1; i<lx; i++) gel(y,i) = mkcolcopy(gel(x,i));
     903         308 :       break;
     904             :     }
     905       26307 :     case t_COL:
     906       26307 :       lx = lg(x);
     907       26307 :       if (lx == 1) return cgetg(1, t_MAT);
     908       26293 :       if (typ(gel(x,1)) == t_VEC) {
     909           7 :         long j, h = lg(gel(x,1));
     910          14 :         for (i=2; i<lx; i++) {
     911           7 :           if (typ(gel(x,i)) != t_VEC || lg(gel(x,i)) != h) break;
     912             :         }
     913           7 :         if (i == lx) { /* matrix with h cols */
     914           7 :           y = cgetg(h, t_MAT);
     915          28 :           for (j=1 ; j<h; j++) {
     916          21 :             gel(y,j) = cgetg(lx, t_COL);
     917          63 :             for (i=1; i<lx; i++) gcoeff(y,i,j) = gcopy(gmael(x,i,j));
     918             :           }
     919           7 :           return y;
     920             :         }
     921             :       }
     922       26286 :       y = mkmatcopy(x); break;
     923       55204 :     case t_MAT:
     924       55204 :       y = gcopy(x); break;
     925        8422 :     case t_QFI: case t_QFR: {
     926             :       GEN b;
     927        8422 :       y = cgetg(3,t_MAT); b = gmul2n(gel(x,2),-1);
     928        8422 :       gel(y,1) = mkcol2(icopy(gel(x,1)), b);
     929        8422 :       gel(y,2) = mkcol2(b, icopy(gel(x,3)));
     930        8422 :       break;
     931             :     }
     932         756 :     default:
     933         756 :       y = cgetg(2,t_MAT); gel(y,1) = mkcolcopy(x);
     934         756 :       break;
     935             :   }
     936       90976 :   return y;
     937             : }
     938             : 
     939             : /* create the diagonal matrix, whose diagonal is given by x */
     940             : GEN
     941         476 : diagonal(GEN x)
     942             : {
     943         476 :   long j, lx, tx = typ(x);
     944             :   GEN y;
     945             : 
     946         476 :   if (! is_matvec_t(tx)) return scalarmat(x,1);
     947         469 :   if (tx==t_MAT)
     948             :   {
     949          14 :     if (RgM_isdiagonal(x)) return gcopy(x);
     950           7 :     pari_err_TYPE("diagonal",x);
     951             :   }
     952         455 :   lx=lg(x); y=cgetg(lx,t_MAT);
     953        1428 :   for (j=1; j<lx; j++)
     954             :   {
     955         973 :     gel(y,j) = zerocol(lx-1);
     956         973 :     gcoeff(y,j,j) = gcopy(gel(x,j));
     957             :   }
     958         455 :   return y;
     959             : }
     960             : /* same, assuming x is a t_VEC/t_COL. Not memory clean. */
     961             : GEN
     962       61083 : diagonal_shallow(GEN x)
     963             : {
     964       61083 :   long j, lx = lg(x);
     965       61083 :   GEN y = cgetg(lx,t_MAT);
     966             : 
     967      163138 :   for (j=1; j<lx; j++)
     968             :   {
     969      102055 :     gel(y,j) = zerocol(lx-1);
     970      102055 :     gcoeff(y,j,j) = gel(x,j);
     971             :   }
     972       61083 :   return y;
     973             : }
     974             : 
     975             : GEN
     976         399 : zv_diagonal(GEN x)
     977             : {
     978         399 :   long j, l = lg(x), n = l-1;
     979         399 :   GEN y = cgetg(l,t_MAT);
     980             : 
     981        1344 :   for (j = 1; j < l; j++)
     982             :   {
     983         945 :     gel(y,j) = zero_Flv(n);
     984         945 :     ucoeff(y,j,j) = uel(x,j);
     985             :   }
     986         399 :   return y;
     987             : }
     988             : 
     989             : /* compute m*diagonal(d) */
     990             : GEN
     991          70 : matmuldiagonal(GEN m, GEN d)
     992             : {
     993             :   long j, lx;
     994          70 :   GEN y = cgetg_copy(m, &lx);
     995             : 
     996          70 :   if (typ(m)!=t_MAT) pari_err_TYPE("matmuldiagonal",m);
     997          70 :   if (! is_vec_t(typ(d))) pari_err_TYPE("matmuldiagonal",d);
     998          70 :   if (lg(d) != lx) pari_err_OP("operation 'matmuldiagonal'", m,d);
     999         350 :   for (j=1; j<lx; j++) gel(y,j) = RgC_Rg_mul(gel(m,j), gel(d,j));
    1000          70 :   return y;
    1001             : }
    1002             : 
    1003             : /* compute A*B assuming the result is a diagonal matrix */
    1004             : GEN
    1005           7 : matmultodiagonal(GEN A, GEN B)
    1006             : {
    1007           7 :   long i, j, hA, hB, lA = lg(A), lB = lg(B);
    1008           7 :   GEN y = matid(lB-1);
    1009             : 
    1010           7 :   if (typ(A) != t_MAT) pari_err_TYPE("matmultodiagonal",A);
    1011           7 :   if (typ(B) != t_MAT) pari_err_TYPE("matmultodiagonal",B);
    1012           7 :   hA = (lA == 1)? lB: lgcols(A);
    1013           7 :   hB = (lB == 1)? lA: lgcols(B);
    1014           7 :   if (lA != hB || lB != hA) pari_err_OP("operation 'matmultodiagonal'", A,B);
    1015          56 :   for (i=1; i<lB; i++)
    1016             :   {
    1017          49 :     GEN z = gen_0;
    1018         392 :     for (j=1; j<lA; j++) z = gadd(z, gmul(gcoeff(A,i,j),gcoeff(B,j,i)));
    1019          49 :     gcoeff(y,i,i) = z;
    1020             :   }
    1021           7 :   return y;
    1022             : }
    1023             : 
    1024             : /* [m[1,1], ..., m[l,l]], internal */
    1025             : GEN
    1026      126985 : RgM_diagonal_shallow(GEN m)
    1027             : {
    1028      126985 :   long i, lx = lg(m);
    1029      126985 :   GEN y = cgetg(lx,t_VEC);
    1030      441068 :   for (i=1; i<lx; i++) gel(y, i) = gcoeff(m,i,i);
    1031      126985 :   return y;
    1032             : }
    1033             : 
    1034             : /* same, public function */
    1035             : GEN
    1036           0 : RgM_diagonal(GEN m)
    1037             : {
    1038           0 :   long i, lx = lg(m);
    1039           0 :   GEN y = cgetg(lx,t_VEC);
    1040           0 :   for (i=1; i<lx; i++) gel(y,i) = gcopy(gcoeff(m,i,i));
    1041           0 :   return y;
    1042             : }

Generated by: LCOV version 1.13