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 - dirichlet.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.18.1 lcov report (development 30504-e6f76a2795) Lines: 483 506 95.5 %
Date: 2025-09-29 09:23:05 Functions: 42 45 93.3 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2015  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; either version 2 of the License, or (at your option) any later
       8             : version. It is distributed in the hope that it will be useful, but WITHOUT
       9             : ANY WARRANTY WHATSOEVER.
      10             : 
      11             : Check the License for details. You should have received a copy of it, along
      12             : with the package; see the file 'COPYING'. If not, write to the Free Software
      13             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      14             : 
      15             : /********************************************************************/
      16             : /**                                                                **/
      17             : /**           Dirichlet series through Euler product               **/
      18             : /**                                                                **/
      19             : /********************************************************************/
      20             : #include "pari.h"
      21             : #include "paripriv.h"
      22             : 
      23             : static void
      24          28 : err_direuler(GEN x)
      25          28 : { pari_err_DOMAIN("direuler","constant term","!=", gen_1,x); }
      26             : 
      27             : /* s = t_POL (tolerate t_SER of valuation 0) of constant term = 1
      28             :  * d = minimal such that p^d > X
      29             :  * V indexed by 1..X will contain the a_n
      30             :  * v[1..n] contains the indices nj such that V[nj] != 0 */
      31             : static long
      32       28728 : dirmuleuler_small(GEN V, GEN v, long n, ulong p, GEN s, long d)
      33             : {
      34       28728 :   long i, j, m = n, D = minss(d+2, lg(s));
      35       28728 :   ulong q = 1, X = lg(V)-1;
      36             : 
      37       94780 :   for (i = 3, q = p; i < D; i++, q *= p) /* q*p does not overflow */
      38             :   {
      39       66052 :     GEN aq = gel(s,i);
      40       66052 :     if (gequal0(aq)) continue;
      41             :     /* j = 1 */
      42       53767 :     gel(V,q) = aq;
      43       53767 :     v[++n] = q;
      44     3268034 :     for (j = 2; j <= m; j++)
      45             :     {
      46     3214267 :       ulong nj = umuluu_le(uel(v,j), q, X);
      47     3214267 :       if (!nj) continue;
      48      192017 :       gel(V,nj) = gmul(aq, gel(V,v[j]));
      49      192017 :       v[++n] = nj;
      50             :     }
      51             :   }
      52       28728 :   return n;
      53             : }
      54             : 
      55             : /* ap != 0 for efficiency, p > sqrt(X) */
      56             : static void
      57      308861 : dirmuleuler_large(GEN V, ulong p, GEN ap)
      58             : {
      59      308861 :   long j, jp, X = lg(V)-1;
      60      308861 :   gel(V,p) = ap;
      61     1506694 :   for (j = 2, jp = 2*p; jp <= X; j++, jp += p) gel(V,jp) = gmul(ap, gel(V,j));
      62      308861 : }
      63             : 
      64             : static ulong
      65       10283 : direulertou(GEN a, GEN fl(GEN))
      66             : {
      67       10283 :   if (typ(a) != t_INT)
      68             :   {
      69          49 :     a = fl(a);
      70          28 :     if (typ(a) != t_INT) pari_err_TYPE("direuler", a);
      71             :   }
      72       10262 :   return signe(a)<=0 ? 0: itou(a);
      73             : }
      74             : 
      75             : static GEN
      76        3731 : direuler_Sbad(GEN V, GEN v, GEN Sbad, ulong *n)
      77             : {
      78        3731 :   long i, l = lg(Sbad);
      79        3731 :   ulong X = lg(V)-1;
      80        3731 :   GEN pbad = gen_1;
      81        9674 :   for (i = 1; i < l; i++)
      82             :   {
      83        5978 :     GEN ai = gel(Sbad,i), gq;
      84             :     ulong q;
      85        5978 :     if (typ(ai) != t_VEC || lg(ai) != 3)
      86          14 :       pari_err_TYPE("direuler [bad primes]",ai);
      87        5964 :     gq = gel(ai,1);
      88        5964 :     if (typ(gq) != t_INT || signe(gq) <= 0)
      89           7 :       pari_err_TYPE("direuler [prime expected]",gq);
      90        5957 :     if (lgefint(gq) == 3 && (q = uel(gq,2)) <= X)
      91             :     {
      92        4823 :       long d = ulogint(X, q) + 1;
      93        4823 :       GEN s = direuler_factor(gel(ai,2), d);
      94        4809 :       *n = dirmuleuler_small(V, v, *n, q, s, d);
      95        4809 :       pbad = muliu(pbad, q);
      96             :     }
      97             :   }
      98        3696 :   return pbad;
      99             : }
     100             : 
     101             : GEN
     102         672 : direuler_bad(void *E, GEN (*eval)(void *,GEN,long), GEN a,GEN b,GEN c, GEN Sbad)
     103             : {
     104             :   ulong au, bu, X, sqrtX, n, p;
     105         672 :   pari_sp av0 = avma;
     106             :   GEN gp, v, V;
     107             :   forprime_t T;
     108         672 :   au = direulertou(a, gceil);
     109         665 :   bu = direulertou(b, gfloor);
     110         658 :   X = c ? direulertou(c, gfloor): bu;
     111         651 :   if (X == 0) return cgetg(1,t_VEC);
     112         644 :   if (bu > X) bu = X;
     113         644 :   if (!u_forprime_init(&T, au, bu)) { set_avma(av0); return mkvec(gen_1); }
     114         630 :   v = vecsmall_ei(X, 1);
     115         630 :   V = vec_ei(X, 1);
     116         630 :   n = 1;
     117         630 :   if (Sbad) Sbad = direuler_Sbad(V, v, Sbad, &n);
     118         595 :   p = 1; gp = cgetipos(3); sqrtX = usqrt(X);
     119        8316 :   while (p <= sqrtX && (p = u_forprime_next(&T)))
     120        7742 :     if (!Sbad || umodiu(Sbad, p))
     121             :     {
     122        7637 :       long d = ulogint(X, p) + 1; /* minimal d such that p^d > X */
     123             :       GEN s;
     124        7637 :       gp[2] = p; s = eval(E, gp, d);
     125        7616 :       n = dirmuleuler_small(V, v, n, p, s, d);
     126             :     }
     127      740201 :   while ((p = u_forprime_next(&T))) /* sqrt(X) < p <= X */
     128      739627 :     if (!Sbad || umodiu(Sbad, p))
     129             :     {
     130             :       GEN s;
     131      739620 :       gp[2] = p; s = eval(E, gp, 2); /* s either t_POL or t_SER of val 0 */
     132      739620 :       if (lg(s) > 3 && !gequal0(gel(s,3)))
     133      139153 :         dirmuleuler_large(V, p, gel(s,3));
     134             :     }
     135         574 :   return gc_GEN(av0,V);
     136             : }
     137             : 
     138             : /* return a t_SER or a truncated t_POL to precision n */
     139             : GEN
     140      752080 : direuler_factor(GEN s, long n)
     141             : {
     142      752080 :   long t = typ(s);
     143      752080 :   if (is_scalar_t(t))
     144             :   {
     145       33194 :     if (!gequal1(s)) err_direuler(s);
     146       33180 :     return scalarpol_shallow(s,0);
     147             :   }
     148      718886 :   switch(t)
     149             :   {
     150        5726 :     case t_POL: break; /* no need to RgXn_red */
     151      712845 :     case t_RFRAC:
     152             :     {
     153      712845 :       GEN p = gel(s,1), q = gel(s,2);
     154      712845 :       q = RgXn_red_shallow(q,n);
     155      712845 :       s = RgXn_inv(q, n);
     156      712845 :       if (typ(p) == t_POL && varn(p) == varn(q))
     157             :       {
     158          28 :         p = RgXn_red_shallow(p, n);
     159          28 :         s = RgXn_mul(s, p, n);
     160             :       }
     161             :       else
     162      712817 :         if (!gequal1(p)) s = RgX_Rg_mul(s, p);
     163      712845 :       if (!signe(s) || !gequal1(gel(s,2))) err_direuler(s);
     164      712831 :       break;
     165             :     }
     166         308 :     case t_SER:
     167         308 :       if (!signe(s) || valser(s) || !gequal1(gel(s,2))) err_direuler(s);
     168         308 :       break;
     169           7 :     default: pari_err_TYPE("direuler", s);
     170             :   }
     171      718865 :   return s;
     172             : }
     173             : 
     174             : struct eval_bad
     175             : {
     176             :   void *E;
     177             :   GEN (*eval)(void *, GEN);
     178             : };
     179             : static GEN
     180      688303 : eval_bad(void *E, GEN p, long n)
     181             : {
     182      688303 :   struct eval_bad *d = (struct eval_bad*) E;
     183      688303 :   return direuler_factor(d->eval(d->E, p), n);
     184             : }
     185             : GEN
     186         301 : direuler(void *E, GEN (*eval)(void *, GEN), GEN a, GEN b, GEN c)
     187             : {
     188             :   struct eval_bad d;
     189         301 :   d.E= E; d.eval = eval;
     190         301 :   return direuler_bad((void*)&d, eval_bad, a, b, c, NULL);
     191             : }
     192             : 
     193             : static GEN
     194       31206 : primelist(forprime_t *T, GEN Sbad, long n, long *running)
     195             : {
     196       31206 :   GEN P = cgetg(n+1, t_VECSMALL);
     197             :   long i, j;
     198      302855 :   for (i = 1, j = 1; i <= n; i++)
     199             :   {
     200      275730 :     ulong p = u_forprime_next(T);
     201      275730 :     if (!p) { *running = 0; break; }
     202      271649 :     if (Sbad && umodiu(Sbad, p)==0) continue;
     203      266952 :     uel(P,j++) = p;
     204             :   }
     205       31206 :   setlg(P, j);
     206       31206 :   return P;
     207             : }
     208             : 
     209             : GEN
     210        4088 : pardireuler(GEN worker, GEN a, GEN b, GEN c, GEN Sbad)
     211             : {
     212             :   ulong au, bu, X, sqrtX, n, snX, nX;
     213        4088 :   pari_sp av0 = avma;
     214             :   GEN v, V;
     215             :   forprime_t T;
     216             :   struct pari_mt pt;
     217        4088 :   long running = 1, pending = 0;
     218        4088 :   au = direulertou(a, gceil);
     219        4088 :   bu = direulertou(b, gfloor);
     220        4088 :   X = c ? direulertou(c, gfloor): bu;
     221        4088 :   if (X == 0) return cgetg(1,t_VEC);
     222        4088 :   if (bu > X) bu = X;
     223        4088 :   if (!u_forprime_init(&T, au, bu)) { set_avma(av0); return mkvec(gen_1); }
     224        4081 :   v = vecsmall_ei(X, 1);
     225        4081 :   V = vec_ei(X, 1);
     226        4081 :   n = 1;
     227        4081 :   if (Sbad) Sbad = direuler_Sbad(V, v, Sbad, &n);
     228        4081 :   sqrtX = usqrt(X); snX = uprimepi(sqrtX); nX = uprimepi(X);
     229        4081 :   if (snX)
     230             :   {
     231        4067 :     GEN P = primelist(&T, Sbad, snX, &running);
     232        4067 :     GEN R = gel(closure_callgenvec(worker, mkvec2(P, utoi(X))), 2);
     233        4067 :     long i, l = lg(P);
     234       20370 :     for (i = 1; i < l; i++)
     235             :     {
     236       16303 :       GEN s = gel(R,i);
     237       16303 :       n = dirmuleuler_small(V, v, n, uel(P,i), s, lg(s));
     238             :     }
     239          14 :   } else snX = 1;
     240        4081 :   mt_queue_start_lim(&pt, worker, (nX+snX-1)/snX);
     241       34267 :   while (running || pending)
     242             :   {
     243             :     GEN done;
     244       30186 :     GEN P = running? primelist(&T, Sbad, snX, &running): NULL;
     245       30186 :     mt_queue_submit(&pt, 0, P ? mkvec2(P, utoi(X)): NULL);
     246       30186 :     done = mt_queue_get(&pt, NULL, &pending);
     247       30186 :     if (done)
     248             :     {
     249       27139 :       GEN P = gel(done,1), R = gel(done,2);
     250       27139 :       long j, l = lg(P);
     251      277788 :       for (j=1; j<l; j++)
     252             :       {
     253      250649 :         GEN F = gel(R,j);
     254      250649 :         if (degpol(F) && !gequal0(gel(F,3)))
     255      169708 :           dirmuleuler_large(V, uel(P,j), gel(F,3));
     256             :       }
     257             :     }
     258             :   }
     259        4081 :   mt_queue_end(&pt);
     260        4081 :   return gc_GEN(av0,V);
     261             : }
     262             : 
     263             : /********************************************************************/
     264             : /**                                                                **/
     265             : /**                 DIRPOWERS and DIRPOWERSSUM                     **/
     266             : /**                                                                **/
     267             : /********************************************************************/
     268             : 
     269             : /* [1^B,...,N^B] */
     270             : GEN
     271         686 : vecpowuu(long N, ulong B)
     272             : {
     273             :   GEN v;
     274             :   long p, i;
     275             :   forprime_t T;
     276             : 
     277         686 :   if (B <= 8000)
     278             :   {
     279         686 :     if (!B) return const_vec(N,gen_1);
     280         679 :     v = cgetg(N+1, t_VEC); if (N == 0) return v;
     281         679 :     gel(v,1) = gen_1;
     282         679 :     if (B == 1)
     283       92736 :       for (i = 2; i <= N; i++) gel(v,i) = utoipos(i);
     284         469 :     else if (B == 2)
     285             :     {
     286             :       ulong o, s;
     287         273 :       if (N & HIGHMASK)
     288           0 :         for (i = 2, o = 3; i <= N; i++, o += 2)
     289           0 :           gel(v,i) = addiu(gel(v,i-1), o);
     290             :       else
     291       31073 :         for (i = 2, s = 1, o = 3; i <= N; i++, s += o, o += 2)
     292       30800 :           gel(v,i) = utoipos(s + o);
     293             :     }
     294         196 :     else if (B == 3)
     295         840 :       for (i = 2; i <= N; i++) gel(v,i) = powuu(i, B);
     296             :     else
     297             :     {
     298         182 :       long k, Bk, e = expu(N);
     299        7553 :       for (i = 3; i <= N; i += 2) gel(v,i) = powuu(i, B);
     300        1239 :       for (k = 1; k <= e; k++)
     301             :       {
     302        1057 :         N >>= 1; Bk = B * k;
     303        8498 :         for (i = 1; i <= N; i += 2) gel(v, i << k) = shifti(gel(v, i), Bk);
     304             :       }
     305             :     }
     306         679 :     return v;
     307             :   }
     308           0 :   v = const_vec(N, NULL);
     309           0 :   u_forprime_init(&T, 3, N);
     310           0 :   while ((p = u_forprime_next(&T)))
     311             :   {
     312             :     long m, pk, oldpk;
     313           0 :     gel(v,p) = powuu(p, B);
     314           0 :     for (pk = p, oldpk = p; pk; oldpk = pk, pk = umuluu_le(pk,p,N))
     315             :     {
     316           0 :       if (pk != p) gel(v,pk) = mulii(gel(v,oldpk), gel(v,p));
     317           0 :       for (m = N/pk; m > 1; m--)
     318           0 :         if (gel(v,m) && m%p) gel(v, m*pk) = mulii(gel(v,m), gel(v,pk));
     319             :     }
     320             :   }
     321           0 :   gel(v,1) = gen_1;
     322           0 :   for (i = 2; i <= N; i+=2)
     323             :   {
     324           0 :     long vi = vals(i);
     325           0 :     gel(v,i) = shifti(gel(v,i >> vi), B * vi);
     326             :   }
     327           0 :   return v;
     328             : }
     329             : 
     330             : /* does x^s require log(x) ? */
     331             : static long
     332       22797 : get_needlog(GEN s)
     333             : {
     334       22797 :   switch(typ(s))
     335             :   {
     336       10157 :     case t_REAL: return 2; /* yes but not powcx */
     337        9492 :     case t_COMPLEX: return 1; /* yes using powcx */
     338        3148 :     default: return 0; /* no */
     339             :   }
     340             : }
     341             : /* [1^B,...,N^B] */
     342             : GEN
     343       21936 : vecpowug(long N, GEN B, long prec)
     344             : {
     345       21936 :   GEN v, logp = NULL;
     346       21936 :   long gp[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
     347       21936 :   long p, precp = 2, prec0, prec1, needlog;
     348             :   forprime_t T;
     349       21936 :   if (N == 1) return mkvec(gen_1);
     350       21915 :   if (typ(B) == t_INT && lgefint(B) <= 3 && signe(B) >= 0)
     351         168 :     return vecpowuu(N, itou(B));
     352       21747 :   needlog = get_needlog(B);
     353       21747 :   prec1 = prec0 = prec;
     354       21747 :   if (needlog == 1) prec1 = powcx_prec(log2((double)N), B, prec);
     355       21747 :   u_forprime_init(&T, 2, N);
     356       21747 :   v = const_vec(N, NULL);
     357       21747 :   gel(v,1) = gen_1;
     358     1688053 :   while ((p = u_forprime_next(&T)))
     359             :   {
     360             :     long m, pk, oldpk;
     361             :     GEN u;
     362     1666306 :     gp[2] = p;
     363     1666306 :     if (needlog)
     364             :     {
     365      205983 :       if (!logp)
     366       37240 :         logp = logr_abs(utor(p, prec1));
     367             :       else
     368             :       { /* Assuming p and precp are odd,
     369             :          * log p = log(precp) + 2 atanh((p - precp) / (p + precp)) */
     370      168743 :         ulong a = p >> 1, b = precp >> 1; /* p = 2a + 1, precp = 2b + 1 */
     371      168743 :         GEN z = atanhuu(a - b, a + b + 1, prec1); /* avoid overflow */
     372      168743 :         shiftr_inplace(z, 1); logp = addrr(logp, z);
     373             :       }
     374       92829 :       u = needlog == 1? powcx(gp, logp, B, prec0)
     375      205983 :                       : mpexp(gmul(B, logp));
     376      205983 :       if (p == 2) logp = NULL; /* reset: precp must be odd */
     377             :     }
     378             :     else
     379     1460323 :       u = gpow(gp, B, prec0);
     380     1666306 :     precp = p;
     381     1666306 :     gel(v,p) = u; /* p^B */
     382     1666306 :     if (prec0 != prec) gel(v,p) = gprec_wtrunc(gel(v,p), prec);
     383     3497443 :     for (pk = p, oldpk = p; pk; oldpk = pk, pk = umuluu_le(pk,p,N))
     384             :     {
     385     1831137 :       if (pk != p) gel(v,pk) = gmul(gel(v,oldpk), gel(v,p));
     386    47170846 :       for (m = N/pk; m > 1; m--)
     387    45339709 :         if (gel(v,m) && m%p) gel(v, m*pk) = gmul(gel(v,m), gel(v,pk));
     388             :     }
     389             :   }
     390       21747 :   return v;
     391             : }
     392             : 
     393             : GEN
     394         665 : dirpowers(long n, GEN x, long prec)
     395             : {
     396             :   pari_sp av;
     397             :   GEN v;
     398         665 :   if (n <= 0) return cgetg(1, t_VEC);
     399         651 :   av = avma; v = vecpowug(n, x, prec);
     400         651 :   if (typ(x) == t_INT && lgefint(x) <= 3 && signe(x) >= 0 && cmpiu(x, 2) <= 0)
     401         133 :     return v;
     402         518 :   return gc_GEN(av, v);
     403             : }
     404             : 
     405             : /* f is a totally multiplicative function of modulus 0 or 1
     406             :  * (essentially a Dirichlet character). Compute simultaneously
     407             :  * sum_{0 < n <= N} f(n)n^s and sum_{0 < n <= N} f(n)n^{-1-conj(s)}
     408             :  * Warning: s is conjugated, but not f. Main application for Riemann-Siegel,
     409             :  * where we need R(chi,s) and conj(R(chi,1-conj(s))). */
     410             : 
     411             : static GEN
     412        3234 : vecmulsqlv(GEN Q, GEN V)
     413             : {
     414             :   long l, i;
     415             :   GEN W;
     416        3234 :   if (typ(V) != t_VEC) return RgV_Rg_mul(Q, V);
     417        1239 :   l = lg(Q); W = cgetg(l, t_VEC);
     418       62244 :   for (i = 1; i < l; i++) gel(W, i) = vecmul(gel(Q, i), V);
     419        1239 :   return W;
     420             : }
     421             : 
     422             : /* P = prime divisors of (squarefree) n, V[i] = i^s for i <= sq.
     423             :  * Return NULL if n is not sq-smooth, else f(n)n^s */
     424             : static GEN
     425    18254653 : smallfact(ulong n, GEN P, ulong sq, GEN V)
     426             : {
     427             :   long i, l;
     428             :   ulong p, m, o;
     429             :   GEN c;
     430    18254653 :   if (n <= sq) return gel(V,n);
     431    18205674 :   l = lg(P); m = p = uel(P, l-1); if (p > sq) return NULL;
     432     3862093 :   for (i = l-2; i > 1; i--, m = o) { p = uel(P,i); o = m*p; if (o > sq) break; }
     433     3828786 :   c = gel(V,m); n /= m; /* m <= sq, o = m * p > sq */
     434     3828786 :   if (n > sq) { c = vecmul(c, gel(V,p)); n /= p; }
     435     3821457 :   return vecmul(c, gel(V,n));
     436             : }
     437             : 
     438             : static GEN
     439         973 : _Qtor(GEN x, long prec)
     440         973 : { return typ(x) == t_FRAC? fractor(x, prec): x; }
     441             : static GEN
     442        2457 : Qtor(GEN x, long prec)
     443             : {
     444        2457 :   long tx = typ(x);
     445        3430 :   if (tx == t_VEC || tx == t_COL) pari_APPLY_same(_Qtor(gel(x, i), prec));
     446        1589 :   return tx == t_FRAC? fractor(x, prec): x;
     447             : }
     448             : 
     449             : static GEN
     450         238 : vecf(long N, void *E, GEN (*f)(void *, ulong, long), long prec)
     451             : {
     452             :   GEN v;
     453             :   long n;
     454         238 :   if (!f) return NULL;
     455         140 :   v = cgetg(N + 1, t_VEC);
     456       33481 :   for (n = 1; n <= N; n++) gel(v,n) = f(E, n, prec);
     457         140 :   return v;
     458             : }
     459             : 
     460             : /* Here #V = #F > 0 is small. Analogous to dotproduct, with following
     461             :  * semantic differences: uses V[1] = 1; V has scalar values but F may have
     462             :  * vector values */
     463             : static GEN
     464         301 : naivedirpowerssum(GEN V, GEN F, long prec)
     465             : {
     466             :   GEN S;
     467         301 :   if (!F) S = RgV_sum(V);
     468             :   else
     469             :   {
     470         189 :     long n, l = lg(V);
     471         189 :     S = gel(F,1); /* V[1] = 1 */
     472       34447 :     for (n = 2; n < l; n++) S = gadd(S, gmul(gel(V, n), gel(F, n)));
     473             :   }
     474         301 :   return Qtor(S, prec);
     475             : }
     476             : 
     477             : /* set *p0 and *p1 to 0 and 1 in the algebra where f takes its values */
     478             : static int
     479        1092 : mk01(void *E, GEN (*f)(void*,ulong,long), long prec, GEN *p0, GEN *p1)
     480             : {
     481        1092 :   *p0 = gen_0; *p1 = gen_1;
     482        1092 :   if (!f) return 1;
     483        1057 :   *p1 = f(E, 1, prec);
     484        1057 :   if (is_vec_t(typ(*p1)))
     485             :   {
     486         406 :     long l = lg(*p1);
     487         406 :     if (l == 1) { *p0 = *p1 = NULL; return 0; }
     488         406 :     *p0 = const_vec(l-1, gen_0);
     489             :   }
     490        1057 :   return 1;
     491             : }
     492             : static GEN
     493           0 : mktrivial(long both)
     494             : {
     495           0 :   if (!both) return cgetg(1, t_VEC);
     496           0 :   return mkvec2(cgetg(1,t_VEC), cgetg(1,t_VEC));
     497             : }
     498             : 
     499             : static GEN
     500         280 : smalldirpowerssum(long N, GEN s, void *E, GEN (*f)(void *, ulong, long),
     501             :                   long both, long prec)
     502             : {
     503             :   GEN F, V, S, SB;
     504         280 :   if (!N)
     505             :   {
     506             :     GEN zerf, onef;
     507          42 :     if (!mk01(E, f, prec, &zerf, &onef)) return mktrivial(both);
     508          42 :     return zerf;
     509             :   }
     510         238 :   V = vecpowug(N, s, prec);
     511         238 :   F = vecf(N, E, f, prec);
     512         238 :   S = naivedirpowerssum(V, F, prec);
     513         238 :   if (!both) return S;
     514          84 :   if ((both==2 || !f) && gequalm1(gmul2n(real_i(s), 1)))
     515          21 :     SB = S;
     516             :   else
     517             :   {
     518          63 :     GEN FB = (both == 1 && F)? conj_i(F): F;
     519             :     long n;
     520        1876 :     for (n = 2; n <= N; n++) gel(V, n) = ginv(gmulsg(n, gel(V, n)));
     521          63 :     SB = naivedirpowerssum(V, FB, prec);
     522             :   }
     523          84 :   return mkvec2(S, SB);
     524             : }
     525             : 
     526             : static void
     527       68543 : v2unpack(GEN v, GEN *pV, GEN *pVB)
     528             : {
     529       68543 :   if (typ(v) == t_COL) { *pV = gel(v,1); *pVB = gel(v,2); }
     530       68410 :   else { *pV = v; *pVB = NULL; }
     531       68543 : }
     532             : static GEN
     533       69527 : v2pack(GEN V, GEN VB) { return VB? mkcol2(V,VB): V; }
     534             : 
     535             : static GEN
     536        1050 : dirpowsuminit(GEN s, GEN onef, GEN zerf, void *E, GEN (*f)(void *, ulong, long),              GEN data, long both)
     537             : {
     538        1050 :   long needlog = data[1], prec0 = data[2], prec1 = data[3];
     539        1050 :   ulong a, b, c, e, q, n, sq = usqrt(data[4]);
     540        1050 :   GEN V = cgetg(sq+1, t_VEC), W = cgetg(sq+1, t_VEC), VB = NULL, WB = NULL;
     541        1050 :   GEN Q = cgetg(sq+1, t_VEC), Z = cgetg(sq+1, t_VEC), QB = NULL, ZB = NULL;
     542        1050 :   GEN logp, c2, Q2, Q3, Q6, c2B = NULL, Q2B = NULL, Q3B = NULL, Q6B = NULL;
     543        1050 :   long gp[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
     544             : 
     545        1050 :   if (both == 1 || (both == 2 && !gequal(real_i(s), gneg(ghalf))))
     546          28 :   { VB = cgetg(sq+1, t_VEC); WB = cgetg(sq+1, t_VEC); QB = cgetg(sq+1, t_VEC);}
     547        1050 :   gel(V, 1) = gel(W, 1) = gel(Q, 1) = onef;
     548        1050 :   if (VB) { gel(VB, 1) = gel(WB, 1) = gel(QB, 1) = onef; }
     549        1050 :   c2 = gpow(gen_2, s, prec0); if (VB) c2B = ginv(gmul2n(conj_i(c2), 1));
     550        1050 :   if (f)
     551             :   {
     552        1029 :     GEN tmp2 = f(E, 2, prec0);
     553        1029 :     c2 = gmul(c2, tmp2); if (VB) c2B = gmul(c2B, tmp2);
     554             :   }
     555        1050 :   gel(V,2) = c2; /* f(2) 2^s */
     556        1050 :   gel(W,2) = Qtor(gadd(c2, onef), prec0);
     557        1050 :   gel(Q,2) = Qtor(gadd(vecsqr(c2), onef), prec0);
     558        1050 :   if (VB)
     559             :   {
     560          28 :     gel(VB, 2) = c2B; gel(WB, 2) = Qtor(gadd(c2B, onef), prec0);
     561          28 :     gel(QB, 2) = Qtor(gadd(vecsqr(c2B), onef), prec0);
     562             :   }
     563        1050 :   logp = NULL;
     564      158795 :   for (n = 3; n <= sq; n++)
     565             :   {
     566      157745 :     GEN u = NULL, uB = NULL, ks = f ? f(E, n, prec0) : gen_1;
     567      157745 :     if (gequal0(ks)) ks = NULL;
     568      157745 :     if (odd(n))
     569             :     {
     570       79352 :       gp[2] = n;
     571       79352 :       if (needlog)
     572             :       {
     573       77805 :         if (!logp)
     574        1029 :           logp = logr_abs(utor(n, prec1));
     575             :         else
     576             :         { /* log n = log(n-2) + 2 atanh(1 / (n - 1)) */
     577       76776 :           GEN z = atanhuu(1, n - 1, prec1);
     578       76776 :           shiftr_inplace(z, 1); logp = addrr(logp, z);
     579             :         }
     580       77805 :         if (ks)
     581       76230 :           u = needlog == 1? powcx(gp, logp, s, prec0) : mpexp(gmul(s, logp));
     582             :       }
     583        1547 :       else if (ks) u = gpow(gp, s, prec0);
     584       79352 :       if (ks)
     585             :       {
     586       77763 :         if (VB) uB = gmul(ginv(gmulsg(n, conj_i(u))), ks);
     587       77763 :         u = gmul(u, ks); /* f(n) n^s */
     588             :       }
     589             :     }
     590             :     else
     591             :     {
     592       78393 :       u = vecmul(c2, gel(V, n >> 1));
     593       78393 :       if (VB) uB = vecmul(c2B, gel(VB, n >> 1));
     594             :     }
     595      157745 :     if (ks)
     596             :     {
     597      154693 :       gel(V,n) = u; /* f(n) n^s */
     598      154693 :       gel(W,n) = gadd(gel(W, n-1), gel(V,n));        /*= sum_{i<=n} f(i)i^s */
     599      154693 :       gel(Q,n) = gadd(gel(Q, n-1), vecsqr(gel(V,n)));/*= sum_{i<=n} f(i^2)i^2s*/
     600      154693 :       if (VB)
     601             :       {
     602         483 :         gel(VB,n) = uB;
     603         483 :         gel(WB,n) = gadd(gel(WB,n-1), gel(VB,n));
     604         483 :         gel(QB,n) = gadd(gel(QB,n-1), vecsqr(gel(VB,n)));
     605             :       }
     606             :     }
     607             :     else
     608             :     {
     609        3052 :       gel(V,n) = zerf; gel(W,n) = gel(W, n-1); gel(Q,n) = gel(Q, n-1);
     610        3052 :       if (VB)
     611          21 :       { gel(VB,n) = zerf; gel(WB,n) = gel(WB, n-1); gel(QB,n) = gel(QB, n-1); }
     612             :     }
     613             :   }
     614        1050 :   Q2 = vecmulsqlv(Q, gel(V,2));
     615        1050 :   Q3 = vecmulsqlv(Q, gel(V,3));
     616        1050 :   Q6 = vecmulsqlv(Q, gel(V,6));
     617        1050 :   if (VB)
     618             :   {
     619          28 :     Q2B = vecmulsqlv(QB, gel(VB,2));
     620          28 :     Q3B = vecmulsqlv(QB, gel(VB,3));
     621          28 :     Q6B = vecmulsqlv(QB, gel(VB,6));
     622             :   }
     623             :   /* a,b,c,e = sqrt(q), sqrt(q/2), sqrt(q/3), sqrt(q/6)
     624             :    * Z[q] = Q[a] + 2^s Q[b] + 3^s Q[c] + 6^s Q[e], with Q[0] = 0 */
     625        1050 :   gel(Z, 1) = onef;
     626        1050 :   gel(Z, 2) = gel(W, 2);
     627        1050 :   gel(Z, 3) = gel(W, 3);
     628        1050 :   gel(Z, 4) = gel(Z, 5) = gel(W, 4);
     629        1050 :   gel(Z, 6) = gel(Z, 7) = gadd(gel(W, 4), gel(V, 6));
     630        1050 :   if (VB)
     631             :   {
     632          28 :     ZB = cgetg(sq+1, t_VEC);
     633          28 :     gel(ZB, 1) = onef;
     634          28 :     gel(ZB, 2) = gel(WB, 2);
     635          28 :     gel(ZB, 3) = gel(WB, 3);
     636          28 :     gel(ZB, 4) = gel(ZB, 5) = gel(WB, 4);
     637          28 :     gel(ZB, 6) = gel(ZB, 7) = gadd(gel(WB, 4), gel(VB, 6));
     638             :   }
     639        1050 :   a = 2; b = c = e = 1;
     640      153545 :   for (q = 8; q <= sq; q++)
     641             :   { /* Gray code: at most one of a,b,c,d differs (by 1) from previous value */
     642      152495 :     GEN z = gel(Z, q - 1), zB = NULL;
     643             :     ulong na, nb, nc, ne, na2, nb2, nc2, ne2;
     644      152495 :     if (VB) zB = gel(ZB, q - 1);
     645      152495 :     if ((na = usqrt(q)) != a)
     646        9191 :     { a = na; na2 = na * na; z = gadd(z, gel(V, na2));
     647        9191 :       if (VB) zB = gadd(zB, gel(VB, na2)); }
     648      143304 :     else if ((nb = usqrt(q / 2)) != b)
     649        6216 :     { b = nb; nb2 = 2 * nb * nb; z = gadd(z, gel(V, nb2));
     650        6216 :       if (VB) zB = gadd(zB, gel(VB, nb2)); }
     651      137088 :     else if ((nc = usqrt(q / 3)) != c)
     652        5418 :     { c = nc; nc2 = 3 * nc * nc; z = gadd(z, gel(V, nc2));
     653        5418 :       if (VB) zB = gadd(zB, gel(VB, nc2)); }
     654      131670 :     else if ((ne = usqrt(q / 6)) != e)
     655        3108 :     { e = ne; ne2 = 6 * ne * ne; z = gadd(z, gel(V, ne2));
     656        3108 :       if (VB) zB = gadd(zB, gel(VB, ne2)); }
     657      152495 :     gel(Z, q) = z; if (VB) gel(ZB, q) = zB;
     658             :   }
     659        1078 :   return v2pack(mkvecn(7, V, W, Q, Q2, Q3, Q6, Z),
     660          28 :                 VB? mkvecn(7, VB, WB, QB, Q2B, Q3B, Q6B, ZB): NULL);
     661             : }
     662             : 
     663             : static GEN
     664       38338 : sumprimeloop(forprime_t *pT, GEN s, long N, GEN data, GEN S, GEN W, GEN WB,
     665             :              void *E, GEN (*f)(void *, ulong, long))
     666             : {
     667       38338 :   pari_sp av = avma;
     668       38338 :   long needlog = data[1], prec0 = data[2], prec1 = data[3];
     669       38338 :   long gp[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
     670       38338 :   ulong p, precp = 0;
     671       38338 :   GEN logp = NULL, SB = WB? S: NULL;
     672             : 
     673     5281428 :   while ((p = u_forprime_next(pT)))
     674             :   {
     675     5242460 :     GEN u = NULL, ks;
     676     5242460 :     if (!f) ks = gen_1;
     677             :     else
     678             :     {
     679     5173573 :       ks = f(E, p, prec1);
     680     5173431 :       if (gequal0(ks)) ks = NULL;
     681     5173894 :       if (ks && gequal1(ks)) ks = gen_1;
     682             :     }
     683     5243073 :     gp[2] = p;
     684     5243073 :     if (needlog)
     685             :     {
     686     5166857 :       if (!logp)
     687       30419 :         logp = logr_abs(utor(p, prec1));
     688             :       else
     689             :       { /* Assuming p and precp are odd,
     690             :          * log p = log(precp) + 2 atanh((p - precp) / (p + precp)) */
     691     5136438 :         ulong a = p >> 1, b = precp >> 1; /* p = 2a + 1, precp = 2b + 1 */
     692     5136438 :         GEN z = atanhuu(a - b, a + b + 1, prec1); /* avoid overflow */
     693     5138379 :         shiftr_inplace(z, 1); logp = addrr(logp, z);
     694             :       }
     695    10332480 :       if (ks) u = needlog == 1? powcx(gp, logp, s, prec0)
     696     5166189 :                               : mpexp(gmul(s, logp));
     697             :     }
     698       76216 :     else if (ks) u = gpow(gp, s, prec0);
     699     5245069 :     if (ks)
     700             :     {
     701     5244961 :       S = gadd(S, vecmul(gel(W, N / p), ks == gen_1? u: gmul(ks, u)));
     702     5242982 :       if (WB)
     703             :       {
     704        2457 :         GEN w = gel(WB, N / p);
     705        2457 :         if (ks != gen_1) w = vecmul(ks, w);
     706        2457 :         SB = gadd(SB, gdiv(w, gmulsg(p, conj_i(u))));
     707             :       }
     708             :     }
     709     5243090 :     precp = p;
     710     5243090 :     if ((p & 0x1ff) == 1)
     711             :     {
     712       19565 :       if (!logp) (void)gc_all(av, SB? 2: 1, &S, &SB);
     713       19285 :       else (void)gc_all(av, SB? 3: 2, &S, &logp, &SB);
     714             :     }
     715             :   }
     716       38195 :   return v2pack(S, SB);
     717             : }
     718             : 
     719             : static GEN
     720       48475 : add4(GEN a, GEN b, GEN c, GEN d) { return gadd(gadd(a,b), gadd(c,d)); }
     721             : 
     722             : static const long step = 2048;
     723             : static int
     724       29557 : mksqfloop(long N, long x1, GEN R, GEN RB, GEN *pS, GEN *pSB)
     725             : {
     726       29557 :   GEN V = gel(R,1), Q = gel(R,3), Q2 = gel(R,4);
     727       29557 :   GEN Q3 = gel(R,5), Q6 = gel(R,6), Z = gel(R,7);
     728       29557 :   GEN v, VB = NULL, QB = NULL, Q2B = NULL, Q3B = NULL, Q6B = NULL, ZB = NULL;
     729       29557 :   long x2, j, lv, sq = lg(V)-1;
     730             : 
     731       29557 :   if (RB)
     732          28 :   { VB = gel(RB,1); QB = gel(RB,3); Q2B = gel(RB,4);
     733          28 :     Q3B = gel(RB,5), Q6B = gel(RB,6); ZB = gel(RB,7); }
     734             :   /* beware overflow, fuse last two bins (avoid a tiny remainder) */
     735       29557 :   x2 = (N >= 2*step && N - 2*step >= x1)? x1-1 + step: N;
     736       29557 :   v = vecfactorsquarefreeu_coprime(x1, x2, mkvecsmall2(2, 3));
     737       30043 :   lv = lg(v);
     738    60100928 :   for (j = 1; j < lv; j++)
     739    60071368 :     if (gel(v,j))
     740             :     {
     741    18256286 :       ulong d = x1 - 1 + j; /* squarefree, coprime to 6 */
     742    18256286 :       GEN t = smallfact(d, gel(v,j), sq, V), u;
     743    18226179 :       GEN tB = NULL, uB = NULL; /* = f(d) d^s */
     744             :       long a, b, c, e, q;
     745    18226179 :       if (!t || gequal0(t)) continue;
     746     3808119 :       if (VB) tB = vecinv(gmulsg(d, conj_i(t)));
     747             :       /* warning: gives 1/conj(f(d)) d^(-1-conj(s)), equal to
     748             :          f(d) d^(-1-conj(s)) only if |f(d)|=1. */
     749             :       /* S += f(d) * d^s * Z[q] */
     750     3875411 :       q = N / d;
     751     3875411 :       if (q == 1)
     752             :       {
     753     1477595 :         *pS = gadd(*pS, t); if (VB) *pSB = gadd(*pSB, tB);
     754     1479434 :         continue;
     755             :       }
     756     2397816 :       if (q <= sq) { u = gel(Z, q); if (VB) uB = gel(ZB, q); }
     757             :       else
     758             :       {
     759       48250 :         a = usqrt(q); b = usqrt(q / 2); c = usqrt(q / 3); e = usqrt(q / 6);
     760       48293 :         u = add4(gel(Q,a), gel(Q2,b), gel(Q3,c), gel(Q6,e));
     761       48293 :         if (VB) uB = add4(gel(QB,a), gel(Q2B,b), gel(Q3B,c), gel(Q6B,e));
     762             :       }
     763     2397859 :       *pS = gadd(*pS, vecmul(t, u)); if (VB) *pSB = gadd(*pSB, vecmul(tB, uB));
     764             :     }
     765       29560 :   return x2 == N;
     766             : }
     767             : 
     768             : static GEN
     769        1050 : mkdata(long N, GEN s, long prec)
     770             : {
     771        1050 :   long needlog, prec0, prec1, m = mt_nbthreads(), STEP = maxss(N / (m * m), 1);
     772        1050 :   prec1 = prec0 = prec + EXTRAPRECWORD;
     773        1050 :   needlog = get_needlog(s);
     774        1050 :   if (needlog == 1) prec1 = powcx_prec(log2((double)N), s, prec);
     775        1050 :   return mkvecsmalln(5, needlog, prec0, prec1, N, STEP);
     776             : }
     777             : 
     778             : static GEN
     779       13748 : gp_callUp(void *E, ulong x, long prec)
     780             : {
     781       13748 :   long n[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0};
     782       13748 :   n[2] = x; return gp_callprec(E, n, prec);
     783             : }
     784             : 
     785             : /* both is
     786             :  * 0: sum_{n<=N}f(n)n^s
     787             :  * 1: sum for (f,s) and (conj(f),-1-s)
     788             :  * 2: sum for (f,s) and (f,-1-s), assuming |f(n)| in {0,1} */
     789             : static GEN
     790         203 : dirpowerssumfun_i(ulong N, GEN s, void *E, GEN (*f)(void *, ulong, long),
     791             :                   long both, long prec)
     792             : {
     793             :   forprime_t T;
     794             :   pari_sp av;
     795             :   GEN onef, zerf, R, RB, W, WB, S, SB, data;
     796             :   ulong x1;
     797             : 
     798         203 :   if ((f && N < 49) || (!f && N < 1000))
     799         140 :     return smalldirpowerssum(N, s, E, f, both, prec);
     800          63 :   if (!mk01(E, f, prec, &zerf, &onef)) return mktrivial(both);
     801          63 :   data = mkdata(N, s, prec); s = gprec_w(s, prec + EXTRAPRECWORD);
     802          63 :   v2unpack(dirpowsuminit(s, onef, zerf, E, f, data, both), &R, &RB);
     803          63 :   W = gel(R,2); WB = RB? gel(RB,2): NULL;
     804          63 :   av = avma; u_forprime_init(&T, lg(W), N);
     805          63 :   v2unpack(sumprimeloop(&T, s, N, data, zerf, W, WB, E, f), &S, &SB);
     806         413 :   for(x1 = 1;; x1 += step)
     807             :   {
     808         413 :     if (mksqfloop(N, x1, R, RB, &S, &SB))
     809          63 :       return both? mkvec2(S, conj_i(SB? SB: S)): S;
     810         350 :     (void)gc_all(av, SB? 2: 1, &S, &SB);
     811             :   }
     812             : }
     813             : GEN
     814         203 : dirpowerssumfun(ulong N, GEN s, void *E, GEN (*f)(void *, ulong, long),
     815             :                 long both, long prec)
     816             : {
     817         203 :   pari_sp av = avma;
     818         203 :   return gc_GEN(av, dirpowerssumfun_i(N, s, E, f, both, prec));
     819             : }
     820             : 
     821             : GEN
     822         133 : dirpowerssum(ulong N, GEN s, long both, long prec)
     823         133 : { return dirpowerssumfun(N, s, NULL, NULL, both, prec); }
     824             : 
     825             : GEN
     826         154 : dirpowerssum0(GEN N, GEN s, GEN f, long both, long prec)
     827             : {
     828         154 :   if (typ(N) != t_INT) pari_err_TYPE("dirpowerssum", N);
     829         147 :   if (signe(N) <= 0) N = gen_0;
     830         147 :   if (!f) return dirpowerssum(itou(N), s, both, prec);
     831          70 :   if (typ(f) != t_CLOSURE) pari_err_TYPE("dirpowerssum", f);
     832          70 :   return dirpowerssumfun(itou(N), s, (void*)f, gp_callUp, both, prec);
     833             : }
     834             : 
     835             : GEN
     836       29144 : parsqf_worker(GEN gk, GEN vR, GEN gN)
     837             : {
     838       29144 :   pari_sp av = avma;
     839             :   GEN R, RB, onef, S, SB;
     840       29144 :   long k = itou(gk), N = itou(gN), x1 = 1 + step * k;
     841       29144 :   v2unpack(vR, &R, &RB); onef = gmael(R,1,1);
     842       29144 :   S = SB = is_vec_t(typ(onef)) ? zerovec(lg(onef) - 1): gen_0;
     843       29144 :   (void)mksqfloop(N, x1, R, RB, &S, &SB);
     844       29147 :   return gc_GEN(av, v2pack(S, RB? SB: NULL));
     845             : }
     846             : 
     847             : static GEN
     848     5350299 : mycallvec(void *f, ulong n, long prec)
     849             : {
     850     5350299 :   GEN F = (GEN)f;
     851     5350299 :   if (!f) return gen_1;
     852      390608 :   if (typ(F) == t_CLOSURE) return gp_callUp(f, n, prec);
     853      390608 :   return gel(F, (n-1) % (lg(F)-1) + 1);
     854             : }
     855             : 
     856             : GEN
     857       38288 : parsumprimefun_worker(GEN gk, GEN s, GEN zerf, GEN data, GEN vW, GEN f)
     858             : {
     859       38288 :   pari_sp av = avma;
     860             :   forprime_t T;
     861             :   GEN W, WB, S;
     862       38288 :   long k = itou(gk), sq, N = data[4], STEP = data[5];
     863             : 
     864       38287 :   v2unpack(vW, &W, &WB); sq = lg(W)-1;
     865       38283 :   if (isintzero(f)) f = NULL;
     866       38277 :   u_forprime_init(&T, k * STEP + sq + 1, minss(N, (k + 1) * STEP + sq));
     867       38265 :   S = sumprimeloop(&T, s, N, data, zerf, W, WB, (void*)f, mycallvec);
     868       38278 :   return gc_GEN(av, S);
     869             : }
     870             : 
     871             : static GEN
     872         987 : vR_get_vW(GEN vR)
     873             : {
     874             :   GEN R, RB, W, WB;
     875         987 :   v2unpack(vR, &R, &RB); W = gel(R,2); WB = RB? gel(RB,2): NULL;
     876         987 :   return v2pack(W, WB);
     877             : }
     878             : 
     879             : static GEN
     880         987 : halfconj(long both, GEN V)
     881         987 : { return both ? mkvec2(gel(V, 1), gconj(gel(V, 2))) : V; }
     882             : 
     883             : static GEN
     884        1127 : pardirpowerssumfun_i(GEN f, ulong N, GEN s, long both, long prec)
     885             : {
     886             :   GEN worker, worker2, data, vR, onef, zerf;
     887             : 
     888        1127 :   if ((f && N < 49) || (!f && N < 10000UL))
     889         140 :     return smalldirpowerssum(N, s, (void*)f, mycallvec, both, prec);
     890         987 :   if (!mk01((void*)f, mycallvec, prec, &zerf, &onef)) return mktrivial(both);
     891         987 :   data = mkdata(N, s, prec); s = gprec_w(s, prec + EXTRAPRECWORD);
     892         987 :   vR = dirpowsuminit(s, onef, zerf, (void*)f, mycallvec, data, both);
     893         987 :   worker = snm_closure(is_entry("_parsumprimefun_worker"),
     894             :                        mkvecn(5, s, zerf, data, vR_get_vW(vR), f? f: gen_0));
     895         987 :   worker2 = snm_closure(is_entry("_parsqf_worker"), mkvec2(vR, utoi(N)));
     896        1974 :   return halfconj(both,
     897         987 :                   gadd(parsum(gen_0, utoipos((N-1) / data[5]), worker),
     898         987 :                        parsum(gen_0, utoipos(maxss((N-1) / step - 1, 0)), worker2)));
     899             : }
     900             : 
     901             : GEN
     902        1127 : pardirpowerssumfun(GEN f, ulong N, GEN s, long both, long prec)
     903             : {
     904        1127 :   pari_sp av = avma;
     905        1127 :   return gc_GEN(av, pardirpowerssumfun_i(f, N, s, both, prec));
     906             : }
     907             : GEN
     908           0 : pardirpowerssum0(GEN N, GEN s, GEN f, long both, long prec)
     909             : {
     910           0 :   if (typ(N) != t_INT) pari_err_TYPE("pardirpowerssum", N);
     911           0 :   return pardirpowerssumfun(f, itou(N), s, both, prec);
     912             : }
     913             : GEN
     914           0 : pardirpowerssum(ulong N, GEN s, long prec)
     915           0 : { return pardirpowerssumfun(NULL, N, s, 0, prec); }

Generated by: LCOV version 1.16