mirror of
https://github.com/vale981/arb
synced 2025-03-05 09:21:38 -05:00
acb_dirichlet_l: use functional equation
This commit is contained in:
Fredrik Johansson
parent
0be71f1294
commit
ec7d9783e6
6 changed files with 270 additions and 116 deletions
|
|
@ -130,7 +130,7 @@ void acb_dirichlet_hardy_theta(acb_ptr res, const acb_t t,
|
|||
const dirichlet_group_t G, const dirichlet_char_t chi,
|
||||
slong len, slong prec);
|
||||
|
||||
void acb_dirichlet_hardy_z(acb_t res, const acb_t t,
|
||||
void acb_dirichlet_hardy_z(acb_ptr res, const acb_t t,
|
||||
const dirichlet_group_t G, const dirichlet_char_t chi,
|
||||
slong len, slong prec);
|
||||
|
||||
|
|
|
|||
|
|
@ -13,7 +13,7 @@
|
|||
#include "acb_poly.h"
|
||||
|
||||
void
|
||||
acb_dirichlet_hardy_z(acb_t res, const acb_t t,
|
||||
acb_dirichlet_hardy_z(acb_ptr res, const acb_t t,
|
||||
const dirichlet_group_t G, const dirichlet_char_t chi,
|
||||
slong len, slong prec)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -13,7 +13,7 @@
|
|||
#include "acb_dirichlet.h"
|
||||
|
||||
void
|
||||
acb_dirichlet_l(acb_t res, const acb_t s,
|
||||
acb_dirichlet_l_general(acb_t res, const acb_t s,
|
||||
const dirichlet_group_t G, const dirichlet_char_t chi, slong prec)
|
||||
{
|
||||
/* this cutoff is probably too conservative when q is large */
|
||||
|
|
@ -27,3 +27,81 @@ acb_dirichlet_l(acb_t res, const acb_t s,
|
|||
}
|
||||
}
|
||||
|
||||
void
|
||||
acb_dirichlet_l(acb_t res, const acb_t s,
|
||||
const dirichlet_group_t G, const dirichlet_char_t chi, slong prec)
|
||||
{
|
||||
if (!acb_is_finite(s))
|
||||
{
|
||||
acb_indeterminate(res);
|
||||
}
|
||||
else if (G->q == 1)
|
||||
{
|
||||
acb_dirichlet_zeta(res, s, prec);
|
||||
}
|
||||
else if (dirichlet_conductor_char(G, chi) == G->q &&
|
||||
(arf_cmp_d(arb_midref(acb_realref(s)), -0.5) < 0 ||
|
||||
(G->q != 1 && dirichlet_parity_char(G, chi) == 0 &&
|
||||
arf_cmpabs_d(arb_midref(acb_imagref(s)), 0.125) < 0 &&
|
||||
arf_cmp_d(arb_midref(acb_realref(s)), 0.125) < 0)))
|
||||
{
|
||||
/* use functional equation */
|
||||
acb_t t, u, v;
|
||||
int parity;
|
||||
ulong q;
|
||||
|
||||
parity = dirichlet_parity_char(G, chi);
|
||||
q = G->q;
|
||||
|
||||
acb_init(t);
|
||||
acb_init(u);
|
||||
acb_init(v);
|
||||
|
||||
/* gamma((1-s+p)/2) / gamma((s+p)/2) */
|
||||
acb_add_ui(t, s, parity, prec);
|
||||
acb_mul_2exp_si(t, t, -1);
|
||||
acb_rgamma(t, t, prec);
|
||||
|
||||
if (!acb_is_zero(t)) /* assumes q != 1 when s = 0 */
|
||||
{
|
||||
acb_neg(u, s);
|
||||
acb_add_ui(u, u, 1 + parity, prec);
|
||||
acb_mul_2exp_si(u, u, -1);
|
||||
acb_gamma(u, u, prec);
|
||||
acb_mul(t, t, u, prec);
|
||||
|
||||
/* epsilon */
|
||||
acb_dirichlet_root_number(u, G, chi, prec);
|
||||
acb_mul(t, t, u, prec);
|
||||
|
||||
/* (pi/q)^(s-1/2) */
|
||||
acb_const_pi(u, prec);
|
||||
acb_div_ui(u, u, q, prec);
|
||||
acb_set_d(v, -0.5);
|
||||
acb_add(v, v, s, prec);
|
||||
acb_pow(u, u, v, prec);
|
||||
acb_mul(t, t, u, prec);
|
||||
|
||||
acb_sub_ui(u, s, 1, prec);
|
||||
acb_neg(u, u);
|
||||
acb_conj(u, u);
|
||||
acb_dirichlet_l_general(u, u, G, chi, prec);
|
||||
acb_conj(u, u);
|
||||
acb_mul(t, t, u, prec);
|
||||
|
||||
if (dirichlet_char_is_real(G, chi) && acb_is_real(s))
|
||||
arb_zero(acb_imagref(t));
|
||||
}
|
||||
|
||||
acb_set(res, t);
|
||||
|
||||
acb_clear(t);
|
||||
acb_clear(u);
|
||||
acb_clear(v);
|
||||
}
|
||||
else
|
||||
{
|
||||
acb_dirichlet_l_general(res, s, G, chi, prec);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
|
|
@ -1,5 +1,5 @@
|
|||
/*
|
||||
Copyright (C) 2016 Pascal Molin
|
||||
Copyright (C) 2016 Fredrik Johansson
|
||||
|
||||
This file is part of Arb.
|
||||
|
||||
|
|
@ -11,135 +11,64 @@
|
|||
|
||||
#include "acb_dirichlet.h"
|
||||
|
||||
#define nq 5
|
||||
#define nx 3
|
||||
|
||||
int main()
|
||||
{
|
||||
|
||||
slong i, j;
|
||||
ulong q[nq] = {3, 5, 61, 91, 800};
|
||||
ulong m[nq] = {2, 4, 11, 2, 3};
|
||||
slong prec = 150;
|
||||
|
||||
acb_ptr x;
|
||||
/* cannot test at s = 1 with hurwitz */
|
||||
const char * x_r[nx] = { "1", "0.5", "0.5" };
|
||||
const char * x_i[nx] = { "1", "0", "6" };
|
||||
|
||||
acb_t ref, res;
|
||||
/*
|
||||
default(realprecision, 54)
|
||||
X = [ 1 + I, 1/2, 1/2 + 6 * I ]
|
||||
C = [Mod(2,3),Mod(4,5),Mod(11,61),Mod(2,91),Mod(3,800)]
|
||||
v = concat([ [lfun(c,x) | x<-X] | c<-C])
|
||||
apply(z->printf("\"%s\",\n",real(z)),v)
|
||||
apply(z->printf("\"%s\",\n",imag(z)),v)
|
||||
*/
|
||||
const char * ref_r[nq * nx] = {
|
||||
"0.655527984002548033786648216345221087939439503905627469",
|
||||
"0.480867557696828626181220063235589877776829730832371863",
|
||||
"1.56831301727577320813799211138797101541772722814204757",
|
||||
"0.521271244517346991221550773660594765406476858135844321",
|
||||
"0.231750947504015755883383661760877226427886696409005898",
|
||||
"0.275543455389521803395512886745330595086898302178508437",
|
||||
"0.598221809458540554839300433683735304093606595684903281",
|
||||
"0.489264190003695740292779374874163221990017067040417393",
|
||||
"0.573331076412428980263984182365365715292560207445592018",
|
||||
"0.510279695870740409778738767334484809708615155539404548",
|
||||
"0.635626509594367380604827545000418331455019188562281349",
|
||||
"0.129304857274642475564179442785425797926079767522671163",
|
||||
"1.18088858810025653590356481638012816019876881487868657",
|
||||
"2.17175778983760437737667738885959799183430688287297767",
|
||||
"3.41568550810774629867945639900431994221065497147578087"
|
||||
};
|
||||
const char * ref_i[nq * nx] = {
|
||||
"0.220206044893215842652155131408935133596486560067476474",
|
||||
"0",
|
||||
"-0.969458654385732175077973304161399773236957587792986099",
|
||||
"0.354614573267731219242838516784605303288232150760467423",
|
||||
"0",
|
||||
"-0.995392028773643947872231871832838543767598301887892796",
|
||||
"1.04370497561090171487193145841005574472705644411957863",
|
||||
"-0.108902811943905225853677097712717212629591264759957602",
|
||||
"-0.232114369998608907337769019848201890558327186146689311",
|
||||
"-0.133300066189980774635445078240315148544665020358019145",
|
||||
"0.0119464572932630291870372694406253796888930803905106876",
|
||||
"-0.567660589679294457801153713636532209809112025502518666",
|
||||
"-0.654079942571300523223917775358845549990877148918886474",
|
||||
"0.970337207245832214408380510793679653538607483205616894",
|
||||
"-1.43652482351673593824956935036654893593947145947637807"
|
||||
};
|
||||
slong iter;
|
||||
flint_rand_t state;
|
||||
|
||||
flint_printf("l....");
|
||||
fflush(stdout);
|
||||
|
||||
x = _acb_vec_init(nx);
|
||||
flint_randinit(state);
|
||||
|
||||
for (j = 0; j < nx; j++)
|
||||
for (iter = 0; iter < 1000 * arb_test_multiplier(); iter++)
|
||||
{
|
||||
if (arb_set_str(acb_realref(x + j), x_r[j], prec) ||
|
||||
arb_set_str(acb_imagref(x + j), x_i[j], prec) )
|
||||
acb_t s, t, u;
|
||||
dirichlet_group_t G;
|
||||
dirichlet_char_t chi;
|
||||
ulong q, k;
|
||||
slong prec;
|
||||
|
||||
acb_init(s);
|
||||
acb_init(t);
|
||||
acb_init(u);
|
||||
|
||||
q = 1 + n_randint(state, 50);
|
||||
prec = 2 + n_randint(state, 100);
|
||||
k = n_randint(state, n_euler_phi(q));
|
||||
|
||||
dirichlet_group_init(G, q);
|
||||
dirichlet_char_init(chi, G);
|
||||
dirichlet_char_index(chi, G, k);
|
||||
|
||||
if (n_randint(state, 2))
|
||||
acb_set_si(s, n_randint(state, 50) - 25);
|
||||
else
|
||||
acb_randtest(s, state, 2 + n_randint(state, 200), 2);
|
||||
|
||||
acb_dirichlet_l_hurwitz(t, s, G, chi, prec);
|
||||
acb_dirichlet_l(u, s, G, chi, prec);
|
||||
|
||||
if (!acb_overlaps(t, u))
|
||||
{
|
||||
flint_printf("error while setting x[%ld] <- %s+I*%s\n",
|
||||
j, x_r[j], x_i[j]);
|
||||
flint_printf("FAIL: overlap\n\n");
|
||||
flint_printf("iter = %ld q = %lu k = %lu prec = %ld\n\n", iter, q, k, prec);
|
||||
flint_printf("s = "); acb_printn(s, 100, 0); flint_printf("\n\n");
|
||||
flint_printf("t = "); acb_printn(t, 100, 0); flint_printf("\n\n");
|
||||
flint_printf("u = "); acb_printn(u, 100, 0); flint_printf("\n\n");
|
||||
abort();
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
acb_init(ref);
|
||||
acb_init(res);
|
||||
|
||||
for (i = 0; i < nq; i++)
|
||||
{
|
||||
dirichlet_group_t G;
|
||||
dirichlet_char_t chi;
|
||||
|
||||
dirichlet_group_init(G, q[i]);
|
||||
dirichlet_char_init(chi, G);
|
||||
|
||||
dirichlet_char_log(chi, G, m[i]);
|
||||
|
||||
for (j = 0; j < nx; j++)
|
||||
{
|
||||
|
||||
if (arb_set_str(acb_realref(ref), ref_r[i * nx + j], prec - 10) ||
|
||||
arb_set_str(acb_imagref(ref), ref_i[i * nx + j], prec - 10) )
|
||||
{
|
||||
flint_printf("error while setting ref <- %s+I*%s\n",
|
||||
ref_r[i * nx + j], ref_i[i * nx + j]);
|
||||
abort();
|
||||
}
|
||||
|
||||
acb_dirichlet_l_hurwitz(res, x + j, G, chi, prec + 10);
|
||||
|
||||
if (!acb_contains(ref, res))
|
||||
{
|
||||
flint_printf("FAIL:\n\n");
|
||||
flint_printf("q = %wu\n", q[i]);
|
||||
flint_printf("m = %wu\n", m[i]);
|
||||
flint_printf("x = ");
|
||||
acb_printd(x, 54);
|
||||
flint_printf("\nref = ");
|
||||
acb_printd(ref, 54);
|
||||
flint_printf("\nl(chi,x) = ");
|
||||
acb_printd(res, 54);
|
||||
flint_printf("\n\n");
|
||||
abort();
|
||||
}
|
||||
|
||||
}
|
||||
dirichlet_char_clear(chi);
|
||||
dirichlet_group_clear(G);
|
||||
|
||||
acb_clear(s);
|
||||
acb_clear(t);
|
||||
acb_clear(u);
|
||||
}
|
||||
acb_clear(ref);
|
||||
acb_clear(res);
|
||||
_acb_vec_clear(x, nx);
|
||||
|
||||
|
||||
flint_randclear(state);
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
||||
|
|
|
|||
145
acb_dirichlet/test/t-l_hurwitz.c
Normal file
145
acb_dirichlet/test/t-l_hurwitz.c
Normal file
|
|
@ -0,0 +1,145 @@
|
|||
/*
|
||||
Copyright (C) 2016 Pascal Molin
|
||||
|
||||
This file is part of Arb.
|
||||
|
||||
Arb is free software: you can redistribute it and/or modify it under
|
||||
the terms of the GNU Lesser General Public License (LGPL) as published
|
||||
by the Free Software Foundation; either version 2.1 of the License, or
|
||||
(at your option) any later version. See <http://www.gnu.org/licenses/>.
|
||||
*/
|
||||
|
||||
#include "acb_dirichlet.h"
|
||||
|
||||
#define nq 5
|
||||
#define nx 3
|
||||
|
||||
int main()
|
||||
{
|
||||
|
||||
slong i, j;
|
||||
ulong q[nq] = {3, 5, 61, 91, 800};
|
||||
ulong m[nq] = {2, 4, 11, 2, 3};
|
||||
slong prec = 150;
|
||||
|
||||
acb_ptr x;
|
||||
/* cannot test at s = 1 with hurwitz */
|
||||
const char * x_r[nx] = { "1", "0.5", "0.5" };
|
||||
const char * x_i[nx] = { "1", "0", "6" };
|
||||
|
||||
acb_t ref, res;
|
||||
/*
|
||||
default(realprecision, 54)
|
||||
X = [ 1 + I, 1/2, 1/2 + 6 * I ]
|
||||
C = [Mod(2,3),Mod(4,5),Mod(11,61),Mod(2,91),Mod(3,800)]
|
||||
v = concat([ [lfun(c,x) | x<-X] | c<-C])
|
||||
apply(z->printf("\"%s\",\n",real(z)),v)
|
||||
apply(z->printf("\"%s\",\n",imag(z)),v)
|
||||
*/
|
||||
const char * ref_r[nq * nx] = {
|
||||
"0.655527984002548033786648216345221087939439503905627469",
|
||||
"0.480867557696828626181220063235589877776829730832371863",
|
||||
"1.56831301727577320813799211138797101541772722814204757",
|
||||
"0.521271244517346991221550773660594765406476858135844321",
|
||||
"0.231750947504015755883383661760877226427886696409005898",
|
||||
"0.275543455389521803395512886745330595086898302178508437",
|
||||
"0.598221809458540554839300433683735304093606595684903281",
|
||||
"0.489264190003695740292779374874163221990017067040417393",
|
||||
"0.573331076412428980263984182365365715292560207445592018",
|
||||
"0.510279695870740409778738767334484809708615155539404548",
|
||||
"0.635626509594367380604827545000418331455019188562281349",
|
||||
"0.129304857274642475564179442785425797926079767522671163",
|
||||
"1.18088858810025653590356481638012816019876881487868657",
|
||||
"2.17175778983760437737667738885959799183430688287297767",
|
||||
"3.41568550810774629867945639900431994221065497147578087"
|
||||
};
|
||||
const char * ref_i[nq * nx] = {
|
||||
"0.220206044893215842652155131408935133596486560067476474",
|
||||
"0",
|
||||
"-0.969458654385732175077973304161399773236957587792986099",
|
||||
"0.354614573267731219242838516784605303288232150760467423",
|
||||
"0",
|
||||
"-0.995392028773643947872231871832838543767598301887892796",
|
||||
"1.04370497561090171487193145841005574472705644411957863",
|
||||
"-0.108902811943905225853677097712717212629591264759957602",
|
||||
"-0.232114369998608907337769019848201890558327186146689311",
|
||||
"-0.133300066189980774635445078240315148544665020358019145",
|
||||
"0.0119464572932630291870372694406253796888930803905106876",
|
||||
"-0.567660589679294457801153713636532209809112025502518666",
|
||||
"-0.654079942571300523223917775358845549990877148918886474",
|
||||
"0.970337207245832214408380510793679653538607483205616894",
|
||||
"-1.43652482351673593824956935036654893593947145947637807"
|
||||
};
|
||||
|
||||
flint_printf("l....");
|
||||
fflush(stdout);
|
||||
|
||||
x = _acb_vec_init(nx);
|
||||
|
||||
for (j = 0; j < nx; j++)
|
||||
{
|
||||
if (arb_set_str(acb_realref(x + j), x_r[j], prec) ||
|
||||
arb_set_str(acb_imagref(x + j), x_i[j], prec) )
|
||||
{
|
||||
flint_printf("error while setting x[%ld] <- %s+I*%s\n",
|
||||
j, x_r[j], x_i[j]);
|
||||
abort();
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
acb_init(ref);
|
||||
acb_init(res);
|
||||
|
||||
for (i = 0; i < nq; i++)
|
||||
{
|
||||
dirichlet_group_t G;
|
||||
dirichlet_char_t chi;
|
||||
|
||||
dirichlet_group_init(G, q[i]);
|
||||
dirichlet_char_init(chi, G);
|
||||
|
||||
dirichlet_char_log(chi, G, m[i]);
|
||||
|
||||
for (j = 0; j < nx; j++)
|
||||
{
|
||||
|
||||
if (arb_set_str(acb_realref(ref), ref_r[i * nx + j], prec - 10) ||
|
||||
arb_set_str(acb_imagref(ref), ref_i[i * nx + j], prec - 10) )
|
||||
{
|
||||
flint_printf("error while setting ref <- %s+I*%s\n",
|
||||
ref_r[i * nx + j], ref_i[i * nx + j]);
|
||||
abort();
|
||||
}
|
||||
|
||||
acb_dirichlet_l_hurwitz(res, x + j, G, chi, prec + 10);
|
||||
|
||||
if (!acb_contains(ref, res))
|
||||
{
|
||||
flint_printf("FAIL:\n\n");
|
||||
flint_printf("q = %wu\n", q[i]);
|
||||
flint_printf("m = %wu\n", m[i]);
|
||||
flint_printf("x = ");
|
||||
acb_printd(x, 54);
|
||||
flint_printf("\nref = ");
|
||||
acb_printd(ref, 54);
|
||||
flint_printf("\nl(chi,x) = ");
|
||||
acb_printd(res, 54);
|
||||
flint_printf("\n\n");
|
||||
abort();
|
||||
}
|
||||
|
||||
}
|
||||
dirichlet_char_clear(chi);
|
||||
dirichlet_group_clear(G);
|
||||
|
||||
}
|
||||
acb_clear(ref);
|
||||
acb_clear(res);
|
||||
_acb_vec_clear(x, nx);
|
||||
|
||||
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
|
@ -13,6 +13,8 @@
|
|||
|
||||
void acb_zeta_si(acb_t z, slong s, slong prec);
|
||||
|
||||
/* todo: use euler product for complex s, and check efficiency
|
||||
for large negative integers */
|
||||
void
|
||||
acb_dirichlet_zeta(acb_t res, const acb_t s, slong prec)
|
||||
{
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue