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 : #include "pari.h"
15 : #include "paripriv.h"
16 :
17 : /*****************************************************************/
18 : /* Program to compute L(chi,s) */
19 : /* for Im(s) large, chi primitive Dirichlet character */
20 : /* In the present branch, only Tyagi's method is used */
21 : /*****************************************************************/
22 : /* In addition, C can also be a polynomial defining an abelian
23 : * extension of Q. */
24 :
25 : /*****************************************************************/
26 : /* Character programs */
27 : /*****************************************************************/
28 : /* A character, always assumed primitive can be given in the following formats:
29 : * - omitted or 0: special to zetaRS,
30 : * - a t_INT: assumed to be a discriminant,
31 : * - a t_INTMOD: a conrey character,
32 : * - a pair [G,chi] or [bnr,chi],
33 : * - [C1,C2,...]~ where the Ci are characters as above with same moduli. */
34 :
35 : /* Given a list of linit/ldata for chars of same conductor F, return
36 : * [Vecan, F, Parities, Gaussums] */
37 : static GEN
38 1127 : mycharinit(GEN C, long bit)
39 : {
40 : GEN L, LVC, LE, LGA;
41 1127 : long F = 0, i, j, lc = lg(C), prec;
42 :
43 1127 : bit += 64; prec = nbits2prec(bit);
44 1127 : L = cgetg(lc, t_VEC);
45 1127 : LE = cgetg(lc, t_VECSMALL);
46 1127 : LGA = cgetg(lc, t_VEC);
47 2289 : for (i = 1; i < lc; i++)
48 : {
49 1162 : GEN gv, ga, gm, ro, ldata = gel(C, i);
50 : long e;
51 1162 : if (is_linit(ldata)) ldata = linit_get_ldata(ldata);
52 1162 : gv = ldata_get_gammavec(ldata); e = itou(gel(gv, 1));
53 1162 : gm = ldata_get_conductor(ldata);
54 1162 : ro = ldata_get_rootno(ldata);
55 1162 : if (isintzero(ro)) ro = lfunrootno(ldata, bit);
56 1162 : ga = gmul(ro, gsqrt(gm, prec)); if (e) ga = mulcxI(ga);
57 1162 : gel(LGA, i) = ga;
58 1162 : LE[i] = e;
59 1162 : if (i == 1) F = itos(gm); /* constant */
60 1162 : gel(L, i) = lfunan(ldata, F, prec);
61 : }
62 1127 : if (lc == 2 && is_vec_t(typ(gmael(L,1,1))))
63 : { /* multichar */
64 7 : LGA = gel(LGA,1); lc = lg(LGA);
65 7 : LVC = gel(L,1);
66 7 : LE = const_vecsmall(lc-1, LE[1]); /* FIXME: can handle mixed values */
67 : }
68 : else
69 : {
70 1120 : LVC = cgetg(F + 1, t_VEC);
71 4627 : for (j = 1; j <= F; j++)
72 : {
73 3507 : GEN v = cgetg(lc, t_VEC);
74 7301 : for (i = 1; i < lc; i++) gel(v, i) = gmael(L, i, j);
75 3507 : gel(LVC, j) = v;
76 : }
77 : }
78 1127 : return mkvec5(LVC, stoi(F), LE, LGA, grootsof1(2*F, prec));
79 : }
80 :
81 : /* n >= 1 and #VC = F, the conductor of the character or multicharacter X.
82 : * VC contains [X(1),X(2),...X(F)] */
83 : static GEN
84 636923 : mycall(GEN VC, long n)
85 : {
86 636923 : long F = lg(VC) - 1;
87 636923 : GEN z = n <= F ? gel(VC, n) : gel(VC, ((n - 1) % F) + 1);
88 636923 : return gequal0(z)? NULL: z;
89 : }
90 :
91 205578 : static GEN get_chivec(GEN VCALL) { return gel(VCALL, 1); }
92 805484 : static long get_modulus(GEN VCALL) { return itos(gel(VCALL, 2)); }
93 204871 : static GEN get_signat(GEN VCALL) { return gel(VCALL, 3); }
94 21 : static GEN get_gauss(GEN VCALL) { return gel(VCALL, 4); }
95 203968 : static GEN get_chiZ(GEN VCALL) { return gel(VCALL, 5); }
96 :
97 : /* (-1)^A[i] * conj(B[i]) */
98 : static GEN
99 7 : gnegconj(GEN A, GEN B)
100 : {
101 7 : long i, l = lg(A);
102 7 : GEN W = cgetg(l, t_VEC);
103 14 : for (i = 1; i < l; i++)
104 7 : { GEN b = gconj(gel(B,i)); gel(W,i) = A[i]? gneg(b): b; }
105 7 : return W;
106 : }
107 : /* g(conj(CHI)) */
108 : static GEN
109 7 : gaussconj(GEN VCALL)
110 7 : { return gnegconj(get_signat(VCALL), get_gauss(VCALL)); }
111 :
112 : static GEN
113 7 : myinitconj(GEN VCALL)
114 : {
115 7 : GEN CONJ = shallowcopy(VCALL);
116 7 : gel(CONJ, 1) = gconj(get_chivec(VCALL));
117 7 : gel(CONJ, 4) = gaussconj(VCALL); return CONJ;
118 : }
119 :
120 : /********************************************************************/
121 : /* Driver Program */
122 : /********************************************************************/
123 :
124 : /* assume |Im(s)| >> 1, in particular s is not a negative integer */
125 : static GEN
126 14 : applyfuneq(GEN gau, GEN s, GEN z, long odd, long q, long bitprec)
127 : {
128 : GEN t, S, rS;
129 : long prec;
130 14 : if (!gequal0(s)) bitprec += maxss(gexpo(s), 0);
131 14 : prec = nbits2prec(bitprec);
132 14 : if (odd) gau = mulcxmI(gau);
133 14 : S = mulcxI(gmul(Pi2n(-1,prec),gsubgs(s, odd)));
134 14 : rS = real_i(S);
135 14 : if (gsigne(rS) < 0) { S = gneg(S); rS = gneg(rS); }
136 14 : if (gcmpgs(rS, bitprec) <= 0)
137 0 : S = gadd(S,glog1p(gexp(gmulgs(S,-2),prec),prec));
138 14 : t = gexp(gneg(gadd(S, glngamma(s, prec))), prec);
139 14 : t = gmul(gpow(gdivgs(Pi2n(1, prec), q), s, prec), t);
140 14 : return gmul(gmul(gau, t), z);
141 : }
142 :
143 : static GEN RZchi(GEN VCALL, GEN s, long prec);
144 :
145 : /* VCALL already initialized */
146 : static GEN
147 1127 : lfunlarge_char(GEN VCALL, GEN s, long bitprec)
148 : {
149 1127 : pari_sp av = avma;
150 : GEN sig, tau, z;
151 1127 : long funeq = 0, ts = typ(s), stau, flconj, q;
152 1127 : if (!is_real_t(ts) && ts != t_COMPLEX) pari_err_TYPE("lfunlarge_char", s);
153 1127 : sig = real_i(s); tau = imag_i(s);
154 1127 : if (gexpo(tau) < 1) pari_err_DOMAIN("lfun","im(s)", "<", gen_2, tau);
155 1127 : stau = gsigne(tau);
156 1127 : if (stau < 0) { tau = gneg(tau); VCALL = myinitconj(VCALL); }
157 1127 : if (gcmp(sig, ghalf) < 0) { funeq = 1; sig = gsubsg(1, sig); }
158 1127 : flconj = ((stau > 0 && funeq) || (stau < 0 && !funeq));
159 1127 : q = get_modulus(VCALL); bitprec += gexpo(stoi(q));
160 1127 : z = RZchi(VCALL, mkcomplex(sig, tau), nbits2prec(bitprec));
161 1127 : if (flconj) z = gconj(z);
162 1127 : if (funeq)
163 : {
164 14 : GEN odd = get_signat(VCALL), gau = get_gauss(VCALL), Vz;
165 14 : long lC = lg(gau), j;
166 14 : Vz = cgetg(lC, t_VEC);
167 28 : for (j = 1; j < lC; j++)
168 14 : gel(Vz,j) = applyfuneq(gel(gau,j), s, gel(z,j), odd[j], q, bitprec);
169 14 : z = Vz;
170 : }
171 1127 : return gc_GEN(av, z);
172 : }
173 :
174 : static GEN
175 14 : lfungetchars(GEN pol)
176 : {
177 14 : GEN L, F, v, bnf = Buchall(pol_x(1), 0, LOWDEFAULTPREC);
178 : GEN w, condall, bnr;
179 : long i, l, lc;
180 14 : condall = rnfconductor(bnf, pol); bnr = gel(condall, 2);
181 14 : L = bnrchar(bnr, gel(condall, 3), NULL); lc = lg(L);
182 14 : F = cgetg(lc, t_VEC);
183 77 : for (i = 1; i < lc; i++)
184 : {
185 63 : GEN chi = gel(L, i), cond = bnrconductor_raw(bnr, chi);
186 63 : gel(F, i) = gcoeff(gel(cond,1), 1, 1);
187 : }
188 14 : w = vec_equiv(F); l = lg(w); v = cgetg(l, t_COL);
189 42 : for (i = 1; i < l; i++)
190 : {
191 28 : GEN wi = gel(w, i), vi;
192 28 : long j, li = lg(wi);
193 28 : gel(v,i) = vi = cgetg(li, t_VEC);
194 28 : if (li == 2 && equali1(gel(F, wi[1]))) /* common conductor is 1 */
195 14 : gel(vi,1) = lfunmisc_to_ldata_shallow(gen_1);
196 : else
197 : {
198 63 : for (j = 1; j < li; j++)
199 49 : gel(vi,j) = lfunmisc_to_ldata_shallow(mkvec2(bnr, gel(L, wi[j])));
200 : }
201 : }
202 14 : return v;
203 : }
204 :
205 : /********************************************************************/
206 : /* NEW RS IMPLEMENTATION FROM SANDEEP TYAGI */
207 : /********************************************************************/
208 : /* See arXiv:2203.02509v2 */
209 :
210 5775 : static long m_n0(GEN sel) { return itos(gel(sel, 1)); }
211 599514 : static GEN m_r0(GEN sel) { return gel(sel, 2); }
212 1127 : static GEN m_al(GEN sel) { return gel(sel, 3); }
213 599514 : static GEN m_aleps(GEN sel) { return gel(sel, 4); }
214 4606 : static GEN m_h(GEN sel) { return gel(sel, 5); }
215 1176 : static GEN m_lin(GEN sel) { return gel(sel, 6); }
216 1176 : static long m_np(GEN sel) { return itou(gel(sel, 7)); }
217 :
218 : static GEN
219 3444 : phi_hat(GEN x, long prec)
220 : {
221 : GEN y;
222 3444 : if (gsigne(imag_i(x)) > 0)
223 1141 : y = gneg(gexpm1(gneg(gmul(PiI2(prec), x)), prec));
224 : else
225 2303 : y = gexpm1(gmul(PiI2(prec), x), prec);
226 3444 : return ginv(y);
227 : }
228 :
229 : static GEN
230 3388 : phi_hat_h0(GEN sel, long k, long prec)
231 : {
232 3388 : GEN t = gdiv(gsubsg(m_n0(sel) + k, m_r0(sel)), m_aleps(sel));
233 3388 : return phi_hat(gdiv(gasinh(t, prec), m_h(sel)), prec);
234 : }
235 :
236 : /* v[i] = A[i] * (a + (-1)^E[i] b) */
237 : static GEN
238 630952 : mul_addsub(GEN A, GEN a, GEN b, GEN E)
239 : {
240 630952 : long i, l = lg(E);
241 630952 : GEN v = cgetg(l, t_VEC);
242 1345352 : for (i = 1; i < l; i++)
243 714400 : gel(v,i) = gmul(gel(A,i), E[i]? gsub(a, b): gadd(a, b));
244 630952 : return v;
245 : }
246 :
247 : static GEN
248 202841 : wd(GEN VCALL, GEN xpmd, long prec)
249 : {
250 202841 : GEN VC = get_chivec(VCALL), E = get_signat(VCALL), Z = get_chiZ(VCALL);
251 202841 : GEN ex, emx, y = NULL;
252 202841 : long md = get_modulus(VCALL), N = 2*md, k;
253 202841 : ex = gexp(mulcxI(xpmd), prec); emx = ginv(ex);
254 833793 : for (k = 1; k <= (md-1) / 2; k++)
255 : {
256 630952 : GEN xc = mycall(VC, k);
257 630952 : if (xc)
258 : {
259 630952 : GEN p3, p2, p1 = gmul(xc, gel(Z, Fl_neg(Fl_sqr(k,N), N) + 1));
260 630952 : GEN a = gmul(ex, gel(Z, N - k + 1)), b = gmul(emx, gel(Z, k + 1));
261 630952 : GEN c = gmul(ex, gel(Z, k + 1)), d = gmul(emx, gel(Z, N - k + 1));
262 630952 : if (odd(md))
263 : {
264 594831 : p2 = ginv(mulcxmI(gmul2n(gsub(a,b), -1))); /* 1 / sin(xpmd - kpmd) */
265 594831 : p3 = ginv(mulcxmI(gmul2n(gsub(c,d), -1))); /* 1 / sin(xpmd + kpmd) */
266 : }
267 : else
268 : {
269 36121 : p2 = mulcxI(gdiv(gadd(a,b), gsub(a,b))); /* cotan(xpmd - kpmd) */
270 36121 : p3 = mulcxI(gdiv(gadd(c,d), gsub(c,d))); /* cotan(xpmd + kpmd) */
271 : }
272 630952 : p1 = mul_addsub(p1, p2, p3, E);
273 630952 : y = y ? gadd(y, p1) : p1;
274 : }
275 : }
276 202841 : return mulcxmI(gdivgs(y, N));
277 : }
278 :
279 : static GEN
280 1127 : series_h0(long n0, GEN s, GEN VCALL, long fl, long prec)
281 : {
282 1127 : GEN C = get_modulus(VCALL) == 1? NULL: get_chivec(VCALL);
283 1127 : GEN R = pardirpowerssumfun(C, n0, gneg(s), fl, prec);
284 1127 : if (C) return R;
285 700 : if (fl) return mkvec2(mkvec(gel(R,1)), mkvec(gel(R,2)));
286 665 : return mkvec(R);
287 : }
288 :
289 : static GEN
290 1176 : series_residues_h0(GEN sel, GEN s, GEN VCALL, long prec)
291 : {
292 1176 : long n0 = m_n0(sel), np = m_np(sel), k;
293 1176 : GEN val = gen_0, VC = get_chivec(VCALL);
294 4676 : for (k = maxss(1 - np, 1 - n0); k <= 1 + np; k++)
295 : {
296 3500 : GEN nk = mycall(VC, n0 + k); /* n0 + k > 0 */
297 3500 : if (nk) val = gadd(val, gmul(gmul(phi_hat_h0(sel, k, prec), nk),
298 : gpow(stoi(n0 + k), gneg(s), prec)));
299 : }
300 1176 : return val;
301 : }
302 :
303 : static GEN
304 596126 : integrand_h0(GEN sel, GEN s, GEN VCALL, GEN x, long prec)
305 : {
306 596126 : pari_sp av = avma;
307 596126 : long md = get_modulus(VCALL);
308 596126 : GEN r0 = m_r0(sel), aleps = m_aleps(sel), zn, p1;
309 596126 : GEN znpmd, pmd = divru(mppi(prec), md), ix = ginv(x);
310 :
311 596126 : zn = gadd(r0, gdivgs(gmul(aleps, gsub(x, ix)), 2));
312 596126 : znpmd = gmul(pmd, zn);
313 596126 : p1 = gsub(mulcxI(gmul(znpmd, zn)), gmul(s, glog(zn, prec)));
314 596126 : p1 = gmul(gexp(p1, prec), gmul(aleps, gadd(x, ix)));
315 596126 : if (md == 1)
316 393285 : p1 = gdiv(mkvec(mulcxI(p1)), gmul2n(gsin(znpmd, prec), 2));
317 : else
318 202841 : p1 = gdivgs(gmul(p1, wd(VCALL, znpmd, prec)), -2);
319 596126 : return gc_upto(av, p1);
320 : }
321 :
322 : static GEN
323 1176 : integral_h0(GEN sel, GEN s, GEN VCALL, long prec)
324 : {
325 1176 : GEN lin_grid = m_lin(sel), S = gen_0;
326 1176 : pari_sp av = avma;
327 1176 : long j, l = lg(lin_grid);
328 592794 : for (j = 1; j < l; j++)
329 : {
330 591618 : S = gadd(S, integrand_h0(sel, s, VCALL, gel(lin_grid, j), prec));
331 591618 : if ((j & 0xff) == 0) S = gc_upto(av, S);
332 : }
333 1176 : return gc_upto(av, gmul(m_h(sel), S));
334 : }
335 :
336 : /* a + log |x|, a t_REAL, low accuracy */
337 : static GEN
338 4676 : addrlogabs(GEN a, GEN x)
339 : {
340 4676 : long prec = DEFAULTPREC;
341 4676 : if (gequal0(x)) return real_m2n(64, prec); /* -oo */
342 4550 : switch(typ(x))
343 : {
344 4550 : case t_COMPLEX: return gmul2n(glog(cxnorm(x), prec), -1);
345 0 : case t_REAL: break;
346 0 : default: x = gtofp(x, prec);
347 : }
348 0 : return addrr(a, logr_abs(x));
349 : }
350 :
351 : struct fun_q_t { GEN sel, s, VCALL, B; };
352 : static GEN
353 4508 : fun_q(void *E, GEN x)
354 : {
355 4508 : struct fun_q_t *S = (struct fun_q_t *)E;
356 4508 : GEN z = integrand_h0(S->sel,S->s,S->VCALL, gexp(x,DEFAULTPREC), DEFAULTPREC);
357 4508 : if (typ(z) == t_VEC) z = vecsum(z);
358 4508 : return addrlogabs(S->B, z);
359 : }
360 : static GEN
361 2338 : brent_q(void *E, GEN (*f)(void*,GEN), GEN q_low, GEN q_hi)
362 : {
363 2338 : if (gsigne(f(E, q_low)) * gsigne(f(E, q_hi)) >= 0) return NULL;
364 0 : return zbrent(E, f, q_low, q_hi, LOWDEFAULTPREC);
365 : }
366 : static GEN
367 1169 : findq(void *E, GEN (*f)(void*,GEN), GEN lq, long prec)
368 : {
369 1169 : GEN q_low, q_hi, q_right, q_left, q_est = gasinh(lq, LOWDEFAULTPREC);
370 1169 : q_low = gdivgs(gmulsg(4, q_est), 5);
371 1169 : q_hi = gdivgs(gmulsg(3, q_est), 2);
372 1169 : q_right = brent_q(E, f, q_low, q_hi); if (!q_right) q_right = q_est;
373 1169 : q_left = brent_q(E, f, gneg(q_low), gneg(q_hi)); if (!q_left) q_left = q_est;
374 1169 : return gprec_w(gmax(q_right, q_left), prec);
375 : }
376 :
377 : static GEN
378 1127 : set_q_value(GEN sel, GEN s, GEN VCALL, long prec)
379 : {
380 : struct fun_q_t E;
381 1127 : GEN al = m_al(sel), lq;
382 1127 : long md = get_modulus(VCALL), LD = DEFAULTPREC;
383 1127 : E.sel = sel; E.s = s; E.VCALL = VCALL, E.B = mulur(prec, mplog2(LD));
384 1127 : lq = gdiv(gsqrt(gdiv(gmulsg(md, E.B), Pi2n(1, LD)), LD), al);
385 1127 : return findq((void*)&E, &fun_q, lq, prec);
386 : }
387 :
388 : static GEN
389 1141 : setlin_grid_exp(GEN h, long m, long prec)
390 : {
391 1141 : GEN w, vex = gpowers(gexp(h, prec), (m - 1)/2);
392 : long i;
393 1141 : w = cgetg(m+1, t_VEC); gel(w, (m + 1)/2) = gen_1;
394 290838 : for (i = (m + 3)/2; i <= m; i++)
395 : {
396 289697 : GEN t1 = gel(vex, i - ((m - 1)/2));
397 289697 : gel(w, i) = t1; gel(w, (m + 1) - i) = ginv(t1);
398 : }
399 1141 : return w;
400 : }
401 :
402 : static long
403 1141 : get_m(GEN q, long prec)
404 : {
405 1141 : GEN t = divrr(mulur(4 * prec, mplog2(prec)), sqrr(mppi(prec)));
406 1141 : return 2*itos(gfloor(mulrr(q, t))) + 1;
407 : }
408 :
409 : static GEN
410 1127 : RZinit(GEN s, GEN VCALL, GEN numpoles, long prec)
411 : {
412 : GEN sel, al, aleps, n0, r0, q, h;
413 1127 : long md = get_modulus(VCALL), m;
414 1127 : al = gcmpgs(gabs(imag_i(s), prec), 100) < 0 ? ginv(stoi(4)) : gen_1;
415 1127 : r0 = gsqrt(gdiv(gmulgs(s, md), PiI2(prec)), prec);
416 1127 : n0 = gfloor(gsub(real_i(r0), imag_i(r0)));
417 1127 : aleps = gmul(al, gexp(PiI2n(-2, prec), prec));
418 1127 : sel = mkvecn(7, n0, r0, al, aleps, NULL, NULL, numpoles);
419 1127 : q = set_q_value(sel, s, VCALL, prec);
420 1127 : m = get_m(q, prec);
421 1127 : gel(sel,5) = h = divru(q, (m - 1) >> 1);
422 1127 : gel(sel,6) = setlin_grid_exp(h, m, prec);
423 1127 : return sel;
424 : }
425 :
426 : static GEN
427 469 : xpquo_one(GEN s, GEN cs, GEN ga, long odd, long md, long prec)
428 : {
429 469 : GEN rho, a = odd? gen_1: gen_0, z = divsr(md, mppi(prec));
430 469 : rho = gmul(gdiv(gpow(gen_I(), gdivgs(gneg(a), 2), prec), gsqrt(ga, prec)),
431 : gpow(stoi(md), ginv(stoi(4)), prec));
432 469 : return gmul(gdiv(gconj(gmul(rho, gpow(z, gdivgs(cs, 2), prec))),
433 : gmul(rho, gpow(z, gdivgs(s, 2), prec))),
434 : gexp(gsub(gconj(glngamma(gdivgs(gadd(cs, a), 2), prec)),
435 : glngamma(gdivgs(gadd(s, a), 2), prec)), prec));
436 : }
437 :
438 : static GEN
439 1127 : xpquo(GEN s, GEN VCALL, long prec)
440 : {
441 1127 : GEN ve = get_signat(VCALL), VC = get_chivec(VCALL), Z = get_chiZ(VCALL);
442 1127 : GEN R, cs, cd = NULL;
443 1127 : long md = get_modulus(VCALL), lve = lg(ve), n, i;
444 1127 : if (!gequal0(s)) prec = nbits2prec(prec + maxss(gexpo(s), 0));
445 1127 : if (md == 1)
446 700 : return gmul(gpow(mppi(prec), gsub(s, ghalf), prec),
447 : gexp(gsub(glngamma(gdivgs(gsubsg(1, s), 2), prec),
448 : glngamma(gdivgs(s, 2), prec)), prec));
449 2898 : for (n = 1; n < md; n++)
450 : {
451 2471 : GEN xn = mycall(VC, n);
452 2471 : if (xn)
453 : {
454 2408 : xn = gmul(xn, gel(Z, 2*(md - n) + 1));
455 2408 : cd = cd ? gadd(cd, xn) : xn;
456 : }
457 : }
458 427 : cs = gsubsg(1, gconj(s)); R = cgetg(lve, t_VEC);
459 896 : for (i = 1; i < lve; i++)
460 469 : gel(R, i) = xpquo_one(s, cs, gel(cd, i), ve[i], md, prec);
461 427 : if (lve == 2) R = gel(R, 1);
462 427 : return R;
463 : }
464 :
465 : static GEN
466 1176 : total_value(GEN serh0, GEN sel, GEN s, GEN VCALL, long prec)
467 : {
468 1176 : return gadd(integral_h0(sel, s, VCALL, prec),
469 : gsub(serh0, series_residues_h0(sel, s, VCALL, prec)));
470 : }
471 : static GEN
472 1127 : dirichlet_ours(GEN s, GEN VCALL, long prec)
473 : {
474 1127 : int fl = !gequal(real_i(s), ghalf);
475 1127 : GEN sel = RZinit(s, VCALL, gen_1, prec);
476 1127 : GEN S1, S2, serh0 = series_h0(m_n0(sel), s, VCALL, fl, prec);
477 1127 : if (!fl)
478 1078 : S2 = S1 = total_value(serh0, sel, s, VCALL, prec);
479 : else
480 : {
481 49 : S1 = total_value(gel(serh0,1), sel, s, VCALL, prec);
482 49 : S2 = total_value(gconj(gel(serh0,2)), sel, gsubsg(1, gconj(s)), VCALL, prec);
483 : }
484 1127 : return gadd(S1, vecmul(xpquo(s, VCALL, prec), gconj(S2)));
485 : }
486 :
487 : /* assume |Im(s)| > 2^-bitprec */
488 : static GEN
489 1127 : RZchi(GEN VCALL, GEN s, long prec)
490 : {
491 1127 : long prec2 = prec + EXTRAPREC64;
492 1127 : return gprec_wtrunc(dirichlet_ours(gprec_w(s, prec2), VCALL, prec2), prec);
493 : }
494 :
495 : /********************************************************************/
496 : /* Utility Functions */
497 : /********************************************************************/
498 : /* lam = 0, return L(s); else Lambda(s) */
499 : static GEN
500 1127 : lfuncharall(GEN VCALL, GEN s, long lam, long bitprec)
501 : {
502 1127 : GEN ve, P, Q, R, z = lfunlarge_char(VCALL, s, bitprec);
503 : long l, i, q, prec;
504 1127 : if (!lam) return z;
505 882 : ve = get_signat(VCALL); l = lg(ve);
506 882 : q = get_modulus(VCALL); prec = nbits2prec(bitprec);
507 882 : R = cgetg(l, t_VEC);
508 882 : Q = divur(q, mppi(prec));
509 882 : P = (q == 1 || zv_equal0(ve))? NULL: gsqrt(utoipos(q), prec);
510 1764 : for (i = 1; i < l; i++)
511 : {
512 882 : GEN se = gmul2n(gaddgs(s, ve[i]), -1), r;
513 882 : if (lam == 1)
514 : {
515 0 : r = gmul(gpow(Q, se, prec), ggamma(se, prec));
516 0 : if (P && ve[i]) r = gdiv(r, P);
517 : }
518 : else
519 : {
520 882 : r = gadd(gmul(se, glog(Q, prec)), glngamma(se, prec));
521 882 : if (P && ve[i]) r = gsub(r, glog(P, prec));
522 : }
523 882 : gel(R, i) = r;
524 : }
525 882 : return lam == 1 ? vecmul(R, z) : gadd(R, glog(z, prec));
526 : }
527 :
528 : static GEN
529 28 : lfunlargeall_from_chars(GEN v, GEN s, long lam, long bit)
530 : {
531 28 : long i, l = lg(v);
532 84 : for (i = 1; i < l; i++)
533 : {
534 56 : GEN w = mycharinit(gel(v, i), bit), L = lfuncharall(w, s, lam, bit);
535 56 : gel(v, i) = lam==-1 ? vecsum(L): vecprod(L);
536 : }
537 28 : return lam==-1 ? vecsum(v): vecprod(v);
538 : }
539 : static GEN
540 1099 : lfunlargeall(GEN ldata, GEN s, long lam, long bit)
541 : {
542 : GEN w, an;
543 1099 : if (lg(ldata) == 2)
544 : { /* HACK: ldata[1] a t_DESC_PRODUCT from lfunabelianrelinit / Q */
545 14 : GEN v = lfunprod_get_fact(linit_get_tech(gel(ldata,1)));
546 : long i, l;
547 14 : v = shallowcopy(gel(v,1)); l = lg(v);
548 42 : for (i = 1; i < l; i++) gel(v,i) = mkvec(gel(v,i));
549 14 : return lfunlargeall_from_chars(v, s, lam, bit);
550 : }
551 1085 : an = gel(ldata_get_an(ldata), 2);
552 1085 : switch(ldata_get_type(ldata))
553 : {
554 14 : case t_LFUN_NF:
555 : {
556 14 : GEN v = lfungetchars(nf_get_pol(an));
557 14 : return lfunlargeall_from_chars(v, s, lam, bit);
558 : }
559 7 : case t_LFUN_CHIGEN:
560 : {
561 7 : GEN chi = gmael(an, 2, 2);
562 7 : if (lg(chi) > 1 && is_vec_t(typ(gel(chi,1))))
563 : { /* multi char */
564 7 : w = mycharinit(mkcol(ldata), bit);
565 7 : return lfuncharall(w, s, lam, bit);
566 : }
567 : }
568 : default: /* single char */
569 1064 : w = mycharinit(mkcol(ldata), bit);
570 1064 : return gel(lfuncharall(w, s, lam, bit), 1);
571 : }
572 : }
573 :
574 : GEN
575 231 : lfunlarge(GEN CHI, GEN s, long bit)
576 231 : { return lfunlargeall(CHI, s, 0, bit); }
577 :
578 : GEN
579 0 : lfunlambdalarge(GEN CHI, GEN s, long bit)
580 0 : { return lfunlargeall(CHI, s, 1, bit); }
581 :
582 : GEN
583 868 : lfunloglambdalarge(GEN CHI, GEN s, long bit)
584 868 : { return lfunlargeall(CHI, s, -1, bit); }
585 :
586 : /********************************************************************/
587 : /* LERCH RS IMPLEMENTATION FROM SANDEEP TYAGI */
588 : /********************************************************************/
589 :
590 : static GEN
591 56 : double_exp_residue_pos_h(GEN selsm, long k, long ind, long prec)
592 : {
593 56 : long nk = itos(gel(selsm, 1)) + k;
594 56 : GEN r = gel(selsm, 2), ale = gel(selsm, 3), aor = gel(selsm, 4);
595 56 : GEN h = gel(selsm, 5), t = gen_0;
596 56 : switch(ind)
597 : {
598 28 : case 0: t = gaddsg(nk, aor); break;
599 0 : case 1: t = gneg(gaddsg(nk, aor)); break;
600 28 : case 2: t = gsubsg(nk, aor); break;
601 : }
602 56 : return gdiv(gasinh(gdiv(gsub(t, r), ale), prec), h);
603 : }
604 :
605 : static GEN
606 56 : phi_hat_h(GEN selsm, long m, long ind, long prec)
607 56 : { return phi_hat(double_exp_residue_pos_h(selsm, m, ind, prec), prec); }
608 :
609 : static long
610 56 : myex(GEN is) { return gequal0(is) ? 0 : maxss(0, 2 + gexpo(is)); }
611 : static GEN
612 56 : gaminus(GEN s, long prec)
613 : {
614 56 : GEN is = imag_i(s), tmp;
615 : long prec2;
616 56 : if (gcmpgs(is, -5*prec) < 0) return gen_0;
617 56 : prec2 = nbits2prec(prec + myex(is));
618 56 : tmp = gexp(gsub(glngamma(s, prec2), gmul(PiI2n(-1, prec2), s)), prec2);
619 56 : return gprec_w(tmp, prec);
620 : }
621 : static GEN
622 28 : gaplus(GEN s, long prec)
623 : {
624 28 : GEN is = imag_i(s), tmp;
625 : long prec2;
626 28 : if (gcmpgs(is, 5*prec) > 0) return gen_0;
627 0 : prec2 = nbits2prec(prec + myex(is));
628 0 : tmp = gexp(gadd(glngamma(s, prec2), gmul(PiI2n(-1, prec2), s)), prec2);
629 0 : return gprec_w(tmp, prec);
630 : }
631 :
632 : GEN
633 3513 : serh_worker(GEN k, GEN z, GEN a, GEN ns, GEN gprec)
634 : {
635 3513 : long prec = itou(gprec);
636 3513 : return gmul(gpow(z, k, prec), gpow(gadd(a, k), ns, prec));
637 : }
638 :
639 : static void
640 28 : set_arg(GEN worker, GEN z, GEN a, GEN ns, long prec)
641 28 : { gel(worker, 7) = mkvec4(z, a, ns, utoi(prec)); }
642 :
643 :
644 : static GEN
645 14 : series_h0l(GEN worker, long n0, GEN s, GEN a, GEN lam, long prec)
646 : {
647 14 : GEN z = typ(lam) == t_INT ? gen_1 : gexp(gmul(PiI2(prec), lam), prec);
648 14 : set_arg(worker, z, a, gneg(s), prec);
649 14 : return parsum(gen_0, utoi(n0), worker);
650 : }
651 :
652 : static GEN
653 14 : series_h1(GEN worker, long n1, GEN s, GEN a, GEN lam, long prec)
654 : {
655 14 : GEN mP, pre_factor, z, sn = gsubgs(s, 1);
656 14 : GEN ini = gequal0(lam) ? gen_1 : gen_0;
657 14 : pre_factor = gaplus(gneg(sn), prec);
658 14 : if (gequal0(pre_factor)) return gen_0;
659 0 : mP = gneg(PiI2(prec));
660 0 : pre_factor = gmul(gmul(pre_factor, gexp(gmul(mP, gmul(a, lam)), prec)),
661 : gpow(Pi2n(1, prec), sn, prec));
662 0 : z = typ(a) == t_INT ? gen_1 : gexp(gmul(mP, a), prec);
663 0 : set_arg(worker, z, lam, sn, prec);
664 0 : return gmul(pre_factor, parsum(ini, stoi(n1 - 1), worker));
665 : }
666 :
667 : static GEN
668 14 : series_h2(GEN worker, long n2, GEN s, GEN a, GEN lam, long prec)
669 : {
670 14 : GEN P, pre_factor, z, sn = gsubgs(s, 1);
671 14 : pre_factor = gaminus(gneg(sn), prec);
672 14 : if (gequal0(pre_factor)) return gen_0;
673 14 : P = PiI2(prec);
674 14 : pre_factor = gmul(gmul(pre_factor, gexp(gmul(gneg(P), gmul(a, lam)), prec)),
675 : gpow(Pi2n(1, prec), sn, prec));
676 14 : z = typ(a) == t_INT ? gen_1 : gexp(gmul(P, a), prec);
677 14 : set_arg(worker, z, gneg(lam), sn, prec);
678 14 : return gmul(pre_factor, parsum(gen_1, stoi(n2), worker));
679 : }
680 :
681 : static GEN
682 14 : series_residues_h0l(long numpoles, GEN selsm0, GEN s, GEN a, GEN lam, long prec)
683 : {
684 14 : GEN P, val = gen_0, ra = real_i(a);
685 14 : long n0 = m_n0(selsm0), k;
686 14 : P = PiI2(prec);
687 42 : for (k = -numpoles + 1; k <= numpoles; k++)
688 28 : if (gsigne(gaddsg(n0 + k, ra)) > 0)
689 28 : val = gadd(val, gmul(gmul(phi_hat_h(selsm0, k, 0, prec),
690 : gexp(gmul(P, gmulgs(lam, n0 + k)), prec)),
691 : gpow(gaddsg(n0 + k, a), gneg(s), prec)));
692 14 : return val;
693 : }
694 :
695 : static GEN
696 14 : series_residues_h1(long numpoles, GEN selsm1, GEN s, GEN a, GEN lam, long prec)
697 : {
698 14 : GEN mP, val = gen_0, rlam = real_i(lam), pre_factor, sn = gsubgs(s, 1);
699 14 : long n1 = m_n0(selsm1), k;
700 14 : pre_factor = gaplus(gneg(sn), prec);
701 14 : if (gequal0(pre_factor)) return gen_0;
702 0 : mP = gneg(PiI2(prec));
703 0 : pre_factor = gmul(gmul(pre_factor, gexp(gmul(mP, gmul(a, lam)), prec)),
704 : gpow(Pi2n(1, prec), sn, prec));
705 0 : for (k = -numpoles; k <= numpoles - 1; k++)
706 0 : if (gsigne(gaddsg(n1 + k, rlam)) > 0)
707 0 : val = gadd(val, gmul(gmul(phi_hat_h(selsm1, k, 1, prec),
708 : gexp(gmul(mP, gmulgs(a, n1 + k)), prec)),
709 : gpow(gaddsg(n1 + k, lam), sn, prec)));
710 0 : return gmul(pre_factor, val);
711 : }
712 :
713 : static GEN
714 14 : series_residues_h2(long numpoles, GEN selsm2, GEN s, GEN a, GEN lam, long prec)
715 : {
716 14 : GEN P, val = gen_0, rlam = real_i(lam), pre_factor, sn = gsubgs(s, 1);
717 14 : long n2 = m_n0(selsm2), k;
718 14 : pre_factor = gaminus(gneg(sn), prec);
719 14 : if (gequal0(pre_factor)) return gen_0;
720 14 : P = PiI2(prec);
721 14 : pre_factor = gmul(gmul(pre_factor, gexp(gmul(gneg(P), gmul(a, lam)), prec)),
722 : gpow(Pi2n(1, prec), sn, prec));
723 42 : for (k = -numpoles + 1; k <= numpoles; k++)
724 28 : if (gsigne(gsubsg(n2 + k, rlam)) > 0)
725 28 : val = gsub(val, gmul(gmul(phi_hat_h(selsm2, k, 2, prec),
726 : gexp(gmul(P, gmulgs(a, n2 + k)), prec)),
727 : gpow(gsubsg(n2 + k, lam), sn, prec)));
728 14 : return gmul(pre_factor, val);
729 : }
730 :
731 : static GEN
732 3430 : integrand_h0l(GEN selsm0, GEN s, GEN alam1, GEN x, long prec)
733 : {
734 3430 : GEN r0 = gel(selsm0, 2), ale = gel(selsm0, 3), a = gel(selsm0, 4);
735 3430 : GEN ix = ginv(x), zn = gadd(r0, gmul2n(gmul(ale, gsub(x, ix)), -1));
736 3430 : GEN P = PiI2n(0, prec), den, num;
737 3430 : den = gexpm1(gmul(P, gmul2n(gsub(zn,a), 1)), prec);
738 3430 : num = gexp(gmul(gmul(P, zn), gsub(alam1, zn)), prec);
739 3430 : num = gmul(gmul(gmul(num, ale), gmul2n(gadd(x, ix), -1)),
740 : gpow(zn, gneg(s), prec));
741 3430 : return gdiv(num, den);
742 : }
743 :
744 : static GEN
745 6860 : integrand_h12(GEN selsm1, GEN s, GEN alam1, GEN x, long prec)
746 : {
747 6860 : GEN r1 = gel(selsm1, 2), ale = gel(selsm1, 3), lam = gel(selsm1, 4);
748 6860 : GEN ix = ginv(x), zn = gadd(r1, gmul2n(gmul(ale, gsub(x, ix)), -1));
749 6860 : GEN P = PiI2n(0, prec), den, num, y;
750 6860 : den = gexpm1(gmul(P, gmul2n(gadd(zn,lam), 1)), prec);
751 6860 : num = gexp(gmul(gmul(P, zn), gadd(alam1, zn)), prec);
752 6860 : num = gmul(gmul(gmul(num, ale), gmul2n(gadd(x, ix), -1)),
753 : gpow(zn, gsubgs(s, 1), prec));
754 6860 : y = gdiv(num, den);
755 6860 : if (gcmp(garg(zn, prec), Pi2n(-2, prec)) > 0)
756 1722 : y = gmul(y, gexp(gmul(PiI2(prec), gsubsg(1, s)), prec));
757 6860 : return y;
758 : }
759 :
760 : static GEN
761 14 : integral_h0l(GEN lin_grid, GEN selsm0, GEN s, GEN a, GEN lam, long prec)
762 : {
763 14 : GEN A = gaddgs(gmul2n(gadd(a, lam),1), 1), S = gen_0;
764 14 : pari_sp av = avma;
765 14 : long j, l = lg(lin_grid);
766 :
767 3388 : for (j = 1; j < l; j++)
768 : {
769 3374 : S = gadd(S, integrand_h0l(selsm0, s, A, gel(lin_grid, j), prec));
770 3374 : if ((j & 0xff) == 0) S = gc_upto(av, S);
771 : }
772 14 : S = gmul(m_h(selsm0), S);
773 14 : A = gmul(a, gaddsg(1, gadd(a, gmul2n(lam, 1))));
774 14 : return gmul(S, gexp(gneg(gmul(PiI2n(0, prec), A)), prec));
775 : }
776 :
777 : /* do not forget a minus sign for index 2 */
778 : static GEN
779 28 : integral_h12(GEN lin_grid, GEN selsm1, GEN s, GEN a, GEN lam, long prec)
780 : {
781 28 : GEN A, E, S = gen_0, ga = gaminus(gsubsg(1, s), prec);
782 28 : pari_sp av = avma;
783 28 : long j, l = lg(lin_grid);
784 :
785 28 : if (gequal0(ga)) return S;
786 28 : A = gaddgs(gmul2n(gadd(a,lam), 1), 1);
787 6776 : for (j = 1; j < l; j++)
788 : {
789 6748 : S = gadd(S, integrand_h12(selsm1, s, A, gel(lin_grid, j), prec));
790 6748 : if ((j & 0xff) == 0) S = gc_upto(av, S);
791 : }
792 28 : if (gequal0(S)) return gen_0;
793 28 : S = gmul(m_h(selsm1), S);
794 28 : E = gexp(gmul(PiI2n(0, prec), gmul(lam, gaddgs(lam, 1))), prec);
795 28 : return gmul(gmul(gmul(S, ga), E), gpow(Pi2n(1, prec), gsubgs(s, 1), prec));
796 : }
797 :
798 : struct _fun_q0_t { GEN sel, s, alam1, B; };
799 : static GEN
800 56 : _fun_q0(void *E, GEN x)
801 : {
802 56 : struct _fun_q0_t *S = (struct _fun_q0_t*)E;
803 56 : GEN z = integrand_h0l(S->sel, S->s, S->alam1, x, DEFAULTPREC);
804 56 : return addrlogabs(S->B, z);
805 : }
806 : static GEN
807 112 : _fun_q12(void *E, GEN x)
808 : {
809 112 : struct _fun_q0_t *S = (struct _fun_q0_t*)E;
810 112 : GEN z = integrand_h12(S->sel, S->s, S->alam1, x, DEFAULTPREC);
811 112 : return addrlogabs(S->B, z);
812 : }
813 :
814 : static GEN
815 14 : RZLERinit(GEN s, GEN a, GEN lam, GEN al, GEN numpoles, long prec)
816 : {
817 : GEN eps, r0, r1, r2, h, lin_grid, q, q0, q1, q2, sel0, sel1, sel2, lq;
818 14 : GEN pinv = ginv(PiI2(prec)), c = gmul2n(gadd(a, lam), -1), n0, n1, n2, c2;
819 : long m;
820 : struct _fun_q0_t E;
821 :
822 14 : if (!al || gequal0(al))
823 0 : al = gcmpgs(gabs(imag_i(s), prec), 100) < 0 ? ginv(stoi(4)) : gen_1;
824 14 : c2 = gsub(gsqr(c), gmul(s, pinv));
825 14 : r0 = gadd(c, gsqrt(c2, prec));
826 14 : r1 = gsqrt(gadd(c2, pinv), prec);
827 14 : r2 = gsub(r1, c);
828 14 : r1 = gneg(gadd(r1, c));
829 14 : n0 = gfloor(gsub(gadd(real_i(r0), imag_i(r0)), a));
830 14 : n1 = gneg(gfloor(gadd(gsub(real_i(r1), imag_i(r1)), real_i(lam))));
831 14 : n2 = gfloor(gadd(gsub(real_i(r2), imag_i(r2)), real_i(lam)));
832 :
833 14 : E.s = s; E.alam1 = gaddgs(gmul2n(gadd(a, lam), 1), 1);
834 14 : E.B = mulur(prec, mplog2(prec));
835 14 : lq = gmul(al, sqrtr_abs(mulrr(divsr(prec, Pi2n(1, DEFAULTPREC)),
836 : mplog2(DEFAULTPREC))));
837 14 : eps = gexp(PiI2n(-2, prec), prec);
838 14 : E.sel = sel0 = mkvec5(n0, r0, gdiv(al, eps), a, gen_0);
839 14 : q0 = findq(&E, &_fun_q0, lq, prec);
840 :
841 14 : if (!gequal1(al)) lq = gdiv(lq, gsqr(al));
842 14 : E.sel = sel1 = mkvec5(n1, r1, gmul(al, eps), lam, gen_0);
843 14 : q1 = findq(&E, &_fun_q12, lq, prec);
844 14 : E.sel = sel2 = mkvec5(n2, r2, gmul(al, eps), lam, gen_0);
845 14 : q2 = findq(&E, &_fun_q12, lq, prec);
846 14 : q = vecmax(mkvec3(q0, q1, q2)); m = get_m(q, prec);
847 14 : gel(sel0, 5) = gel(sel1, 5) = gel(sel2, 5) = h = divru(q, (m-1) >> 1);
848 14 : lin_grid = setlin_grid_exp(h, m, prec);
849 14 : if (!numpoles) numpoles = gen_1;
850 14 : return mkvec5(sel0, sel1, sel2, lin_grid, numpoles);
851 : }
852 :
853 42 : static GEN add3(GEN x, GEN y, GEN z) { return gadd(x, gadd(y,z)); }
854 14 : static GEN addsub(GEN x, GEN y, GEN z) { return gadd(x, gsub(y,z)); }
855 :
856 : static GEN
857 14 : lerch_ours(GEN sel, GEN s, GEN a, GEN lam, long prec)
858 : {
859 14 : GEN selsm0 = gel(sel, 1), selsm1 = gel(sel, 2), selsm2 = gel(sel, 3);
860 14 : GEN lin_grid = gel(sel, 4), v0, v1, v2;
861 14 : long numpoles = itos(gel(sel, 5));
862 14 : GEN worker = snm_closure(is_entry("_serh_worker"),
863 : mkvec4(NULL, NULL, NULL, NULL));
864 14 : v0 = add3(series_h0l(worker, m_n0(selsm0), s, a, lam, prec),
865 : series_residues_h0l(numpoles, selsm0, s, a, lam, prec),
866 : integral_h0l(lin_grid, selsm0, s, a, lam, prec));
867 14 : v1 = add3(series_h1(worker, m_n0(selsm1), s, a, lam, prec),
868 : series_residues_h1(numpoles, selsm1, s, a, lam, prec),
869 : integral_h12(lin_grid, selsm1, s, a, lam, prec));
870 14 : v2 = addsub(series_h2(worker, m_n0(selsm2), s, a, lam, prec),
871 : series_residues_h2(numpoles, selsm2, s, a, lam, prec),
872 : integral_h12(lin_grid, selsm2, s, a, lam, prec));
873 14 : return add3(v0, v1, v2);
874 : }
875 :
876 : static GEN
877 0 : RZlerch_easy(GEN s, GEN a, GEN lam, long prec)
878 : {
879 0 : pari_sp av = avma;
880 : GEN z, y, N;
881 0 : N = gdiv(gmulsg(prec + 5, mplog2(LOWDEFAULTPREC)),
882 : gmul(Pi2n(1, LOWDEFAULTPREC), imag_i(lam)));
883 0 : if (gexpo(N) > 40) pari_err_IMPL("precision too large in lerchzeta");
884 0 : N = gceil(N); prec += EXTRAPREC64;
885 0 : z = typ(lam) == t_INT ? gen_1 : gexp(gmul(PiI2(prec), lam), prec);
886 0 : y = parsum(gen_0, N, snm_closure(is_entry("_serh_worker"),
887 : mkvec4(z, a, gneg(s), stoi(prec))));
888 0 : return gc_GEN(av, gprec_wtrunc(y, prec));
889 : }
890 :
891 : static GEN
892 14 : mygfrac(GEN z)
893 0 : { return typ(z) == t_COMPLEX ? mkcomplex(gfrac(real_i(z)), imag_i(z))
894 14 : : gfrac(z); }
895 :
896 : static GEN
897 42 : lerchlarge(GEN s, GEN a, GEN lam, GEN al, GEN numpoles, long prec)
898 : {
899 42 : pari_sp av = avma;
900 42 : GEN val, sel, imlam = imag_i(lam);
901 : long prec2;
902 42 : switch(gsigne(imlam))
903 : {
904 0 : case -1: pari_err_IMPL("imag(lam) < 0");
905 0 : case 1: if (gexpo(imlam) >= -16) return RZlerch_easy(s, a, lam, prec);
906 : }
907 42 : if (gcmpgs(real_i(a), 1) < 0)
908 : {
909 14 : GEN P = gexp(gmul(PiI2(prec), lam), prec);
910 14 : GEN L = lerchlarge(s, gaddgs(a, 1), lam, al, numpoles, prec);
911 14 : return gc_upto(av, gadd(gpow(a, gneg(s), prec), gmul(P, L)));
912 : }
913 28 : if (gcmpgs(real_i(a), 2) >= 0)
914 : {
915 0 : GEN L, P = gexp(gneg(gmul(PiI2(prec), lam)), prec);
916 0 : a = gsubgs(a, 1); L = lerchlarge(s, a, lam, al, numpoles, prec);
917 0 : return gc_upto(av, gmul(P, gsub(L, gpow(a, gneg(s), prec))));
918 : }
919 28 : if (gsigne(imag_i(s)) > 0)
920 : {
921 : GEN L;
922 14 : lam = mygfrac(gneg(gconj(lam)));
923 14 : L = lerchlarge(gconj(s), a, lam, al, numpoles, prec);
924 14 : return gc_upto(av, gconj(L));
925 : }
926 14 : prec2 = prec + EXTRAPREC64;
927 14 : a = gprec_w(a, prec2);
928 14 : s = gprec_w(s, prec2);
929 14 : lam = gprec_w(lam, prec2);
930 14 : sel = RZLERinit(s, a, lam, al, numpoles, prec2);
931 14 : val = lerch_ours(sel, s, a, lam, prec2);
932 14 : return gc_GEN(av, gprec_wtrunc(val, prec));
933 : }
934 :
935 : GEN
936 7 : zetahurwitzlarge(GEN s, GEN a, long prec)
937 7 : { return lerchlarge(s, a, gen_0, gen_1, gen_1, prec); }
938 :
939 : GEN
940 7 : lerchzetalarge(GEN s, GEN a, GEN lam, long prec)
941 7 : { return lerchlarge(s, a, lam, gen_1, gen_1, prec); }
|
