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-bordeaux1.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 - modules - DedekZeta.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 16912-212c0f0) Lines: 353 363 97.2 %
Date: 2014-10-20 Functions: 12 14 85.7 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 189 226 83.6 %

           Branch data     Line data    Source code
       1                 :            : /* Copyright (C) 2000  The PARI group.
       2                 :            : 
       3                 :            : This file is part of the PARI/GP package.
       4                 :            : 
       5                 :            : PARI/GP is free software; you can redistribute it and/or modify it under the
       6                 :            : terms of the GNU General Public License as published by the Free Software
       7                 :            : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8                 :            : ANY WARRANTY WHATSOEVER.
       9                 :            : 
      10                 :            : Check the License for details. You should have received a copy of it, along
      11                 :            : with the package; see the file 'COPYING'. If not, write to the Free Software
      12                 :            : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13                 :            : /*******************************************************************/
      14                 :            : /*                                                                 */
      15                 :            : /*                     DEDEKIND ZETA FUNCTION                      */
      16                 :            : /*                                                                 */
      17                 :            : /*******************************************************************/
      18                 :            : #include "pari.h"
      19                 :            : #include "paripriv.h"
      20                 :            : static GEN
      21                 :         56 : dirzetak0(GEN nf, ulong N)
      22                 :            : {
      23                 :         56 :   GEN vect, c, c2, T = nf_get_pol(nf), index = nf_get_index(nf);
      24                 :         56 :   pari_sp av = avma, av2;
      25                 :         56 :   const ulong SQRTN = (ulong)(sqrt(N) + 1e-3);
      26                 :            :   ulong i, p, lx;
      27                 :         56 :   long court[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
      28                 :            :   forprime_t S;
      29                 :            : 
      30                 :         56 :   (void)evallg(N+1);
      31                 :         56 :   c  = cgetalloc(t_VECSMALL, N+1);
      32                 :         56 :   c2 = cgetalloc(t_VECSMALL, N+1);
      33         [ +  + ]:     159264 :   c2[1] = c[1] = 1; for (i=2; i<=N; i++) c[i] = 0;
      34                 :         56 :   u_forprime_init(&S, 2, N);
      35                 :         56 :   av2 = avma;
      36         [ +  + ]:      18536 :   while ( (p = u_forprime_next(&S)) )
      37                 :            :   {
      38                 :      18480 :     avma = av2;
      39         [ +  + ]:      18480 :     if (umodiu(index, p)) /* p does not divide index */
      40                 :            :     {
      41                 :      18466 :       vect = gel(Flx_degfact(ZX_to_Flx(T,p), p),1);
      42                 :      18466 :       lx = lg(vect);
      43                 :            :     }
      44                 :            :     else
      45                 :            :     {
      46                 :            :       GEN P;
      47                 :         14 :       court[2] = p; P = idealprimedec(nf,court);
      48                 :         14 :       lx = lg(P); vect = cgetg(lx,t_VECSMALL);
      49         [ +  + ]:         35 :       for (i=1; i<lx; i++) vect[i] = pr_get_f(gel(P,i));
      50                 :            :     }
      51         [ +  + ]:      18480 :     if (p <= SQRTN)
      52         [ +  + ]:        987 :       for (i=1; i<lx; i++)
      53                 :            :       {
      54                 :        735 :         ulong qn, q = upowuu(p, vect[i]); /* Norm P[i] */
      55 [ +  + ][ +  + ]:        735 :         if (!q || q > N) break;
      56                 :        567 :         memcpy(c2 + 2, c + 2, (N-1)*sizeof(long));
      57                 :            :         /* c2[i] <- c[i] + sum_{k = 1}^{v_q(i)} c[i/q^k] for all i <= N */
      58         [ +  - ]:       1295 :         for (qn = q; qn <= N; qn *= q)
      59                 :            :         {
      60                 :       1295 :           ulong k0 = N/qn, k, k2; /* k2 = k*qn */
      61         [ +  + ]:     376012 :           for (k = k0, k2 = k*qn; k > 0; k--, k2 -=qn) c2[k2] += c[k];
      62         [ +  + ]:       1295 :           if (q > k0) break; /* <=> q*qn > N */
      63                 :            :         }
      64                 :        567 :         swap(c, c2);
      65                 :            :       }
      66                 :            :     else /* p > sqrt(N): simpler */
      67         [ +  + ]:      35623 :       for (i=1; i<lx; i++)
      68                 :            :       {
      69                 :            :         ulong k, k2; /* k2 = k*p */
      70         [ +  + ]:      35308 :         if (vect[i] > 1) break;
      71                 :            :         /* c2[i] <- c[i] + sum_{k = 1}^{v_q(i)} c[i/q^k] for all i <= N */
      72         [ +  + ]:     116907 :         for (k = N/p, k2 = k*p; k > 0; k--, k2 -= p) c[k2] += c[k];
      73                 :            :       }
      74                 :            :   }
      75                 :         56 :   avma = av;
      76                 :         56 :   pari_free(c2); return c;
      77                 :            : }
      78                 :            : 
      79                 :            : GEN
      80                 :          7 : dirzetak(GEN nf, GEN b)
      81                 :            : {
      82                 :            :   GEN z, c;
      83                 :            :   long n;
      84                 :            : 
      85         [ -  + ]:          7 :   if (typ(b) != t_INT) pari_err_TYPE("dirzetak",b);
      86         [ -  + ]:          7 :   if (signe(b) <= 0) return cgetg(1,t_VEC);
      87                 :          7 :   nf = checknf(nf);
      88         [ -  + ]:          7 :   n = itou_or_0(b); if (!n) pari_err_OVERFLOW("dirzetak");
      89                 :          7 :   c = dirzetak0(nf, n);
      90                 :          7 :   z = vecsmall_to_vec(c); pari_free(c); return z;
      91                 :            : }
      92                 :            : 
      93                 :            : /* return a t_REAL */
      94                 :            : GEN
      95                 :        539 : zeta_get_limx(long r1, long r2, long bit)
      96                 :            : {
      97                 :        539 :   pari_sp av = avma;
      98                 :            :   GEN p1, p2, c0, c1, A0;
      99                 :        539 :   long r = r1 + r2, N = r + r2;
     100                 :            : 
     101                 :            :   /* c1 = N 2^(-2r2 / N) */
     102                 :        539 :   c1 = mulrs(powrfrac(real2n(1, DEFAULTPREC), -2*r2, N), N);
     103                 :            : 
     104                 :        539 :   p1 = powru(Pi2n(1, DEFAULTPREC), r - 1);
     105                 :        539 :   p2 = mulir(powuu(N,r), p1); shiftr_inplace(p2, -r2);
     106                 :        539 :   c0 = sqrtr( divrr(p2, powru(c1, r+1)) );
     107                 :            : 
     108                 :        539 :   A0 = logr_abs( gmul2n(c0, bit) ); p2 = divrr(A0, c1);
     109                 :        539 :   p1 = divrr(mulur(N*(r+1), logr_abs(p2)), addsr(2*(r+1), gmul2n(A0,2)));
     110                 :        539 :   return gerepileuptoleaf(av, divrr(addrs(p1, 1), powruhalf(p2, N)));
     111                 :            : }
     112                 :            : /* N_0 = floor( C_K / limx ). Large */
     113                 :            : long
     114                 :        651 : zeta_get_N0(GEN C,  GEN limx)
     115                 :            : {
     116                 :            :   long e;
     117                 :        651 :   pari_sp av = avma;
     118                 :        651 :   GEN z = gcvtoi(gdiv(C, limx), &e); /* avoid truncation error */
     119 [ +  - ][ -  + ]:        651 :   if (e >= 0 || is_bigint(z))
     120                 :          0 :     pari_err_OVERFLOW("zeta_get_N0 [need too many primes]");
     121         [ -  + ]:        651 :   if (DEBUGLEVEL>1) err_printf("\ninitzeta: N0 = %Ps\n", z);
     122                 :        651 :   avma = av; return itos(z);
     123                 :            : }
     124                 :            : 
     125                 :            : /* even i such that limx^i ( (i\2)! )^r1 ( i! )^r2 ~ B. Small. */
     126                 :            : static long
     127                 :        224 : get_i0(long r1, long r2, GEN B, GEN limx)
     128                 :            : {
     129                 :        224 :   long imin = 1, imax = 1400;
     130         [ +  + ]:       2240 :   while(imax - imin >= 4)
     131                 :            :   {
     132                 :       2016 :     long i = (imax + imin) >> 1;
     133                 :       2016 :     GEN t = powru(limx, i);
     134         [ +  + ]:       2016 :     if (!r1)      t = mulrr(t, powru(mpfactr(i  , DEFAULTPREC), r2));
     135         [ +  + ]:       1953 :     else if (!r2) t = mulrr(t, powru(mpfactr(i/2, DEFAULTPREC), r1));
     136                 :            :     else {
     137                 :        126 :       GEN u1 = mpfactr(i/2, DEFAULTPREC);
     138                 :        126 :       GEN u2 = mpfactr(i,   DEFAULTPREC);
     139         [ +  + ]:        126 :       if (r1 == r2) t = mulrr(t, powru(mulrr(u1,u2), r1));
     140                 :         63 :       else t = mulrr(t, mulrr(powru(u1,r1), powru(u2,r2)));
     141                 :            :     }
     142         [ +  + ]:       2016 :     if (mpcmp(t, B) >= 0) imax = i; else imin = i;
     143                 :            :   }
     144                 :        224 :   return imax & ~1; /* make it even */
     145                 :            : }
     146                 :            : 
     147                 :            : /* assume limx = zeta_get_limx(r1, r2, bit), a t_REAL */
     148                 :            : long
     149                 :        224 : zeta_get_i0(long r1, long r2, long bit, GEN limx)
     150                 :            : {
     151                 :        224 :   pari_sp av = avma;
     152                 :        224 :   GEN B = gmul(sqrtr( divrr(powrs(mppi(DEFAULTPREC), r2-3), limx) ),
     153                 :            :                gmul2n(powuu(5, r1), bit + r2));
     154                 :        224 :   long i0 = get_i0(r1, r2, B, limx);
     155         [ -  + ]:        224 :   if (DEBUGLEVEL>1) { err_printf("i0 = %ld\n",i0); err_flush(); }
     156                 :        224 :   avma = av; return i0;
     157                 :            : }
     158                 :            : 
     159                 :            : /* sum(j=1, r-k+1, A[j] * B[j]), assumes k <= r and A[j],B[j] are 'mp' */
     160                 :            : INLINE GEN
     161                 :   26201042 : sumprod(GEN A, GEN B, long r, long k)
     162                 :            : {
     163         [ +  + ]:   26201042 :   GEN s = signe(gel(A,1))? mpmul(gel(A,1), gel(B,1)): gen_0;
     164                 :            :   long j;
     165         [ +  + ]:   52236716 :   for (j=2; j<=r-k+1; j++)
     166         [ +  + ]:   26035674 :     if (signe(gel(A,j))) s = mpadd(s, mpmul(gel(A,j), gel(B,j)));
     167                 :   26201042 :   return s;
     168                 :            : }
     169                 :            : 
     170                 :            : GEN
     171                 :         49 : initzeta(GEN T, long prec)
     172                 :            : {
     173                 :            :   GEN znf, nf, bnf, gr2, gru, p1, p2, cst, coef;
     174                 :            :   GEN limx, resi,zet,C,coeflog,racpi,aij,tabj,colzero, tabcstn, tabcstni;
     175                 :            :   GEN c_even, ck_even, c_odd, ck_odd, serie_even, serie_odd, serie_exp;
     176                 :            :   GEN VOID;
     177                 :         49 :   long N0, i0, r1, r2, r, R, N, i, j, k, n, bit = prec2nbits(prec) + 6;
     178                 :            :   pari_timer ti;
     179                 :            : 
     180                 :            :   pari_sp av, av2;
     181                 :         49 :   long court[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
     182                 :            : 
     183                 :            :   /*************** residue & constants ***************/
     184                 :         49 :   T = get_bnfpol(T, &bnf, &nf);
     185         [ +  - ]:         49 :   if (!nf) {
     186                 :         49 :     bnf = Buchall(T, 0, prec);
     187                 :         49 :     nf = bnf_get_nf(bnf);
     188                 :            :   }
     189         [ #  # ]:          0 :   else if (!bnf) {
     190         [ #  # ]:          0 :     if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
     191                 :          0 :     bnf = Buchall(nf, 0, prec);
     192         [ #  # ]:          0 :   } else if (nf_get_prec(bnf) < prec) {
     193                 :          0 :     bnf = bnfnewprec(bnf, prec);
     194                 :          0 :     nf = bnf_get_nf(bnf);
     195                 :            :   }
     196                 :            : 
     197                 :         49 :   prec = precdbl(prec);
     198                 :         49 :   racpi = sqrtr(mppi(prec));
     199                 :            :   /* all the above left on stack for efficiency */
     200                 :            : 
     201                 :            :   /* class number & regulator */
     202                 :         49 :   N = nf_get_degree(nf);
     203                 :         49 :   nf_get_sign(nf, &r1, &r2);
     204                 :         49 :   gr2 = gmael(nf,2,2);
     205                 :         49 :   r = r1 + r2; R = r+2;
     206                 :            : 
     207                 :         49 :   znf = cgetg(10,t_VEC);
     208                 :         49 :   VOID = cgetg(1, t_STR); /* unused value */
     209                 :         49 :   av = avma;
     210                 :         49 :   resi = gerepileupto(av,
     211                 :            :            gdivgs(mpmul(shifti(bnf_get_no(bnf),r1), bnf_get_reg(bnf)),
     212                 :            :                   bnf_get_tuN(bnf))); /* hr 2^r1 / w*/
     213                 :         49 :   av = avma;
     214                 :         49 :   p1 = sqrtr_abs(itor(nf_get_disc(nf), prec));
     215                 :         49 :   p2 = gmul2n(powru(racpi,N), r2);
     216                 :         49 :   cst = gerepileuptoleaf(av, divrr(p1,p2));
     217                 :            : 
     218                 :            :   /* N0, i0 */
     219                 :         49 :   limx = zeta_get_limx(r1, r2, bit);
     220                 :         49 :   N0 = zeta_get_N0(cst, limx);
     221                 :         49 :   i0 = zeta_get_i0(r1, r2, bit + 4, limx);
     222                 :            : 
     223                 :            :   /* Array of i/cst (i=1..N0) */
     224                 :         49 :   av = avma;
     225                 :         49 :   tabcstn  = cgetg(N0+1,t_VEC);
     226                 :         49 :   tabcstni = cgetg(N0+1,t_VEC);
     227                 :         49 :   p1 = invr(cst);
     228         [ +  + ]:     159103 :   for (i=1; i<=N0; i++) gel(tabcstni,i) = gel(tabcstn,i) = mulur(i,p1);
     229                 :         49 :   tabcstn  = gclone(tabcstn);
     230                 :         49 :   tabcstni = gclone(tabcstni);
     231                 :            : 
     232                 :            :   /********** compute a(i,j) **********/
     233         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_start(&ti);
     234                 :         49 :   zet = cgetg(R,t_VEC); gel(zet,1) = mpeuler(prec);
     235         [ +  + ]:        147 :   for (i=2; i<R; i++) gel(zet,i) = szeta(i, prec);
     236                 :            : 
     237                 :         49 :   aij = cgetg(i0+1,t_VEC);
     238         [ +  + ]:      15309 :   for (i=1; i<=i0; i++) gel(aij,i) = cgetg(R,t_VEC);
     239                 :            : 
     240                 :         49 :   c_even = real2n(r1, prec);
     241                 :         49 :   c_odd = gmul(c_even, powru(racpi,r1));
     242         [ +  + ]:         49 :   if (r&1) c_odd = gneg_i(c_odd);
     243                 :         49 :   ck_even = cgetg(R,t_VEC); ck_odd = cgetg(r2+2,t_VEC);
     244         [ +  + ]:        196 :   for (k=1; k<R; k++)
     245                 :            :   {
     246                 :        147 :     GEN t = mulri(gel(zet,k), addis(shifti(gr2, k), r1));
     247                 :        147 :     shiftr_inplace(t, -k);
     248         [ +  + ]:        147 :     if (k&1) togglesign(t);
     249                 :        147 :     gel(ck_even,k) = t;
     250                 :            :   }
     251                 :         49 :   gru = utoipos(r);
     252         [ +  + ]:        126 :   for (k = 1; k <= r2+1; k++)
     253                 :            :   {
     254                 :         77 :     GEN t = mulri(gel(zet,k), subis(shifti(gru,k), r1));
     255                 :         77 :     shiftr_inplace(t, -k);
     256         [ +  + ]:         77 :     if (k&1) togglesign(t);
     257                 :         77 :     gel(ck_odd,k) = addsr(r, t);
     258                 :            :   }
     259         [ +  + ]:         49 :   if (r1) gel(ck_odd,1) = subrr(gel(ck_odd,1), mulur(r1, mplog2(prec)));
     260                 :         49 :   serie_even = cgetg(r+3,t_SER);
     261                 :         49 :   serie_odd = cgetg(r2+3,t_SER);
     262                 :         49 :   serie_even[1] = serie_odd[1] = evalsigne(1)|_evalvalp(1);
     263                 :         49 :   i = 0;
     264         [ +  + ]:       7679 :   while (i < i0/2)
     265                 :            :   {
     266         [ +  + ]:      28329 :     for (k=1; k<R; k++) gel(serie_even,k+1) = gdivgs(gel(ck_even,k),k);
     267                 :       7630 :     serie_exp = gmul(c_even, gexp(serie_even,0));
     268                 :       7630 :     p1 = gel(aij,2*i+1);
     269         [ +  + ]:      28329 :     for (j=1; j<R; j++) p1[j] = serie_exp[r+3-j];
     270                 :            : 
     271         [ +  + ]:      18053 :     for (k=1; k<=r2+1; k++) gel(serie_odd,k+1) = gdivgs(gel(ck_odd,k),k);
     272                 :       7630 :     serie_exp = gmul(c_odd, gexp(serie_odd,0));
     273                 :       7630 :     p1 = gel(aij,2*i+2);
     274         [ +  + ]:      18053 :     for (j=1; j<=r2+1; j++) p1[j] = serie_exp[r2+3-j];
     275         [ +  + ]:      17906 :     for (   ; j<R; j++) gel(p1,j) = gen_0;
     276                 :       7630 :     i++;
     277                 :            : 
     278                 :       7630 :     c_even = gdiv(c_even,gmul(powuu(i,r),powuu(2*i-1,r2)));
     279                 :       7630 :     c_odd  = gdiv(c_odd, gmul(powuu(i,r2),powuu(2*i+1,r)));
     280                 :       7630 :     c_even = gmul2n(c_even,-r2);
     281                 :       7630 :     c_odd  = gmul2n(c_odd,r1-r2);
     282         [ +  + ]:       7630 :     if (r1 & 1) { c_even = gneg_i(c_even); c_odd = gneg_i(c_odd); }
     283                 :       7630 :     p1 = gr2; p2 = gru;
     284         [ +  + ]:      28329 :     for (k=1; k<R; k++)
     285                 :            :     {
     286                 :      20699 :       p1 = gdivgs(p1,2*i-1); p2 = gdivgs(p2,2*i);
     287                 :      20699 :       gel(ck_even,k) = gadd(gel(ck_even,k), gadd(p1,p2));
     288                 :            :     }
     289                 :       7630 :     p1 = gru; p2 = gr2;
     290         [ +  + ]:      18053 :     for (k=1; k<=r2+1; k++)
     291                 :            :     {
     292                 :      10423 :       p1 = gdivgs(p1,2*i+1); p2 = gdivgs(p2,2*i);
     293                 :      10423 :       gel(ck_odd,k) = gadd(gel(ck_odd,k), gadd(p1,p2));
     294                 :            :     }
     295                 :            :   }
     296                 :         49 :   aij = gerepilecopy(av, aij);
     297         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_printf(&ti,"a(i,j)");
     298                 :         49 :   gel(znf,1) = mkvecsmall2(r1,r2);
     299                 :         49 :   gel(znf,2) = resi;
     300                 :         49 :   gel(znf,5) = cst;
     301                 :         49 :   gel(znf,6) = logr_abs(cst);
     302                 :         49 :   gel(znf,7) = aij;
     303                 :            : 
     304                 :            :   /************* Calcul du nombre d'ideaux de norme donnee *************/
     305                 :         49 :   coef = dirzetak0(nf,N0); tabj = cgetg(N0+1,t_MAT);
     306         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_printf(&ti,"coef");
     307                 :         49 :   colzero = zerocol(r+1);
     308         [ +  + ]:     159103 :   for (i=1; i<=N0; i++)
     309         [ +  + ]:     159054 :     if (coef[i])
     310                 :            :     {
     311                 :      65814 :       GEN t = cgetg(r+2,t_COL);
     312                 :      65814 :       p1 = logr_abs(gel(tabcstn,i)); togglesign(p1);
     313                 :      65814 :       gel(t,1) = utoi(coef[i]);
     314                 :      65814 :       gel(t,2) = mulur(coef[i], p1);
     315         [ +  + ]:     195615 :       for (j=2; j<=r; j++)
     316                 :            :       {
     317                 :     129801 :         pari_sp av2 = avma;
     318                 :     129801 :         gel(t,j+1) = gerepileuptoleaf(av2, divru(mulrr(gel(t,j), p1), j));
     319                 :            :       }
     320                 :      65814 :       gel(tabj,i) = t; /* tabj[n,j] = coef(n)*ln(c/n)^(j-1)/(j-1)! */
     321                 :            :     }
     322                 :            :     else
     323                 :      93240 :       gel(tabj,i) = colzero;
     324         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_printf(&ti,"a(n)");
     325                 :            : 
     326                 :         49 :   coeflog=cgetg(N0+1,t_VEC); gel(coeflog,1) = VOID;
     327         [ +  + ]:     159054 :   for (i=2; i<=N0; i++)
     328         [ +  + ]:     159005 :     if (coef[i])
     329                 :            :     {
     330                 :      65765 :       court[2] = i; p1 = glog(court,prec);
     331                 :      65765 :       setsigne(p1, -1); gel(coeflog,i) = p1;
     332                 :            :     }
     333                 :      93240 :     else gel(coeflog,i) = VOID;
     334         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_printf(&ti,"log(n)");
     335                 :            : 
     336                 :         49 :   gel(znf,3) = tabj;
     337                 :         49 :   gel(znf,8) = vecsmall_copy(coef);
     338                 :         49 :   gel(znf,9) = coeflog;
     339                 :            : 
     340                 :            :   /******************** Calcul des coefficients Cik ********************/
     341                 :         49 :   C = cgetg(r+1,t_MAT);
     342         [ +  + ]:        147 :   for (k=1; k<=r; k++) gel(C,k) = cgetg(i0+1,t_COL);
     343                 :         49 :   av2 = avma;
     344         [ +  + ]:      15309 :   for (i=1; i<=i0; i++)
     345                 :            :   {
     346         [ +  + ]:      41398 :     for (k=1; k<=r; k++)
     347                 :            :     {
     348                 :      26138 :       GEN A = gel(aij,i) + k; /* A[j] = aij[i, j+k] */
     349                 :      26138 :       GEN t = sumprod(A, gel(tabj,1), r, k);
     350                 :            :       /* n = 1 */
     351 [ +  + ][ +  + ]:      26138 :       if (i > 1 && signe(t)) t = mpmul(t, gel(tabcstni,1));
     352         [ +  + ]:   62902168 :       for (n=2; n<=N0; n++)
     353         [ +  + ]:   62876030 :         if (coef[n])
     354                 :            :         { /* sum(j=1, r-k+1, aij[i,j+k] * tabj[n, j]) */
     355                 :   26043584 :           GEN s = sumprod(A, gel(tabj,n), r, k);
     356         [ +  + ]:   26043584 :           if (!signe(s)) continue;
     357         [ +  + ]:   21369502 :           if (i > 1) s = mpmul(s, gel(tabcstni,n));
     358                 :   21369502 :           t = mpadd(t,s);
     359                 :            :         }
     360                 :      26138 :       togglesign(t);
     361                 :      26138 :       gcoeff(C,i,k) = gerepileuptoleaf(av2,t);
     362                 :      26138 :       av2 = avma;
     363                 :            :     }
     364 [ +  + ][ +  + ]:      15260 :     if (i > 1 && i < i0) {
     365         [ +  + ]:   20929566 :       for (n=1; n<=N0; n++)
     366         [ +  + ]:   20914404 :         if (coef[n]) mpmulz(gel(tabcstni,n), gel(tabcstn,n), gel(tabcstni,n));
     367                 :            :     }
     368                 :            :   }
     369                 :         49 :   gel(znf,4) = C;
     370         [ -  + ]:         49 :   if (DEBUGLEVEL>1) timer_printf(&ti,"Cik");
     371                 :         49 :   gunclone(tabcstn); gunclone(tabcstni);
     372                 :         49 :   pari_free((void*)coef); return znf;
     373                 :            : }
     374                 :            : 
     375                 :            : static void
     376                 :        280 : znf_get_sign(GEN znf, long *r1, long *r2)
     377                 :        280 : { GEN v = gel(znf,1); *r1 = v[1]; *r2 = v[2]; }
     378                 :            : 
     379                 :            : /* s != 0,1 */
     380                 :            : static GEN
     381                 :         91 : slambdak(GEN znf, long s, long flag, long prec)
     382                 :            : {
     383                 :            :   GEN resi, C, cst, cstlog, coeflog, cs, coef;
     384                 :            :   GEN lambd, gammas2, gammaunmoins2, var1, var2;
     385                 :            :   GEN gar, val, valm, valk, valkm;
     386                 :            :   long r1, r2, r, i0, i, k, N0;
     387                 :            : 
     388                 :         91 :   znf_get_sign(znf, &r1, &r2); r = r1+r2;
     389                 :         91 :   resi   = gel(znf,2);
     390                 :         91 :   C      = gel(znf,4);
     391                 :         91 :   cst    = gel(znf,5);
     392                 :         91 :   cstlog = gel(znf,6);
     393                 :         91 :   coef   = gel(znf,8);
     394                 :         91 :   coeflog= gel(znf,9);
     395                 :         91 :   i0 = nbrows(C);
     396                 :         91 :   N0 = lg(coef)-1;
     397                 :            : 
     398 [ +  + ][ +  - ]:         91 :   if (s < 0 && (r2 || !odd(s))) s = 1 - s;
                 [ -  + ]
     399                 :            : 
     400                 :            :   /* s >= 2 or s < 0 */
     401                 :         91 :   lambd = gdiv(resi, mulss(s, s-1));
     402                 :         91 :   gammas2 = ggamma(gmul2n(stoi(s),-1),prec);
     403                 :         91 :   gar = gpowgs(gammas2,r1);
     404                 :         91 :   cs = mpexp( mulrs(cstlog,s) );
     405                 :         91 :   val = stoi(s); valm = stoi(1 - s);
     406         [ +  + ]:         91 :   if (s < 0) /* r2 = 0 && odd(s) */
     407                 :            :   {
     408                 :         21 :     gammaunmoins2 = ggamma(gmul2n(valm,-1),prec); /* gamma((1-s) / 2) */
     409                 :         21 :     var1 = var2 = gen_1;
     410         [ +  + ]:       3031 :     for (i=2; i<=N0; i++)
     411         [ +  + ]:       3010 :       if (coef[i])
     412                 :            :       {
     413                 :       1309 :         GEN gexpro = mpexp(mulrs(gel(coeflog,i),s));
     414                 :       1309 :         var1 = gadd(var1, mulsr(coef[i],gexpro));
     415                 :       1309 :         var2 = gadd(var2, divsr(coef[i],mulsr(i,gexpro)));
     416                 :            :       }
     417                 :         21 :     lambd = gadd(lambd,gmul(gmul(var1,cs),gar));
     418                 :         21 :     lambd = gadd(lambd,gmul(gmul(var2,gdiv(cst,cs)),
     419                 :            :                             gpowgs(gammaunmoins2,r1)));
     420                 :         21 :     var1 = gen_0;
     421         [ +  + ]:       3038 :     for (i=1; i<=i0; i+=2)
     422                 :            :     {
     423                 :       3017 :       valk  = val;
     424                 :       3017 :       valkm = valm;
     425         [ +  + ]:       9814 :       for (k=1; k<=r; k++)
     426                 :            :       {
     427                 :       6797 :         GEN c = gcoeff(C,i,k);
     428                 :       6797 :         var1 = mpadd(var1,mpdiv(c,valk )); valk  = mulii(val,valk);
     429                 :       6797 :         var1 = mpadd(var1,mpdiv(c,valkm)); valkm = mulii(valm,valkm);
     430                 :            :       }
     431                 :       3017 :       val  = addis(val, 2);
     432                 :       3017 :       valm = addis(valm,2);
     433                 :            :     }
     434                 :            :   }
     435                 :            :   else
     436                 :            :   {
     437                 :         70 :     GEN tabj = gel(znf,3), aij = gel(znf,7), A = gel(aij,s);
     438                 :            :     long n;
     439                 :            : 
     440                 :         70 :     gar = gmul(gar,gpowgs(mpfactr(s-1,prec),r2)); /* x gamma(s)^r2 */
     441                 :            :     /* n = 1 */
     442                 :         70 :     var1 = gen_1;
     443         [ +  - ]:         70 :     var2 = (s <= i0)? sumprod(A, gel(tabj,1), r, 0): gen_0;
     444         [ +  + ]:     317492 :     for (n=2; n<=N0; n++)
     445         [ +  + ]:     317422 :       if (coef[n])
     446                 :            :       {
     447                 :     131250 :         GEN gexpro = mpexp( mulrs(gel(coeflog,n),s) );
     448                 :     131250 :         var1 = mpadd(var1, mulsr(coef[n],gexpro));
     449         [ +  - ]:     131250 :         if (s <= i0)
     450                 :            :         {
     451                 :     131250 :           GEN t = sumprod(A, gel(tabj,n), r, 0);
     452         [ -  + ]:     131250 :           if (!signe(t)) continue;
     453                 :     131250 :           var2 = mpadd(var2, mpdiv(t, mulsr(n,gexpro)));
     454                 :            :         }
     455                 :            :       }
     456                 :         70 :     lambd = gadd(lambd,gmul(gmul(var1,cs),gar));
     457                 :         70 :     lambd = gadd(lambd,gmul(var2,gdiv(cst,cs)));
     458                 :         70 :     var1 = gen_0;
     459         [ +  + ]:      16982 :     for (n=1; n<=i0; n++)
     460                 :            :     {
     461                 :      16912 :       valk  = val;
     462                 :      16912 :       valkm = valm;
     463         [ +  + ]:      16912 :       if (n == s)
     464         [ +  + ]:        224 :         for (k=1; k<=r; k++)
     465                 :            :         {
     466                 :        154 :           GEN c = gcoeff(C,n,k);
     467                 :        154 :           var1 = mpadd(var1,mpdiv(c,valk)); valk = mulii(val,valk);
     468                 :            :         }
     469                 :            :       else
     470         [ +  + ]:      50848 :       for (k=1; k<=r; k++)
     471                 :            :       {
     472                 :      34006 :           GEN c = gcoeff(C,n,k);
     473                 :      34006 :           var1 = mpadd(var1,mpdiv(c,valk )); valk  = mulii(val,valk);
     474                 :      34006 :           var1 = mpadd(var1,mpdiv(c,valkm)); valkm = mulii(valm,valkm);
     475                 :            :       }
     476                 :      16912 :       val  = addis(val,1);
     477                 :      16912 :       valm = addis(valm,1);
     478                 :            :     }
     479                 :            :   }
     480                 :         91 :   lambd = gadd(lambd, var1);
     481         [ +  - ]:         91 :   if (!flag) lambd = gdiv(lambd,gmul(gar,cs)); /* zetak */
     482                 :         91 :   return lambd;
     483                 :            : }
     484                 :            : 
     485                 :            : /* s not an integer */
     486                 :            : static GEN
     487                 :         42 : cxlambdak(GEN znf, GEN s, long flag, long prec)
     488                 :            : {
     489                 :            :   GEN resi, C, cst, cstlog, coeflog, cs, coef;
     490                 :            :   GEN lambd, gammas, gammaunmoins, gammas2, gammaunmoins2, var1, var2;
     491                 :            :   GEN smoinun, unmoins, gar, val, valm, valk, valkm, Pi, sPi;
     492                 :            :   long r1, r2, r, i0, i, k, N0, bigprec;
     493                 :            : 
     494                 :         42 :   znf_get_sign(znf, &r1, &r2);
     495                 :         42 :   resi   = gel(znf,2);
     496                 :         42 :   C      = gel(znf,4);
     497                 :         42 :   cst    = gel(znf,5);
     498                 :         42 :   cstlog = gel(znf,6);
     499                 :         42 :   coef   = gel(znf,8);
     500                 :         42 :   coeflog= gel(znf,9);
     501                 :         42 :   r1 = mael(znf,1,1);
     502                 :         42 :   r2 = mael(znf,1,2); r = r1+r2;
     503                 :         42 :   i0 = nbrows(C);
     504                 :         42 :   N0 = lg(coef)-1;
     505                 :         42 :   bigprec = precision(cst);
     506                 :            : 
     507                 :         42 :   Pi = mppi(bigprec);
     508                 :         42 :   s = gtofp(s, bigprec); sPi = gmul(s, Pi);
     509                 :         42 :   smoinun = gsubgs(s,1);
     510                 :         42 :   unmoins = gneg(smoinun);
     511                 :         42 :   lambd = gdiv(resi,gmul(s,smoinun));
     512                 :         42 :   gammas = ggamma(s,prec);
     513                 :         42 :   gammas2= ggamma(gmul2n(s,-1),prec);
     514                 :         42 :   gar = gmul(gpowgs(gammas,r2),gpowgs(gammas2,r1));
     515                 :         42 :   cs = gexp(gmul(cstlog,s),prec);
     516                 :         42 :   gammaunmoins = gdiv(Pi, gmul(gsin(sPi,prec),gammas));
     517                 :         42 :   gammaunmoins2= gdiv(gmul(gmul(sqrtr(Pi),gpow(gen_2,smoinun,prec)),
     518                 :            :                            gammas2),
     519                 :            :                       gmul(gcos(gmul2n(sPi,-1),prec),gammas));
     520                 :         42 :   var1 = var2 = gen_1;
     521         [ +  + ]:     159019 :   for (i=2; i<=N0; i++)
     522         [ +  + ]:     158977 :     if (coef[i])
     523                 :            :     {
     524                 :      65737 :       GEN gexpro = gexp(gmul(gel(coeflog,i),s),bigprec);
     525                 :      65737 :       var1 = gadd(var1,gmulsg(coef[i], gexpro));
     526                 :      65737 :       var2 = gadd(var2,gdivsg(coef[i], gmulsg(i,gexpro)));
     527                 :            :     }
     528                 :         42 :   lambd = gadd(lambd,gmul(gmul(var1,cs),gar));
     529                 :         42 :   lambd = gadd(lambd,gmul(gmul(gmul(var2,gdiv(cst,cs)),
     530                 :            :                                gpowgs(gammaunmoins,r2)),
     531                 :            :                           gpowgs(gammaunmoins2,r1)));
     532                 :         42 :   val  = s;
     533                 :         42 :   valm = unmoins;
     534                 :         42 :   var1 = gen_0;
     535         [ +  + ]:       7735 :   for (i=1; i<=i0; i++)
     536                 :            :   {
     537                 :       7693 :     valk  = val;
     538                 :       7693 :     valkm = valm;
     539         [ +  + ]:      22484 :     for (k=1; k<=r; k++)
     540                 :            :     {
     541                 :      14791 :       GEN c = gcoeff(C,i,k);
     542                 :      14791 :       var1 = gadd(var1,gdiv(c,valk )); valk  = gmul(val, valk);
     543                 :      14791 :       var1 = gadd(var1,gdiv(c,valkm)); valkm = gmul(valm,valkm);
     544                 :            :     }
     545         [ +  + ]:       7693 :     if (r2)
     546                 :            :     {
     547                 :       4676 :       val  = gaddgs(val, 1);
     548                 :       4676 :       valm = gaddgs(valm,1);
     549                 :            :     }
     550                 :            :     else
     551                 :            :     {
     552                 :       3017 :       val  = gaddgs(val, 2);
     553                 :       3017 :       valm = gaddgs(valm,2); i++;
     554                 :            :     }
     555                 :            :   }
     556                 :         42 :   lambd = gadd(lambd, var1);
     557         [ +  - ]:         42 :   if (!flag) lambd = gdiv(lambd,gmul(gar,cs)); /* zetak */
     558                 :         42 :   return lambd;
     559                 :            : }
     560                 :            : 
     561                 :            : GEN
     562                 :        189 : gzetakall(GEN znf, GEN s, long flag, long prec)
     563                 :            : {
     564                 :        189 :   pari_sp av = avma;
     565                 :            :   GEN z;
     566                 :            : 
     567 [ +  - ][ +  - ]:        189 :   if (typ(znf)!=t_VEC || lg(znf)!=10 || typ(gel(znf,1)) != t_VECSMALL)
                 [ -  + ]
     568                 :          0 :     pari_err_TYPE("zetakall", znf);
     569         [ +  + ]:        189 :   if (isint(s, &s))
     570                 :            :   {
     571                 :        147 :     long ss = itos(s), r1, r2;
     572         [ -  + ]:        147 :     if (ss==1) pari_err_DOMAIN("zetak", "argument", "=", gen_1, s);
     573                 :        147 :     znf_get_sign(znf, &r1, &r2);
     574         [ +  + ]:        147 :     if (ss==0)
     575                 :            :     {
     576                 :         35 :       avma = av;
     577         [ -  + ]:         35 :       if (flag) pari_err_DOMAIN("zetak", "argument", "=", gen_0, s);
     578         [ +  + ]:         35 :       if (r1 + r2 > 1) return gen_0;
     579         [ -  + ]:          7 :       return r1? mkfrac(gen_m1, gen_2): gneg(gel(znf, 2));
     580                 :            :     }
     581 [ +  - ][ +  + ]:        112 :     if (!flag && ss < 0 && (r2 || !odd(ss))) return gen_0;
         [ +  + ][ -  + ]
     582                 :        147 :     z = slambdak(znf, itos(s), flag, prec+EXTRAPRECWORD);
     583                 :            :   }
     584                 :            :   else
     585                 :         42 :     z = cxlambdak(znf, s, flag, prec+EXTRAPRECWORD);
     586         [ +  - ]:        133 :   if (gprecision(z) > prec) z = gprec_w(z, prec);
     587                 :        189 :   return gerepileupto(av, z);
     588                 :            : }
     589                 :            : GEN
     590                 :          0 : gzetak(GEN nfz, GEN s, long prec) { return gzetakall(nfz,s,0,prec); }
     591                 :            : GEN
     592                 :          0 : glambdak(GEN nfz, GEN s, long prec) { return gzetakall(nfz,s,1,prec); }

Generated by: LCOV version 1.9