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 - Ser.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 24988-2584e74448) Lines: 127 127 100.0 %
Date: 2020-01-26 05:57:03 Functions: 18 18 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000, 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             : #include "pari.h"
      15             : #include "paripriv.h"
      16             : /*******************************************************************/
      17             : /*                                                                 */
      18             : /*                     Conversion --> t_SER                        */
      19             : /*                                                                 */
      20             : /*******************************************************************/
      21             : static GEN
      22     1450833 : RgX_to_ser_i(GEN x, long l, long lx, long v, int copy)
      23             : {
      24             :   GEN y;
      25             :   long i;
      26     1450833 :   if (lx == 2) return zeroser(varn(x), l-2);
      27     1450735 :   if (l < 2) pari_err_BUG("RgX_to_ser (l < 2)");
      28     1450735 :   y = cgetg(l,t_SER); y[1] = x[1];
      29             :   /* e.g. Mod(0,3) * x^0 */
      30     1450735 :   if (v == LONG_MAX) { v = 1; lx = 3; } else { x += v; lx = minss(lx-v, l); }
      31     1450735 :   setvalp(y, v);
      32     1450735 :   if (copy)
      33          35 :     for (i = 2; i <lx; i++) gel(y,i) = gcopy(gel(x,i));
      34             :   else
      35     1450700 :     for (i = 2; i <lx; i++) gel(y,i) = gel(x,i);
      36     1450735 :   for (     ; i < l; i++) gel(y,i) = gen_0;
      37     1450735 :   return normalize(y);
      38             : }
      39             : /* enlarge/truncate t_POL x to a t_SER with lg l */
      40             : GEN
      41     1040129 : RgX_to_ser(GEN x, long l) { return RgX_to_ser_i(x, l, lg(x), RgX_val(x), 0); }
      42             : GEN
      43      410669 : RgX_to_ser_inexact(GEN x, long l)
      44             : {
      45      410669 :   long i, lx = lg(x);
      46      410669 :   int first = 1;
      47      416857 :   for (i = 2; i < lx && gequal0(gel(x,i)); i++) /* ~ RgX_valrem + normalize */
      48        6188 :     if (first && !isexactzero(gel(x,i)))
      49             :     {
      50           7 :       pari_warn(warner,"normalizing a series with 0 leading term");
      51           7 :       first = 0;
      52             :     }
      53      410669 :   return RgX_to_ser_i(x, l, lx, i - 2, 0);
      54             : }
      55             : GEN
      56         259 : rfrac_to_ser(GEN x, long l)
      57             : {
      58         259 :   GEN d = gel(x,2);
      59         259 :   if (l == 2)
      60             :   {
      61          14 :     long v = varn(d);
      62          14 :     return zeroser(varn(d), gvaluation(x, pol_x(v)));
      63             :   }
      64         245 :   return gdiv(gel(x,1), RgX_to_ser(d, l));
      65             : }
      66             : 
      67             : static GEN
      68        4711 : RgV_to_ser_i(GEN x, long v, long l, int copy)
      69             : {
      70        4711 :   long j, lx = minss(lg(x), l-1);
      71             :   GEN y;
      72        4711 :   if (lx == 1) return zeroser(v, l-2);
      73        4704 :   y = cgetg(l, t_SER); y[1] = evalvarn(v)|evalvalp(0);
      74        4704 :   x--;
      75        4704 :   if (copy)
      76         231 :     for (j = 2; j <= lx; j++) gel(y,j) = gcopy(gel(x,j));
      77             :   else
      78        4473 :     for (j = 2; j <= lx; j++) gel(y,j) = gel(x,j);
      79        4704 :   for (     ; j < l;   j++) gel(y,j) = gen_0;
      80        4704 :   return normalize(y);
      81             : }
      82             : GEN
      83        4473 : RgV_to_ser(GEN x, long v, long l) { return RgV_to_ser_i(x, v, l, 0); }
      84             : 
      85             : /* x a t_SER, prec >= 0 */
      86             : GEN
      87        2184 : sertoser(GEN x, long prec)
      88             : {
      89        2184 :   long i, lx = lg(x), l;
      90             :   GEN y;
      91        2184 :   if (lx == 2) return zeroser(varn(x), prec);
      92        2170 :   l = prec+2; lx = minss(lx, l);
      93        2170 :   y = cgetg(l,t_SER); y[1] = x[1];
      94        2170 :   for (i = 2; i < lx; i++) gel(y,i) = gel(x,i);
      95        2170 :   for (     ; i < l;  i++) gel(y,i) = gen_0;
      96        2170 :   return y;
      97             : }
      98             : 
      99             : /* R(1/x) = x^v * n/d, val(n) = val(d) = 0 */
     100             : long
     101          84 : rfracrecip(GEN *pn, GEN *pd)
     102             : {
     103          84 :   long v = degpol(*pd);
     104          84 :   if (typ(*pn) == t_POL && varn(*pn) == varn(*pd))
     105             :   {
     106          42 :     v -= degpol(*pn);
     107          42 :     (void)RgX_valrem(*pn, pn); *pn = RgX_recip(*pn);
     108             :   }
     109          84 :   (void)RgX_valrem(*pd, pd); *pd = RgX_recip(*pd);
     110          84 :   return v;
     111             : }
     112             : 
     113             : /* R(1/x) + O(x^N) */
     114             : GEN
     115          49 : rfracrecip_to_ser_absolute(GEN R, long N)
     116             : {
     117          49 :   GEN n = gel(R,1), d = gel(R,2);
     118          49 :   long v = rfracrecip(&n, &d); /* R(1/x) = x^v * n/d, val(n) = val(d) = 0 */
     119          49 :   if (N <= v) return zeroser(varn(d), N);
     120          49 :   R = gdiv(n, RgX_to_ser(d, N-v+2));
     121          49 :   setvalp(R, v); return R;
     122             : }
     123             : 
     124             : /* assume prec >= 0 */
     125             : GEN
     126       20944 : scalarser(GEN x, long v, long prec)
     127             : {
     128             :   long i, l;
     129             :   GEN y;
     130             : 
     131       20944 :   if (gequal0(x))
     132             :   {
     133        1680 :     if (isrationalzero(x)) return zeroser(v, prec);
     134         448 :     if (!isexactzero(x)) prec--;
     135         448 :     y = cgetg(3, t_SER);
     136         448 :     y[1] = evalsigne(0) | _evalvalp(prec) | evalvarn(v);
     137         448 :     gel(y,2) = gcopy(x); return y;
     138             :   }
     139       19264 :   l = prec + 2; y = cgetg(l, t_SER);
     140       19264 :   y[1] = evalsigne(1) | _evalvalp(0) | evalvarn(v);
     141       19264 :   gel(y,2) = gcopy(x); for (i=3; i<l; i++) gel(y,i) = gen_0;
     142       19264 :   return y;
     143             : }
     144             : 
     145             : /* assume x a t_[SER|POL], apply gtoser to all coeffs */
     146             : static GEN
     147           7 : coefstoser(GEN x, long v, long prec)
     148             : {
     149             :   long i, lx;
     150           7 :   GEN y = cgetg_copy(x, &lx); y[1] = x[1];
     151           7 :   for (i=2; i<lx; i++) gel(y,i) = gtoser(gel(x,i), v, prec);
     152           7 :   return y;
     153             : }
     154             : 
     155             : static void
     156          14 : err_ser_priority(GEN x, long v) { pari_err_PRIORITY("Ser", x, "<", v); }
     157             : /* x a t_POL */
     158             : static GEN
     159          56 : poltoser(GEN x, long v, long prec)
     160             : {
     161          56 :   long s = varncmp(varn(x), v);
     162          56 :   if (s < 0) err_ser_priority(x,v);
     163          49 :   if (s > 0) return scalarser(x, v, prec);
     164          35 :   return RgX_to_ser_i(x, prec+2, lg(x), RgX_val(x), 1);
     165             : }
     166             : /* x a t_RFRAC */
     167             : static GEN
     168          77 : rfractoser(GEN x, long v, long prec)
     169             : {
     170          77 :   long s = varncmp(varn(gel(x,2)), v);
     171             :   pari_sp av;
     172          77 :   if (s < 0) err_ser_priority(x,v);
     173          70 :   if (s > 0) return scalarser(x, v, prec);
     174          35 :   av = avma; return gerepileupto(av, rfrac_to_ser(x, prec+2));
     175             : }
     176             : GEN
     177      980916 : toser_i(GEN x)
     178             : {
     179      980916 :   switch(typ(x))
     180             :   {
     181      102732 :     case t_SER: return x;
     182        1400 :     case t_POL: return RgX_to_ser(x, precdl+2);
     183         140 :     case t_RFRAC: return rfrac_to_ser(x, precdl+2);
     184             :   }
     185      876644 :   return NULL;
     186             : }
     187             : 
     188             : /* conversion: prec ignored if t_VEC or t_SER in variable v */
     189             : GEN
     190         392 : gtoser(GEN x, long v, long prec)
     191             : {
     192         392 :   long tx = typ(x);
     193             : 
     194         392 :   if (v < 0) v = 0;
     195         392 :   if (prec < 0) pari_err_DOMAIN("Ser", "precision", "<", gen_0, stoi(prec));
     196         392 :   if (tx == t_SER)
     197             :   {
     198          35 :     long s = varncmp(varn(x), v);
     199          35 :     if      (s < 0) return coefstoser(x, v, prec);
     200          28 :     else if (s > 0) return scalarser(x, v, prec);
     201          14 :     return gcopy(x);
     202             :   }
     203         357 :   if (is_scalar_t(tx)) return scalarser(x,v,prec);
     204         315 :   switch(tx)
     205             :   {
     206          56 :     case t_POL: return poltoser(x, v, prec);
     207          77 :     case t_RFRAC: return rfractoser(x, v, prec);
     208             :     case t_QFR:
     209          42 :     case t_QFI: return RgV_to_ser_i(x, v, 4+1, 1);
     210          14 :     case t_VECSMALL: x = zv_to_ZV(x);/*fall through*/
     211             :     case t_VEC: case t_COL:
     212         133 :       if (varncmp(gvar(x), v) <= 0) pari_err_PRIORITY("Ser", x, "<=", v);
     213         126 :       return RgV_to_ser_i(x, v, lg(x)+1, 1);
     214             :   }
     215           7 :   pari_err_TYPE("Ser",x);
     216             :   return NULL; /* LCOV_EXCL_LINE */
     217             : }
     218             : /* impose prec */
     219             : GEN
     220         175 : gtoser_prec(GEN x, long v, long prec)
     221             : {
     222         175 :   pari_sp av = avma;
     223         175 :   if (v < 0) v = 0;
     224         175 :   if (prec < 0) pari_err_DOMAIN("Ser", "precision", "<", gen_0, stoi(prec));
     225         168 :   switch(typ(x))
     226             :   {
     227          28 :     case t_SER: if (varn(x) != v) break;
     228          21 :                 return gerepilecopy(av, sertoser(x, prec));
     229             :     case t_QFR:
     230             :     case t_QFI:
     231          28 :       x = RgV_to_ser_i(mkvec3(gel(x,1),gel(x,2),gel(x,3)), v, prec+2, 1);
     232          28 :       return gerepileupto(av, x);
     233          14 :     case t_VECSMALL: x = zv_to_ZV(x);/*fall through*/
     234             :     case t_VEC: case t_COL:
     235          42 :       if (varncmp(gvar(x), v) <= 0) pari_err_PRIORITY("Ser", x, "<=", v);
     236          42 :       return RgV_to_ser_i(x, v, prec+2, 1);
     237             :   }
     238          77 :   return gtoser(x, v, prec);
     239             : }
     240             : GEN
     241         469 : Ser0(GEN x, long v, GEN d, long prec)
     242             : {
     243         469 :   if (!d) return gtoser(x, v, prec);
     244         175 :   if (typ(d) != t_INT)
     245             :   {
     246           7 :     d = gceil(d);
     247           7 :     if (typ(d) != t_INT) pari_err_TYPE("Ser [precision]",d);
     248             :   }
     249         175 :   return gtoser_prec(x, v, itos(d));
     250             : }

Generated by: LCOV version 1.13