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 - alglin3.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 19044-129ab8a) Lines: 462 505 91.5 %
Date: 2016-06-26 Functions: 41 47 87.2 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 353 441 80.0 %

           Branch data     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                 :      83747 : vecsum(GEN v)
      31                 :            : {
      32                 :      83747 :   pari_sp av = avma;
      33                 :            :   long i, l;
      34                 :            :   GEN p;
      35         [ +  + ]:      83747 :   if (!is_vec_t(typ(v)))
      36                 :          7 :     pari_err_TYPE("vecsum", v);
      37                 :      83740 :   l = lg(v);
      38         [ +  + ]:      83740 :   if (l == 1) return gen_0;
      39                 :      83733 :   p = gel(v,1);
      40         [ +  + ]:      83733 :   if (l == 2) return gcopy(p);
      41         [ +  + ]:     193764 :   for (i=2; i<l; i++)
      42                 :            :   {
      43                 :     135388 :     p = gadd(p, gel(v,i));
      44         [ -  + ]:     135388 :     if (gc_needed(av, 2))
      45                 :            :     {
      46         [ #  # ]:          0 :       if (DEBUGMEM>1) pari_warn(warnmem,"sum");
      47                 :          0 :       p = gerepileupto(av, p);
      48                 :            :     }
      49                 :            :   }
      50                 :      83740 :   return gerepileupto(av, p);
      51                 :            : }
      52                 :            : 
      53                 :            : /*******************************************************************/
      54                 :            : /*                                                                 */
      55                 :            : /*                         TRANSPOSE                               */
      56                 :            : /*                                                                 */
      57                 :            : /*******************************************************************/
      58                 :            : /* A[x0,]~ */
      59                 :            : static GEN
      60                 :    4631519 : row_transpose(GEN A, long x0)
      61                 :            : {
      62                 :    4631519 :   long i, lB = lg(A);
      63                 :    4631519 :   GEN B  = cgetg(lB, t_COL);
      64         [ +  + ]:   56717935 :   for (i=1; i<lB; i++) gel(B, i) = gcoeff(A, x0, i);
      65                 :    4631519 :   return B;
      66                 :            : }
      67                 :            : static GEN
      68                 :      16632 : row_transposecopy(GEN A, long x0)
      69                 :            : {
      70                 :      16632 :   long i, lB = lg(A);
      71                 :      16632 :   GEN B  = cgetg(lB, t_COL);
      72         [ +  + ]:     141246 :   for (i=1; i<lB; i++) gel(B, i) = gcopy(gcoeff(A, x0, i));
      73                 :      16632 :   return B;
      74                 :            : }
      75                 :            : 
      76                 :            : /* No copy*/
      77                 :            : GEN
      78                 :     860799 : shallowtrans(GEN x)
      79                 :            : {
      80                 :            :   long i, dx, lx;
      81                 :            :   GEN y;
      82   [ +  +  +  - ]:     860799 :   switch(typ(x))
      83                 :            :   {
      84                 :          7 :     case t_VEC: y = leafcopy(x); settyp(y,t_COL); break;
      85                 :       5250 :     case t_COL: y = leafcopy(x); settyp(y,t_VEC); break;
      86                 :            :     case t_MAT:
      87         [ +  + ]:     855542 :       lx = lg(x); if (lx==1) return cgetg(1,t_MAT);
      88                 :     855437 :       dx = lgcols(x); y = cgetg(dx,t_MAT);
      89         [ +  + ]:    5486956 :       for (i = 1; i < dx; i++) gel(y,i) = row_transpose(x,i);
      90                 :     855437 :       break;
      91                 :          0 :     default: pari_err_TYPE("shallowtrans",x); return NULL;
      92                 :            :   }
      93                 :     860799 :   return y;
      94                 :            : }
      95                 :            : 
      96                 :            : GEN
      97                 :       4921 : gtrans(GEN x)
      98                 :            : {
      99                 :            :   long i, dx, lx;
     100                 :            :   GEN y;
     101   [ +  +  +  + ]:       4921 :   switch(typ(x))
     102                 :            :   {
     103                 :        140 :     case t_VEC: y = gcopy(x); settyp(y,t_COL); break;
     104                 :       2905 :     case t_COL: y = gcopy(x); settyp(y,t_VEC); break;
     105                 :            :     case t_MAT:
     106         [ +  + ]:       1869 :       lx = lg(x); if (lx==1) return cgetg(1,t_MAT);
     107                 :       1862 :       dx = lgcols(x); y = cgetg(dx,t_MAT);
     108         [ +  + ]:      18494 :       for (i = 1; i < dx; i++) gel(y,i) = row_transposecopy(x,i);
     109                 :       1862 :       break;
     110                 :          7 :     default: pari_err_TYPE("gtrans",x); return NULL;
     111                 :            :   }
     112                 :       4914 :   return y;
     113                 :            : }
     114                 :            : 
     115                 :            : /*******************************************************************/
     116                 :            : /*                                                                 */
     117                 :            : /*                           EXTRACTION                            */
     118                 :            : /*                                                                 */
     119                 :            : /*******************************************************************/
     120                 :            : 
     121                 :            : static long
     122                 :        182 : str_to_long(char *s, char **pt)
     123                 :            : {
     124                 :        182 :   long a = atol(s);
     125         [ -  + ]:        182 :   while (isspace((int)*s)) s++;
     126 [ +  + ][ -  + ]:        182 :   if (*s == '-' || *s == '+') s++;
     127 [ +  + ][ -  + ]:        385 :   while (isdigit((int)*s) || isspace((int)*s)) s++;
     128                 :        182 :   *pt = s; return a;
     129                 :            : }
     130                 :            : 
     131                 :            : static int
     132                 :        112 : get_range(char *s, long *a, long *b, long *cmpl, long lx)
     133                 :            : {
     134                 :        112 :   long max = lx - 1;
     135                 :            : 
     136                 :        112 :   *a = 1; *b = max;
     137         [ +  + ]:        112 :   if (*s == '^') { *cmpl = 1; s++; } else *cmpl = 0;
     138         [ -  + ]:        112 :   if (!*s) return 0;
     139         [ +  + ]:        112 :   if (*s != '.')
     140                 :            :   {
     141                 :        105 :     *a = str_to_long(s, &s);
     142         [ +  + ]:        105 :     if (*a < 0) *a += lx;
     143 [ +  - ][ -  + ]:        105 :     if (*a<1 || *a>max) return 0;
     144                 :            :   }
     145         [ +  + ]:        112 :   if (*s == '.')
     146                 :            :   {
     147         [ -  + ]:        105 :     s++; if (*s != '.') return 0;
     148         [ -  + ]:        105 :     do s++; while (isspace((int)*s));
     149         [ +  + ]:        105 :     if (*s)
     150                 :            :     {
     151                 :         77 :       *b = str_to_long(s, &s);
     152         [ +  + ]:         77 :       if (*b < 0) *b += lx;
     153 [ +  - ][ +  + ]:         77 :       if (*b<1 || *b>max || *s) return 0;
                 [ -  + ]
     154                 :            :     }
     155                 :         98 :     return 1;
     156                 :            :   }
     157         [ -  + ]:          7 :   if (*s) return 0;
     158                 :        112 :   *b = *a; return 1;
     159                 :            : }
     160                 :            : 
     161                 :            : static int
     162                 :         35 : extract_selector_ok(long lx, GEN L)
     163                 :            : {
     164                 :            :   long i, l;
     165   [ +  +  +  +  :         35 :   switch (typ(L))
                      - ]
     166                 :            :   {
     167                 :            :     case t_INT: {
     168                 :            :       long maxj;
     169         [ -  + ]:          7 :       if (!signe(L)) return 1;
     170                 :          7 :       l = lgefint(L)-1;
     171 [ -  + ][ -  + ]:          7 :       maxj = BITS_IN_LONG - bfffo(*int_MSW(L));
         [ -  + ][ -  + ]
     172                 :          7 :       return ((l-2) * BITS_IN_LONG + maxj < lx);
     173                 :            :     }
     174                 :            :     case t_STR: {
     175                 :            :       long first, last, cmpl;
     176                 :          7 :       return get_range(GSTR(L), &first, &last, &cmpl, lx);
     177                 :            :     }
     178                 :            :     case t_VEC: case t_COL:
     179                 :         14 :       l = lg(L);
     180         [ +  + ]:         28 :       for (i=1; i<l; i++)
     181                 :            :       {
     182                 :         21 :         long j = itos(gel(L,i));
     183 [ +  + ][ -  + ]:         21 :         if (j>=lx || j<=0) return 0;
     184                 :            :       }
     185                 :          7 :       return 1;
     186                 :            :     case t_VECSMALL:
     187                 :          7 :       l = lg(L);
     188         [ +  + ]:         21 :       for (i=1; i<l; i++)
     189                 :            :       {
     190                 :         14 :         long j = L[i];
     191 [ +  - ][ -  + ]:         14 :         if (j>=lx || j<=0) return 0;
     192                 :            :       }
     193                 :          7 :       return 1;
     194                 :            :   }
     195                 :         35 :   return 0;
     196                 :            : }
     197                 :            : 
     198                 :            : GEN
     199                 :       3255 : shallowextract(GEN x, GEN L)
     200                 :            : {
     201                 :       3255 :   long i,j, tl = typ(L), tx = typ(x), lx = lg(x);
     202                 :            :   GEN y;
     203                 :            : 
     204         [ +  + ]:       3255 :   switch(tx)
     205                 :            :   {
     206                 :            :     case t_VEC:
     207                 :            :     case t_COL:
     208                 :            :     case t_MAT:
     209                 :       3248 :     case t_VECSMALL: break;
     210                 :          7 :     default: pari_err_TYPE("extract",x);
     211                 :            : 
     212                 :            :   }
     213         [ +  + ]:       3248 :   if (tl==t_INT)
     214                 :            :   { /* extract components of x as per the bits of mask L */
     215                 :            :     long k, l, ix, iy, maxj;
     216                 :            :     GEN Ld;
     217         [ +  + ]:       3066 :     if (!signe(L)) return cgetg(1,tx);
     218                 :       3059 :     y = new_chunk(lx);
     219                 :       3059 :     l = lgefint(L)-1; ix = iy = 1;
     220 [ +  + ][ +  + ]:       3059 :     maxj = BITS_IN_LONG - bfffo(*int_MSW(L));
         [ +  + ][ +  + ]
     221 [ +  + ][ +  + ]:       3059 :     if ((l-2) * BITS_IN_LONG + maxj >= lx)
     222                 :          7 :       pari_err_TYPE("vecextract [mask too large]", L);
     223         [ +  + ]:       3427 :     for (k = 2, Ld = int_LSW(L); k < l; k++, Ld = int_nextW(Ld))
     224                 :            :     {
     225                 :        375 :       ulong B = *Ld;
     226         [ +  + ]:      20439 :       for (j = 0; j < BITS_IN_LONG; j++, B >>= 1, ix++)
     227         [ +  + ]:      20064 :         if (B & 1) y[iy++] = x[ix];
     228                 :            :     }
     229                 :            :     { /* k = l */
     230                 :       3052 :       ulong B = *Ld;
     231         [ +  + ]:      27319 :       for (j = 0; j < maxj; j++, B >>= 1, ix++)
     232         [ +  + ]:      24267 :         if (B & 1) y[iy++] = x[ix];
     233                 :            :     }
     234                 :       3052 :     y[0] = evaltyp(tx) | evallg(iy);
     235                 :       3052 :     return y;
     236                 :            :   }
     237         [ +  + ]:        182 :   if (tl==t_STR)
     238                 :            :   {
     239                 :        105 :     char *s = GSTR(L);
     240                 :            :     long first, last, cmpl, d;
     241         [ +  + ]:        105 :     if (! get_range(s, &first, &last, &cmpl, lx))
     242                 :          7 :       pari_err_TYPE("vecextract [incorrect range]", L);
     243         [ -  + ]:         98 :     if (lx == 1) return cgetg(1,tx);
     244                 :         98 :     d = last - first;
     245         [ +  + ]:         98 :     if (cmpl)
     246                 :            :     {
     247         [ +  + ]:         21 :       if (d >= 0)
     248                 :            :       {
     249                 :         14 :         y = cgetg(lx - (1+d),tx);
     250         [ +  + ]:        469 :         for (j=1; j<first; j++) gel(y,j) = gel(x,j);
     251         [ +  + ]:        266 :         for (i=last+1; i<lx; i++,j++) gel(y,j) = gel(x,i);
     252                 :            :       }
     253                 :            :       else
     254                 :            :       {
     255                 :          7 :         y = cgetg(lx - (1-d),tx);
     256         [ +  + ]:         14 :         for (j=1,i=lx-1; i>first; i--,j++) gel(y,j) = gel(x,i);
     257         [ +  + ]:         14 :         for (i=last-1; i>0; i--,j++) gel(y,j) = gel(x,i);
     258                 :            :       }
     259                 :            :     }
     260                 :            :     else
     261                 :            :     {
     262         [ +  + ]:         77 :       if (d >= 0)
     263                 :            :       {
     264                 :         35 :         y = cgetg(d+2,tx);
     265         [ +  + ]:        112 :         for (i=first,j=1; i<=last; i++,j++) gel(y,j) = gel(x,i);
     266                 :            :       }
     267                 :            :       else
     268                 :            :       {
     269                 :         42 :         y = cgetg(2-d,tx);
     270         [ +  + ]:        203 :         for (i=first,j=1; i>=last; i--,j++) gel(y,j) = gel(x,i);
     271                 :            :       }
     272                 :            :     }
     273                 :         98 :     return y;
     274                 :            :   }
     275                 :            : 
     276         [ +  + ]:         77 :   if (is_vec_t(tl))
     277                 :            :   {
     278                 :         63 :     long ll=lg(L); y=cgetg(ll,tx);
     279         [ +  + ]:        140 :     for (i=1; i<ll; i++)
     280                 :            :     {
     281                 :         91 :       j = itos(gel(L,i));
     282         [ +  + ]:         91 :       if (j<=0) pari_err_COMPONENT("vecextract","<=",gen_0,stoi(j));
     283         [ +  + ]:         84 :       if (j>=lx) pari_err_COMPONENT("vecextract",">=",stoi(lx),stoi(j));
     284                 :         77 :       gel(y,i) = gel(x,j);
     285                 :            :     }
     286                 :         49 :     return y;
     287                 :            :   }
     288         [ +  + ]:         14 :   if (tl == t_VECSMALL)
     289                 :            :   {
     290                 :          7 :     long ll=lg(L); y=cgetg(ll,tx);
     291         [ +  + ]:         21 :     for (i=1; i<ll; i++)
     292                 :            :     {
     293                 :         14 :       j = L[i];
     294         [ -  + ]:         14 :       if (j<=0) pari_err_COMPONENT("vecextract","<=",gen_0,stoi(j));
     295         [ -  + ]:         14 :       if (j>=lx) pari_err_COMPONENT("vecextract",">=",stoi(lx),stoi(j));
     296                 :         14 :       gel(y,i) = gel(x,j);
     297                 :            :     }
     298                 :          7 :     return y;
     299                 :            :   }
     300                 :          7 :   pari_err_TYPE("vecextract [mask]", L);
     301                 :       3213 :   return NULL; /* not reached */
     302                 :            : }
     303                 :            : 
     304                 :            : /* does the component selector l select 0 component ? */
     305                 :            : static int
     306                 :         70 : select_0(GEN l)
     307                 :            : {
     308      [ +  +  + ]:         70 :   switch(typ(l))
     309                 :            :   {
     310                 :            :     case t_INT:
     311                 :         14 :       return (!signe(l));
     312                 :            :     case t_VEC: case t_COL: case t_VECSMALL:
     313                 :         35 :       return (lg(l) == 1);
     314                 :            :   }
     315                 :         70 :   return 0;
     316                 :            : }
     317                 :            : 
     318                 :            : GEN
     319                 :        700 : extract0(GEN x, GEN l1, GEN l2)
     320                 :            : {
     321                 :        700 :   pari_sp av = avma, av2;
     322                 :            :   GEN y;
     323         [ +  + ]:        700 :   if (! l2)
     324                 :            :   {
     325                 :        630 :     y = shallowextract(x, l1);
     326 [ +  - ][ +  + ]:        588 :     if (lg(y) == 1 || typ(y) == t_VECSMALL) return y;
     327                 :        581 :     av2 = avma;
     328                 :        581 :     y = gcopy(y);
     329                 :            :   }
     330                 :            :   else
     331                 :            :   {
     332         [ -  + ]:         70 :     if (typ(x) != t_MAT) pari_err_TYPE("extract",x);
     333                 :         70 :     y = shallowextract(x,l2);
     334         [ +  + ]:         70 :     if (select_0(l1)) { avma = av; return zeromat(0, lg(y)-1); }
     335 [ +  + ][ +  - ]:         56 :     if (lg(y) == 1 && lg(x) > 1)
     336                 :            :     {
     337         [ +  + ]:         35 :       if (!extract_selector_ok(lgcols(x), l1))
     338                 :          7 :         pari_err_TYPE("vecextract [incorrect mask]", l1);
     339                 :         28 :       avma = av; return cgetg(1, t_MAT);
     340                 :            :     }
     341                 :         21 :     y = shallowextract(shallowtrans(y), l1);
     342                 :         21 :     av2 = avma;
     343                 :         21 :     y = gtrans(y);
     344                 :            :   }
     345                 :        602 :   stackdummy(av, av2);
     346                 :        651 :   return y;
     347                 :            : }
     348                 :            : 
     349                 :            : static long
     350                 :       1753 : vecslice_parse_arg(long lA, long *y1, long *y2, long *skip)
     351                 :            : {
     352                 :       1753 :   *skip=0;
     353         [ +  + ]:       1753 :   if (*y1==LONG_MAX)
     354                 :            :   {
     355         [ +  + ]:        105 :     if (*y2!=LONG_MAX)
     356                 :            :     {
     357         [ +  + ]:         63 :       if (*y2<0) *y2 += lA;
     358 [ +  - ][ +  - ]:         63 :       if (*y2<0 || *y2==LONG_MAX || *y2>=lA)
                 [ -  + ]
     359                 :          0 :         pari_err_DIM("_[..]");
     360                 :         63 :       *skip=*y2;
     361                 :            :     }
     362                 :        105 :     *y1 = 1; *y2 = lA-1;
     363                 :            :   }
     364         [ +  + ]:       1648 :   else if (*y2==LONG_MAX) *y2 = *y1;
     365         [ +  + ]:       1753 :   if (*y1<=0) *y1 += lA;
     366         [ +  + ]:       1753 :   if (*y2<0) *y2 += lA;
     367 [ +  + ][ +  - ]:       1753 :   if (*y1<=0 || *y1>*y2+1 || *y2>=lA) pari_err_DIM("_[..]");
                 [ +  + ]
     368                 :       1739 :   return *y2 - *y1 + 2 - !!*skip;
     369                 :            : }
     370                 :            : 
     371                 :            : static GEN
     372                 :       1809 : vecslice_i(GEN A, long t, long lB, long y1, long skip)
     373                 :            : {
     374                 :       1809 :   GEN B = cgetg(lB, t);
     375                 :            :   long i;
     376         [ +  + ]:     354094 :   for (i=1; i<lB; i++, y1++)
     377                 :            :   {
     378         [ +  + ]:     352285 :     if (y1 == skip) { i--; continue; }
     379                 :     352222 :     gel(B,i) = gcopy(gel(A,y1));
     380                 :            :   }
     381                 :       1809 :   return B;
     382                 :            : }
     383                 :            : 
     384                 :            : static GEN
     385                 :         14 : rowslice_i(GEN A, long lB, long x1, long y1, long skip)
     386                 :            : {
     387                 :         14 :   GEN B = cgetg(lB, t_VEC);
     388                 :            :   long i;
     389         [ +  + ]:         77 :   for (i=1; i<lB; i++, y1++)
     390                 :            :   {
     391         [ +  + ]:         63 :     if (y1 == skip) { i--; continue; }
     392                 :         56 :     gel(B,i) = gcopy(gcoeff(A,x1,y1));
     393                 :            :   }
     394                 :         14 :   return B;
     395                 :            : }
     396                 :            : 
     397                 :            : static GEN
     398                 :          0 : rowsmallslice_i(GEN A, long lB, long x1, long y1, long skip)
     399                 :            : {
     400                 :          0 :   GEN B = cgetg(lB, t_VECSMALL);
     401                 :            :   long i;
     402         [ #  # ]:          0 :   for (i=1; i<lB; i++, y1++)
     403                 :            :   {
     404         [ #  # ]:          0 :     if (y1 == skip) { i--; continue; }
     405                 :          0 :     B[i] = coeff(A,x1,y1);
     406                 :            :   }
     407                 :          0 :   return B;
     408                 :            : }
     409                 :            : 
     410                 :            : static GEN
     411                 :         28 : vecsmallslice_i(GEN A, long t, long lB, long y1, long skip)
     412                 :            : {
     413                 :         28 :   GEN B = cgetg(lB, t);
     414                 :            :   long i;
     415         [ +  + ]:        126 :   for (i=1; i<lB; i++, y1++)
     416                 :            :   {
     417         [ +  + ]:         98 :     if (y1 == skip) { i--; continue; }
     418                 :         91 :     B[i] = A[y1];
     419                 :            :   }
     420                 :         28 :   return B;
     421                 :            : }
     422                 :            : GEN
     423                 :       1515 : vecslice0(GEN A, long y1, long y2)
     424                 :            : {
     425                 :       1515 :   long skip, lB, t = typ(A);
     426                 :       1515 :   lB = vecslice_parse_arg(lg(A), &y1, &y2, &skip);
     427      [ +  +  + ]:       1501 :   switch(t)
     428                 :            :   {
     429                 :            :     case t_VEC: case t_COL:
     430                 :       1466 :       return vecslice_i(A, t,lB,y1,skip);
     431                 :            :     case t_VECSMALL:
     432                 :         28 :       return vecsmallslice_i(A, t,lB,y1,skip);
     433                 :            :     default:
     434                 :          7 :       pari_err_TYPE("_[_.._]",A);
     435                 :       1494 :       return NULL;
     436                 :            :   }
     437                 :            : }
     438                 :            : 
     439                 :            : GEN
     440                 :        126 : matslice0(GEN A, long x1, long x2, long y1, long y2)
     441                 :            : {
     442                 :            :   GEN B;
     443                 :        126 :   long i, lB, lA = lg(A), t, skip, rskip, rlB;
     444 [ +  + ][ +  + ]:        126 :   long is_col = y1!=LONG_MAX && y2==LONG_MAX;
     445 [ +  + ][ +  + ]:        126 :   long is_row = x1!=LONG_MAX && x2==LONG_MAX;
     446                 :            :   GEN (*slice)(GEN A, long t, long lB, long y1, long skip);
     447         [ -  + ]:        126 :   if (typ(A)!=t_MAT) pari_err_TYPE("_[_.._,_.._]",A);
     448                 :        126 :   lB = vecslice_parse_arg(lA, &y1, &y2, &skip);
     449         [ +  + ]:        126 :   if (is_col) return vecslice0(gel(A, y1), x1, x2);
     450                 :            : 
     451                 :            :   /* lA > 1 */
     452                 :        112 :   rlB = vecslice_parse_arg(lg(gel(A,1)), &x1, &x2, &rskip);
     453                 :        112 :   t = typ(gel(A,1));
     454 [ +  + ][ +  - ]:        112 :   if (is_row) return t == t_COL ? rowslice_i(A, lB, x1, y1, skip):
     455                 :          0 :                                   rowsmallslice_i(A, lB, x1, y1, skip);
     456         [ +  - ]:         98 :   slice = t == t_COL? &vecslice_i: &vecsmallslice_i;
     457                 :            : 
     458                 :         98 :   B = cgetg(lB, t_MAT);
     459         [ +  + ]:        455 :   for (i=1; i<lB; i++, y1++)
     460                 :            :   {
     461         [ +  + ]:        357 :     if (y1 == skip) { i--; continue; }
     462                 :        343 :     gel(B,i) = slice(gel(A,y1),t,rlB, x1, rskip);
     463                 :            :   }
     464                 :        126 :   return B;
     465                 :            : }
     466                 :            : 
     467                 :            : GEN
     468                 :       4819 : vecrange(GEN a, GEN b)
     469                 :            : {
     470                 :            :   GEN y;
     471                 :            :   long i, l;
     472         [ +  + ]:       4819 :   if (typ(a)!=t_INT) pari_err_TYPE("[_.._]",a);
     473         [ +  + ]:       4812 :   if (typ(b)!=t_INT) pari_err_TYPE("[_.._]",b);
     474         [ +  + ]:       4805 :   if (cmpii(a,b)>0) return cgetg(1,t_VEC);
     475                 :       4798 :   l = itos(subii(b,a))+1;
     476                 :       4798 :   a = setloop(a);
     477                 :       4798 :   y = cgetg(l+1, t_VEC);
     478         [ +  + ]:     317628 :   for (i=1; i<=l; a = incloop(a), i++)
     479                 :     312830 :     gel(y,i) = icopy(a);
     480                 :       4805 :   return y;
     481                 :            : }
     482                 :            : 
     483                 :            : GEN
     484                 :          0 : vecrangess(long a, long b)
     485                 :            : {
     486                 :            :   GEN y;
     487                 :            :   long i, l;
     488         [ #  # ]:          0 :   if (a>b) return cgetg(1,t_VEC);
     489                 :          0 :   l = b-a+1;
     490                 :          0 :   y = cgetg(l+1, t_VEC);
     491         [ #  # ]:          0 :   for (i=1; i<=l; a++, i++)
     492                 :          0 :     gel(y,i) = stoi(a);
     493                 :          0 :   return y;
     494                 :            : }
     495                 :            : 
     496                 :            : GEN
     497                 :       4669 : genindexselect(void *E, long (*f)(void* E, GEN x), GEN A)
     498                 :            : {
     499                 :            :   long l, i, lv;
     500                 :            :   GEN v, z;
     501                 :            :   pari_sp av;
     502                 :       4669 :   clone_lock(A);
     503      [ +  +  + ]:       4669 :   switch(typ(A))
     504                 :            :   {
     505                 :            :     case t_LIST:
     506                 :          7 :       z = list_data(A);
     507         [ +  - ]:          7 :       l = z? lg(z): 1;
     508                 :          7 :       break;
     509                 :            :     case t_VEC: case t_COL: case t_MAT:
     510                 :       4655 :       l = lg(A);
     511                 :       4655 :       z = A;
     512                 :       4655 :       break;
     513                 :            :     default:
     514                 :          7 :       pari_err_TYPE("select",A);
     515                 :          0 :       return NULL;/*not reached*/
     516                 :            :   }
     517                 :       4662 :   v = cgetg(l, t_VECSMALL);
     518                 :       4662 :   av = avma;
     519         [ +  + ]:      51296 :   for (i = lv = 1; i < l; i++) {
     520         [ +  + ]:      46634 :     if (f(E, gel(z,i))) v[lv++] = i;
     521                 :      46634 :     avma = av;
     522                 :            :   }
     523                 :       4662 :   clone_unlock(A); fixlg(v, lv); return v;
     524                 :            : }
     525                 :            : static GEN
     526                 :       4657 : extract_copy(GEN A, GEN v)
     527                 :            : {
     528                 :       4657 :   long i, l = lg(v);
     529                 :       4657 :   GEN B = cgetg(l, typ(A));
     530         [ +  + ]:       9437 :   for (i = 1; i < l; i++) gel(B,i) = gcopy(gel(A,v[i]));
     531                 :       4657 :   return B;
     532                 :            : }
     533                 :            : /* as genselect, but treat A [ t_VEC,t_COL, or t_MAT] as a t_VEC */
     534                 :            : GEN
     535                 :          0 : vecselect(void *E, long (*f)(void* E, GEN x), GEN A)
     536                 :            : {
     537                 :            :   GEN v;
     538                 :          0 :   clone_lock(A);
     539                 :          0 :   v = genindexselect(E, f, A);
     540                 :          0 :   A = extract_copy(A, v); settyp(A, t_VEC);
     541                 :          0 :   clone_unlock(A); return A;
     542                 :            : }
     543                 :            : GEN
     544                 :       4662 : genselect(void *E, long (*f)(void* E, GEN x), GEN A)
     545                 :            : {
     546                 :            :   GEN y, z, v;/* v left on stack for efficiency */
     547                 :       4662 :   clone_lock(A);
     548      [ +  +  + ]:       4662 :   switch(typ(A))
     549                 :            :   {
     550                 :            :     case t_LIST:
     551                 :         14 :       z = list_data(A);
     552         [ -  + ]:         14 :       if (!z) y = listcreate();
     553                 :            :       else
     554                 :            :       {
     555                 :            :         GEN B;
     556                 :         14 :         y = cgetg(3, t_LIST);
     557                 :         14 :         v = genindexselect(E, f, z);
     558                 :         14 :         B = extract_copy(z, v);
     559                 :         14 :         y[1] = lg(B)-1;
     560                 :         14 :         list_data(y) = B;
     561                 :            :       }
     562                 :         14 :       break;
     563                 :            :     case t_VEC: case t_COL: case t_MAT:
     564                 :       4641 :       v = genindexselect(E, f, A);
     565                 :       4641 :       y = extract_copy(A, v);
     566                 :       4641 :       break;
     567                 :            :     default:
     568                 :          7 :       pari_err_TYPE("select",A);
     569                 :          0 :       return NULL;/*not reached*/
     570                 :            :   }
     571                 :       4655 :   clone_unlock(A); return y;
     572                 :            : }
     573                 :            : 
     574                 :            : static void
     575                 :       5431 : check_callgen1(GEN f, const char *s)
     576                 :            : {
     577 [ +  - ][ +  - ]:       5431 :   if (typ(f) != t_CLOSURE || closure_is_variadic(f)  || closure_arity(f) < 1)
                 [ -  + ]
     578                 :          0 :     pari_err_TYPE(s, f);
     579                 :       5431 : }
     580                 :            : 
     581                 :            : GEN
     582                 :       4676 : select0(GEN f, GEN x, long flag)
     583                 :            : {
     584                 :       4676 :   check_callgen1(f, "select");
     585      [ +  +  - ]:       4676 :   switch(flag)
     586                 :            :   {
     587                 :       4662 :     case 0: return genselect((void *) f, gp_callbool, x);
     588                 :         14 :     case 1: return genindexselect((void *) f, gp_callbool, x);
     589                 :          0 :     default: pari_err_FLAG("select");
     590                 :       4662 :              return NULL;/*not reached*/
     591                 :            :   }
     592                 :            : }
     593                 :            : 
     594                 :            : GEN
     595                 :          4 : parselect(GEN C, GEN D, long flag)
     596                 :            : {
     597                 :            :   pari_sp av, av2;
     598                 :          4 :   long lv, l = lg(D), i, pending = 0, workid;
     599                 :            :   GEN V, worker, done;
     600                 :            :   struct pari_mt pt;
     601                 :          4 :   check_callgen1(C, "parselect");
     602         [ -  + ]:          4 :   if (!is_vec_t(typ(D))) pari_err_TYPE("parapply",D);
     603                 :          4 :   V = cgetg(l, t_VECSMALL); av = avma;
     604                 :          4 :   worker = strtoclosure("_parapply_worker", 1, C);
     605                 :          4 :   av2 = avma;
     606                 :          4 :   mt_queue_start(&pt, worker);
     607 [ +  + ][ +  + ]:       4068 :   for (i=1; i<l || pending; i++)
     608                 :            :   {
     609         [ +  + ]:       4064 :     mt_queue_submit(&pt, i, i<l? mkvec(gel(D,i)): NULL);
     610                 :       4064 :     done = mt_queue_get(&pt, &workid, &pending);
     611         [ +  + ]:       4064 :     if (done) V[workid] = !gequal0(done);
     612                 :       4064 :     avma = av2;
     613                 :            :   }
     614                 :          4 :   mt_queue_end(&pt);
     615                 :          4 :   avma = av;
     616         [ +  + ]:       4008 :   for (lv=1, i=1; i<l; i++)
     617         [ +  + ]:       4004 :     if (V[i]) V[lv++]=i;
     618                 :          4 :   fixlg(V, lv);
     619         [ +  + ]:          4 :   return flag? V: extract_copy(D, V);
     620                 :            : }
     621                 :            : 
     622                 :            : GEN
     623                 :          0 : veccatapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     624                 :            : {
     625                 :          0 :   pari_sp av = avma;
     626                 :          0 :   GEN v = vecapply(E, f, x);
     627         [ #  # ]:          0 :   return lg(v) == 1? v: gerepilecopy(av, shallowconcat1(v));
     628                 :            : }
     629                 :            : 
     630                 :            : static GEN
     631                 :         14 : vecapply2(void *E, GEN (*f)(void* E, GEN x), GEN x)
     632                 :            : {
     633                 :            :   long i, lx;
     634                 :         14 :   GEN y = cgetg_copy(x, &lx); y[1] = x[1];
     635         [ +  + ]:         56 :   for (i=2; i<lx; i++) gel(y,i) = f(E, gel(x,i));
     636                 :         14 :   return y;
     637                 :            : }
     638                 :            : static GEN
     639                 :     111314 : vecapply1(void *E, GEN (*f)(void* E, GEN x), GEN x)
     640                 :            : {
     641                 :            :   long i, lx;
     642                 :     111314 :   GEN y = cgetg_copy(x, &lx);
     643         [ +  + ]:     745556 :   for (i=1; i<lx; i++) gel(y,i) = f(E, gel(x,i));
     644                 :     111314 :   return y;
     645                 :            : }
     646                 :            : static GEN
     647                 :          7 : mapapply1(void *E, GEN (*f)(void* E, GEN x), GEN x)
     648                 :            : {
     649                 :            :   long i, lx;
     650                 :          7 :   GEN y = cgetg_copy(x, &lx);
     651         [ +  + ]:         28 :   for (i=1; i<lx; i++)
     652                 :            :   {
     653                 :         21 :     GEN xi = gel(x, i);
     654                 :         21 :     gel(y,i) = mkvec2(mkvec2(gcopy(gmael(xi, 1, 1)), f(E,gmael(xi, 1, 2))),
     655                 :         21 :                              gcopy(gel(xi, 2)));
     656                 :            :   }
     657                 :          7 :   return y;
     658                 :            : }
     659                 :            : /* as genapply, but treat A [ t_VEC,t_COL, or t_MAT] as a t_VEC */
     660                 :            : GEN
     661                 :     110565 : vecapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     662                 :            : {
     663                 :            :   GEN y;
     664                 :     110565 :   clone_lock(x); y = vecapply1(E,f,x);
     665                 :     110565 :   clone_unlock(x); settyp(y, t_VEC); return y;
     666                 :            : }
     667                 :            : GEN
     668                 :        749 : genapply(void *E, GEN (*f)(void* E, GEN x), GEN x)
     669                 :            : {
     670                 :        749 :   long i, lx, tx = typ(x);
     671                 :            :   GEN y, z;
     672         [ -  + ]:        749 :   if (is_scalar_t(tx)) return f(E, x);
     673                 :        749 :   clone_lock(x);
     674   [ +  +  +  +  :        749 :   switch(tx) {
                   +  - ]
     675                 :          7 :     case t_POL: y = normalizepol(vecapply2(E,f,x)); break;
     676                 :            :     case t_SER:
     677         [ -  + ]:          7 :       y = ser_isexactzero(x)? gcopy(x): normalize(vecapply2(E,f,x));
     678                 :          7 :       break;
     679                 :            :     case t_LIST:
     680                 :            :       {
     681                 :         28 :         long t = list_typ(x);
     682                 :         28 :         z = list_data(x);
     683         [ +  + ]:         28 :         if (!z)
     684                 :          7 :           y = listcreate_typ(t);
     685                 :            :         else
     686                 :            :         {
     687                 :         21 :           y = cgetg(3, t_LIST);
     688                 :         21 :           y[1] = evaltyp(t)|evallg(lg(z)-1);
     689      [ +  +  - ]:         21 :           switch(t)
     690                 :            :           {
     691                 :            :           case t_LIST_RAW:
     692                 :         14 :             list_data(y) = vecapply1(E,f,z);
     693                 :         14 :             break;
     694                 :            :           case t_LIST_MAP:
     695                 :          7 :             list_data(y) = mapapply1(E,f,z);
     696                 :          7 :             break;
     697                 :            :           }
     698                 :            :         }
     699                 :            :       }
     700                 :         28 :       break;
     701                 :            :     case t_MAT:
     702                 :         28 :       y = cgetg_copy(x, &lx);
     703         [ +  + ]:         84 :       for (i = 1; i < lx; i++) gel(y,i) = vecapply1(E,f,gel(x,i));
     704                 :         28 :       break;
     705                 :            : 
     706                 :        679 :     case t_VEC: case t_COL: y = vecapply1(E,f,x); break;
     707                 :            :     default:
     708                 :          0 :       pari_err_TYPE("apply",x); return NULL;/*not reached*/
     709                 :            :   }
     710                 :        749 :   clone_unlock(x); return y;
     711                 :            : }
     712                 :            : 
     713                 :            : GEN
     714                 :        749 : apply0(GEN f, GEN x)
     715                 :            : {
     716                 :        749 :   check_callgen1(f, "apply");
     717                 :        749 :   return genapply((void *) f, gp_call, x);
     718                 :            : }
     719                 :            : 
     720                 :            : GEN
     721                 :        252 : vecselapply(void *Epred, long (*pred)(void* E, GEN x), void *Efun,
     722                 :            :                          GEN (*fun)(void* E, GEN x), GEN A)
     723                 :            : {
     724                 :            :   GEN y;
     725                 :        252 :   long i, l = lg(A), nb=1;
     726                 :        252 :   clone_lock(A); y = cgetg(l, t_VEC);
     727         [ +  + ]:     167139 :   for (i=1; i<l; i++)
     728         [ +  + ]:     166887 :     if (pred(Epred, gel(A,i))) gel(y,nb++) = fun(Efun, gel(A,i));
     729                 :        252 :   fixlg(y,nb); clone_unlock(A); return y;
     730                 :            : }
     731                 :            : 
     732                 :            : GEN
     733                 :          0 : veccatselapply(void *Epred, long (*pred)(void* E, GEN x), void *Efun,
     734                 :            :                             GEN (*fun)(void* E, GEN x), GEN A)
     735                 :            : {
     736                 :          0 :   pari_sp av = avma;
     737                 :          0 :   GEN v = vecselapply(Epred, pred, Efun, fun, A);
     738         [ #  # ]:          0 :   return lg(v) == 1? v: gerepilecopy(av, shallowconcat1(v));
     739                 :            : }
     740                 :            : 
     741                 :            : GEN
     742                 :       3993 : parapply_worker(GEN d, GEN C)
     743                 :            : {
     744                 :       3993 :   return closure_callgen1(C, d);
     745                 :            : }
     746                 :            : 
     747                 :            : GEN
     748                 :          2 : parapply(GEN C, GEN D)
     749                 :            : {
     750                 :          2 :   pari_sp av = avma;
     751                 :          2 :   long l = lg(D), i, pending = 0, workid;
     752                 :            :   GEN V, worker, done;
     753                 :            :   struct pari_mt pt;
     754                 :          2 :   check_callgen1(C, "parapply");
     755         [ -  + ]:          2 :   if (!is_vec_t(typ(D))) pari_err_TYPE("parapply",D);
     756                 :          2 :   worker = strtoclosure("_parapply_worker", 1, C);
     757                 :          2 :   V = cgetg(l, typ(D));
     758                 :          2 :   mt_queue_start(&pt, worker);
     759 [ +  + ][ +  + ]:         10 :   for (i=1; i<l || pending; i++)
     760                 :            :   {
     761         [ +  + ]:          8 :     mt_queue_submit(&pt, i, i<l? mkvec(gel(D,i)): NULL);
     762                 :          8 :     done = mt_queue_get(&pt, &workid, &pending);
     763         [ +  + ]:          8 :     if (done) gel(V,workid) = done;
     764                 :            :   }
     765                 :          2 :   mt_queue_end(&pt);
     766                 :          2 :   return gerepilecopy(av, V);
     767                 :            : }
     768                 :            : 
     769                 :            : GEN
     770                 :         28 : genfold(void *E, GEN (*f)(void* E, GEN x, GEN y), GEN x)
     771                 :            : {
     772                 :         28 :   pari_sp av = avma;
     773                 :            :   GEN z;
     774                 :         28 :   long i, l = lg(x);
     775 [ +  - ][ -  + ]:         28 :   if (!is_vec_t(typ(x))|| l==1  ) pari_err_TYPE("fold",x);
     776                 :         28 :   clone_lock(x);
     777                 :         28 :   z = gel(x,1);
     778         [ +  + ]:        119 :   for (i=2; i<l; i++)
     779                 :            :   {
     780                 :         91 :     z = f(E,z,gel(x,i));
     781         [ -  + ]:         91 :     if (gc_needed(av, 2))
     782                 :            :     {
     783         [ #  # ]:          0 :       if (DEBUGMEM>1) pari_warn(warnmem,"fold");
     784                 :          0 :       z = gerepilecopy(av, z);
     785                 :            :     }
     786                 :            :   }
     787                 :         28 :   clone_unlock(x);
     788                 :         28 :   return gerepilecopy(av, z);
     789                 :            : }
     790                 :            : 
     791                 :            : GEN
     792                 :         28 : fold0(GEN f, GEN x)
     793                 :            : {
     794 [ +  - ][ -  + ]:         28 :   if (typ(f) != t_CLOSURE || closure_arity(f) < 2) pari_err_TYPE("apply",f);
     795                 :         28 :   return genfold((void *) f, gp_call2, x);
     796                 :            : }
     797                 :            : /*******************************************************************/
     798                 :            : /*                                                                 */
     799                 :            : /*                     SCALAR-MATRIX OPERATIONS                    */
     800                 :            : /*                                                                 */
     801                 :            : /*******************************************************************/
     802                 :            : GEN
     803                 :      37557 : gtomat(GEN x)
     804                 :            : {
     805                 :            :   long lx, i;
     806                 :            :   GEN y;
     807                 :            : 
     808         [ -  + ]:      37557 :   if (!x) return cgetg(1, t_MAT);
     809   [ +  +  +  +  :      37557 :   switch(typ(x))
                   +  + ]
     810                 :            :   {
     811                 :            :     case t_LIST:
     812         [ +  + ]:         21 :       if (list_typ(x)==t_LIST_MAP)
     813                 :          7 :         return maptomat(x);
     814                 :         14 :       x = list_data(x);
     815         [ +  + ]:         14 :       if (!x) return cgetg(1, t_MAT);
     816                 :            :       /* fall through */
     817                 :            :     case t_VEC: {
     818                 :        280 :       lx=lg(x); y=cgetg(lx,t_MAT);
     819         [ -  + ]:        280 :       if (lx == 1) break;
     820         [ +  + ]:        280 :       if (typ(gel(x,1)) == t_COL) {
     821                 :        126 :         long h = lgcols(x);
     822         [ +  + ]:        854 :         for (i=2; i<lx; i++) {
     823 [ +  - ][ +  - ]:        728 :           if (typ(gel(x,i)) != t_COL || lg(gel(x,i)) != h) break;
     824                 :            :         }
     825         [ +  - ]:        126 :         if (i == lx) { /* matrix with h-1 rows */
     826                 :        126 :           y = cgetg(lx, t_MAT);
     827         [ +  + ]:        980 :           for (i=1 ; i<lx; i++) gel(y,i) = gcopy(gel(x,i));
     828                 :        126 :           return y;
     829                 :            :         }
     830                 :            :       }
     831         [ +  + ]:        518 :       for (i=1; i<lx; i++) gel(y,i) = mkcolcopy(gel(x,i));
     832                 :        154 :       break;
     833                 :            :     }
     834                 :            :     case t_COL:
     835                 :       7392 :       lx = lg(x);
     836         [ -  + ]:       7392 :       if (lx == 1) return cgetg(1, t_MAT);
     837         [ +  + ]:       7392 :       if (typ(gel(x,1)) == t_VEC) {
     838                 :          7 :         long j, h = lg(gel(x,1));
     839         [ +  + ]:         14 :         for (i=2; i<lx; i++) {
     840 [ +  - ][ +  - ]:          7 :           if (typ(gel(x,i)) != t_VEC || lg(gel(x,i)) != h) break;
     841                 :            :         }
     842         [ +  - ]:          7 :         if (i == lx) { /* matrix with h cols */
     843                 :          7 :           y = cgetg(h, t_MAT);
     844         [ +  + ]:         28 :           for (j=1 ; j<h; j++) {
     845                 :         21 :             gel(y,j) = cgetg(lx, t_COL);
     846         [ +  + ]:         63 :             for (i=1; i<lx; i++) gcoeff(y,i,j) = gcopy(gmael(x,i,j));
     847                 :            :           }
     848                 :          7 :           return y;
     849                 :            :         }
     850                 :            :       }
     851                 :       7385 :       y = mkmatcopy(x); break;
     852                 :            :     case t_MAT:
     853                 :      20211 :       y = gcopy(x); break;
     854                 :            :     case t_QFI: case t_QFR: {
     855                 :            :       GEN b;
     856                 :       8428 :       y = cgetg(3,t_MAT); b = gmul2n(gel(x,2),-1);
     857                 :       8428 :       gel(y,1) = mkcol2(icopy(gel(x,1)), b);
     858                 :       8428 :       gel(y,2) = mkcol2(b, icopy(gel(x,3)));
     859                 :       8428 :       break;
     860                 :            :     }
     861                 :            :     default:
     862                 :       1232 :       y = cgetg(2,t_MAT); gel(y,1) = mkcolcopy(x);
     863                 :       1232 :       break;
     864                 :            :   }
     865                 :      37557 :   return y;
     866                 :            : }
     867                 :            : 
     868                 :            : /* create the diagonal matrix, whose diagonal is given by x */
     869                 :            : GEN
     870                 :        476 : diagonal(GEN x)
     871                 :            : {
     872                 :        476 :   long j, lx, tx = typ(x);
     873                 :            :   GEN y;
     874                 :            : 
     875         [ +  + ]:        476 :   if (! is_matvec_t(tx)) return scalarmat(x,1);
     876         [ +  + ]:        469 :   if (tx==t_MAT)
     877                 :            :   {
     878         [ +  + ]:         14 :     if (RgM_isdiagonal(x)) return gcopy(x);
     879                 :          7 :     pari_err_TYPE("diagonal",x);
     880                 :            :   }
     881                 :        455 :   lx=lg(x); y=cgetg(lx,t_MAT);
     882         [ +  + ]:       1127 :   for (j=1; j<lx; j++)
     883                 :            :   {
     884                 :        672 :     gel(y,j) = zerocol(lx-1);
     885                 :        672 :     gcoeff(y,j,j) = gcopy(gel(x,j));
     886                 :            :   }
     887                 :        469 :   return y;
     888                 :            : }
     889                 :            : /* same, assuming x is a t_VEC/t_COL. Not memory clean. */
     890                 :            : GEN
     891                 :      18305 : diagonal_shallow(GEN x)
     892                 :            : {
     893                 :      18305 :   long j, lx = lg(x);
     894                 :      18305 :   GEN y = cgetg(lx,t_MAT);
     895                 :            : 
     896         [ +  + ]:      46571 :   for (j=1; j<lx; j++)
     897                 :            :   {
     898                 :      28266 :     gel(y,j) = zerocol(lx-1);
     899                 :      28266 :     gcoeff(y,j,j) = gel(x,j);
     900                 :            :   }
     901                 :      18305 :   return y;
     902                 :            : }
     903                 :            : 
     904                 :            : /* compute m*diagonal(d) */
     905                 :            : GEN
     906                 :          7 : matmuldiagonal(GEN m, GEN d)
     907                 :            : {
     908                 :            :   long j, lx;
     909                 :          7 :   GEN y = cgetg_copy(m, &lx);
     910                 :            : 
     911         [ -  + ]:          7 :   if (typ(m)!=t_MAT) pari_err_TYPE("matmuldiagonal",m);
     912         [ -  + ]:          7 :   if (! is_vec_t(typ(d))) pari_err_TYPE("matmuldiagonal",d);
     913         [ -  + ]:          7 :   if (lg(d) != lx) pari_err_OP("operation 'matmuldiagonal'", m,d);
     914         [ +  + ]:         56 :   for (j=1; j<lx; j++) gel(y,j) = RgC_Rg_mul(gel(m,j), gel(d,j));
     915                 :          7 :   return y;
     916                 :            : }
     917                 :            : 
     918                 :            : /* compute A*B assuming the result is a diagonal matrix */
     919                 :            : GEN
     920                 :          7 : matmultodiagonal(GEN A, GEN B)
     921                 :            : {
     922                 :          7 :   long i, j, hA, hB, lA = lg(A), lB = lg(B);
     923                 :          7 :   GEN y = matid(lB-1);
     924                 :            : 
     925         [ -  + ]:          7 :   if (typ(A) != t_MAT) pari_err_TYPE("matmultodiagonal",A);
     926         [ -  + ]:          7 :   if (typ(B) != t_MAT) pari_err_TYPE("matmultodiagonal",B);
     927         [ +  - ]:          7 :   hA = (lA == 1)? lB: lgcols(A);
     928         [ +  - ]:          7 :   hB = (lB == 1)? lA: lgcols(B);
     929 [ +  - ][ -  + ]:          7 :   if (lA != hB || lB != hA) pari_err_OP("operation 'matmultodiagonal'", A,B);
     930         [ +  + ]:         56 :   for (i=1; i<lB; i++)
     931                 :            :   {
     932                 :         49 :     GEN z = gen_0;
     933         [ +  + ]:        392 :     for (j=1; j<lA; j++) z = gadd(z, gmul(gcoeff(A,i,j),gcoeff(B,j,i)));
     934                 :         49 :     gcoeff(y,i,i) = z;
     935                 :            :   }
     936                 :          7 :   return y;
     937                 :            : }
     938                 :            : 
     939                 :            : /* [m[1,1], ..., m[l,l]], internal */
     940                 :            : GEN
     941                 :     129329 : RgM_diagonal_shallow(GEN m)
     942                 :            : {
     943                 :     129329 :   long i, lx = lg(m);
     944                 :     129329 :   GEN y = cgetg(lx,t_VEC);
     945         [ +  + ]:     583082 :   for (i=1; i<lx; i++) gel(y, i) = gcoeff(m,i,i);
     946                 :     129329 :   return y;
     947                 :            : }
     948                 :            : 
     949                 :            : /* same, public function */
     950                 :            : GEN
     951                 :          0 : RgM_diagonal(GEN m)
     952                 :            : {
     953                 :          0 :   long i, lx = lg(m);
     954                 :          0 :   GEN y = cgetg(lx,t_VEC);
     955         [ #  # ]:          0 :   for (i=1; i<lx; i++) gel(y,i) = gcopy(gcoeff(m,i,i));
     956                 :          0 :   return y;
     957                 :            : }
     958                 :            : 
     959                 :            : 

Generated by: LCOV version 1.9