add some fmprb test code, and fix a couple of bugs

This commit is contained in:
Fredrik Johansson 2012-09-10 12:27:32 +02:00
parent da38e48fe5
commit 621c0fc031
16 changed files with 943 additions and 4 deletions

37
fmprb.h
View file

@ -105,6 +105,8 @@ fmprb_set(fmprb_t x, const fmprb_t y)
fmpr_set(fmprb_radref(x), fmprb_radref(y));
}
void fmprb_set_round(fmprb_t z, const fmprb_t x, long prec);
static __inline__ void
fmprb_neg(fmprb_t x, const fmprb_t y)
{
@ -269,7 +271,42 @@ fmprb_set_fmpq(fmprb_t y, const fmpq_t x, long prec)
int fmprb_contains_fmpr(const fmprb_t x, const fmpr_t y);
int fmprb_contains_fmpq(const fmprb_t x, const fmpq_t y);
int fmprb_contains_fmpz(const fmprb_t x, const fmpz_t y);
int fmprb_contains_mpfr(const fmprb_t x, const mpfr_t y);
int fmprb_contains_zero(const fmprb_t x);
static __inline__ long
fmprb_rel_error_bits(const fmprb_t x)
{
fmpz_t midmag, radmag;
long result;
if (fmpr_is_zero(fmprb_radref(x)))
return -FMPR_PREC_EXACT;
if (fmpr_is_special(fmprb_midref(x)) || fmpr_is_special(fmprb_radref(x)))
return FMPR_PREC_EXACT;
fmpz_init(midmag);
fmpz_init(radmag);
fmpz_add_ui(midmag, fmpr_expref(fmprb_midref(x)),
fmpz_bits(fmpr_manref(fmprb_midref(x))));
fmpz_add_ui(radmag, fmpr_expref(fmprb_radref(x)),
fmpz_bits(fmpr_manref(fmprb_radref(x))));
fmpz_add_ui(radmag, radmag, 1);
result = _fmpz_sub_small(radmag, midmag);
fmpz_clear(midmag);
fmpz_clear(radmag);
return result;
}
static __inline__ long
fmprb_rel_accuracy_bits(const fmprb_t x)
{
return -fmprb_rel_error_bits(x);
}
#endif

View file

@ -42,5 +42,5 @@ fmprb_const_pi(fmprb_t x, long prec)
fmprb_const_pi_cached_prec = prec;
}
fmprb_set(x, fmprb_const_pi_cache);
fmprb_set_round(x, fmprb_const_pi_cache, prec);
}

View file

@ -35,7 +35,19 @@ zeta3_bsplit(fmprb_t P, fmprb_t Q, fmprb_t T, long a, long b, long wp, int cont)
fmprb_ui_pow_ui(Q, 2*a + 1, 5, wp);
fmprb_mul_ui(Q, Q, a ? 32 : 64, wp);
fmprb_set_ui(T, 205*a*a + 250*a + 77); /* xxx: fix overflow */
/* 205 * a^2 + 250 * a + 77 */
{
fmpz_t t;
fmpz_init(t);
fmpz_set_ui(t, a);
fmpz_mul_ui(t, t, 205);
fmpz_add_ui(t, t, 250);
fmpz_mul_ui(t, t, a);
fmpz_add_ui(t, t, 77);
fmprb_set_fmpz(T, t);
fmpz_clear(t);
}
if (a % 2)
fmprb_neg(T, T);
fmprb_mul(T, T, P, wp);

38
fmprb/contains_mpfr.c Normal file
View file

@ -0,0 +1,38 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int
fmprb_contains_mpfr(const fmprb_t x, const mpfr_t y)
{
int ans;
fmpr_t t;
fmpr_init(t);
fmpr_set_mpfr(t, y);
ans = fmprb_contains_fmpr(x, t);
fmpr_clear(t);
return ans;
}

38
fmprb/set_round.c Normal file
View file

@ -0,0 +1,38 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
void
fmprb_set_round(fmprb_t z, const fmprb_t x, long prec)
{
long r;
r = fmpr_set_round(fmprb_midref(z), fmprb_midref(x), prec, FMPR_RND_DOWN);
fmpr_add_error_result(fmprb_radref(z), fmprb_radref(x),
fmprb_midref(z), r, FMPRB_RAD_PREC, FMPR_RND_UP);
fmprb_adjust(z);
}

View file

@ -0,0 +1,78 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("const_euler_brent_mcmillan....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 250; iter++)
{
fmprb_t r;
mpfr_t s;
long accuracy, prec;
prec = 2 + n_randint(state, 1 << n_randint(state, 16));
fmprb_init(r);
mpfr_init2(s, prec + 1000);
fmprb_const_euler_brent_mcmillan(r, prec);
mpfr_const_euler(s, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: poor accuracy\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
mpfr_free_cache();
printf("PASS\n");
return EXIT_SUCCESS;
}

78
fmprb/test/t-const_pi.c Normal file
View file

@ -0,0 +1,78 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("const_pi....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 250; iter++)
{
fmprb_t r;
mpfr_t s;
long accuracy, prec;
prec = 2 + n_randint(state, 1 << n_randint(state, 18));
fmprb_init(r);
mpfr_init2(s, prec + 1000);
fmprb_const_pi(r, prec);
mpfr_const_pi(s, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: poor accuracy\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
mpfr_free_cache();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,78 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("const_pi_chudnovsky....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 250; iter++)
{
fmprb_t r;
mpfr_t s;
long accuracy, prec;
prec = 2 + n_randint(state, 1 << n_randint(state, 18));
fmprb_init(r);
mpfr_init2(s, prec + 1000);
fmprb_const_pi_chudnovsky(r, prec);
mpfr_const_pi(s, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: poor accuracy\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
mpfr_free_cache();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,78 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("const_zeta3_bsplit....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 250; iter++)
{
fmprb_t r;
mpfr_t s;
long accuracy, prec;
prec = 2 + n_randint(state, 1 << n_randint(state, 17));
fmprb_init(r);
mpfr_init2(s, prec + 1000);
fmprb_const_zeta3_bsplit(r, prec);
mpfr_zeta_ui(s, 3, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: poor accuracy\n\n");
printf("prec = %ld\n", prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
mpfr_free_cache();
printf("PASS\n");
return EXIT_SUCCESS;
}

81
fmprb/test/t-zeta_ui.c Normal file
View file

@ -0,0 +1,81 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 1000; iter++)
{
fmprb_t r;
ulong n;
mpfr_t s;
long prec, accuracy;
prec = 2 + n_randint(state, 1 << n_randint(state, 14));
fmprb_init(r);
mpfr_init2(s, prec + 100);
do { n = n_randint(state, 1 << n_randint(state, 10)); } while (n == 1);
fmprb_zeta_ui(r, n, prec);
mpfr_zeta_ui(s, n, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,82 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui_bernoulli....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 1000; iter++)
{
fmprb_t r;
ulong n;
mpfr_t s;
long prec, accuracy;
prec = 2 + n_randint(state, 1 << n_randint(state, 14));
fmprb_init(r);
mpfr_init2(s, prec + 100);
do { n = n_randint(state, 1 << n_randint(state, 10)); }
while (n % 2 || n == 0);
fmprb_zeta_ui_bernoulli(r, n, prec);
mpfr_zeta_ui(s, n, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,81 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui_bsplit....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 1000; iter++)
{
fmprb_t r;
ulong n;
mpfr_t s;
long prec, accuracy;
prec = 2 + n_randint(state, 1 << n_randint(state, 14));
fmprb_init(r);
mpfr_init2(s, prec + 100);
do { n = n_randint(state, 1 << n_randint(state, 10)); } while (n == 1);
fmprb_zeta_ui_bsplit(r, n, prec);
mpfr_zeta_ui(s, n, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,81 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui_euler_product....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 1000; iter++)
{
fmprb_t r;
ulong n;
mpfr_t s;
long prec, accuracy;
do { n = n_randint(state, 1 << n_randint(state, 10)); } while (n < 6);
prec = 2 + n_randint(state, 12 * n);
fmprb_init(r);
mpfr_init2(s, prec + 100);
fmprb_zeta_ui_euler_product(r, n, prec);
mpfr_zeta_ui(s, n, MPFR_RNDN);
if (!fmprb_contains_mpfr(r, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("r = "); fmprb_printd(r, prec / 3.33); printf("\n\n");
abort();
}
fmprb_clear(r);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,88 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui_vec....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 100; iter++)
{
fmprb_struct * r;
ulong n;
long i, num;
mpfr_t s;
long prec, accuracy;
prec = 2 + n_randint(state, 1 << n_randint(state, 13));
num = 1 + n_randint(state, 20);
r = _fmprb_vec_init(num);
mpfr_init2(s, prec + 100);
do { n = n_randint(state, 1 << n_randint(state, 10)); } while (n < 2);
fmprb_zeta_ui_vec(r, n, num, prec);
for (i = 0; i < num; i++)
{
mpfr_zeta_ui(s, n + i, MPFR_RNDN);
if (!fmprb_contains_mpfr(r + i, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n + i);
printf("r = "); fmprb_printd(r + i, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r + i);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("n = %lu\n\n", n + i);
printf("r = "); fmprb_printd(r + i, prec / 3.33); printf("\n\n");
abort();
}
}
_fmprb_vec_clear(r, num);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -0,0 +1,89 @@
/*=============================================================================
This file is part of ARB.
ARB is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
ARB is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with ARB; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
=============================================================================*/
/******************************************************************************
Copyright (C) 2012 Fredrik Johansson
******************************************************************************/
#include "fmprb.h"
int main()
{
long iter;
flint_rand_t state;
printf("zeta_ui_vec_borwein....");
fflush(stdout);
flint_randinit(state);
for (iter = 0; iter < 100; iter++)
{
fmprb_struct * r;
ulong n;
long i, num, step;
mpfr_t s;
long prec, accuracy;
prec = 2 + n_randint(state, 1 << n_randint(state, 13));
num = 1 + n_randint(state, 20);
step = 1 + n_randint(state, 5);
r = _fmprb_vec_init(num);
mpfr_init2(s, prec + 100);
do { n = n_randint(state, 1 << n_randint(state, 10)); } while (n < 2);
fmprb_zeta_ui_vec_borwein(r, n, num, step, prec);
for (i = 0; i < num; i++)
{
mpfr_zeta_ui(s, n + i * step, MPFR_RNDN);
if (!fmprb_contains_mpfr(r + i, s))
{
printf("FAIL: containment\n\n");
printf("n = %lu\n\n", n + i * step);
printf("r = "); fmprb_printd(r + i, prec / 3.33); printf("\n\n");
printf("s = "); mpfr_printf("%.275Rf\n", s); printf("\n\n");
abort();
}
accuracy = fmprb_rel_accuracy_bits(r + i);
if (accuracy < prec - 4)
{
printf("FAIL: accuracy = %ld, prec = %ld\n\n", accuracy, prec);
printf("n = %lu\n\n", n + i * step);
printf("r = "); fmprb_printd(r + i, prec / 3.33); printf("\n\n");
abort();
}
}
_fmprb_vec_clear(r, num);
mpfr_clear(s);
}
flint_randclear(state);
_fmpz_cleanup();
printf("PASS\n");
return EXIT_SUCCESS;
}

View file

@ -60,7 +60,7 @@ fmprb_zeta_ui(fmprb_t x, ulong n, long prec)
fmprb_zeta_ui_bsplit(x, n, prec);
}
/* large n */
else if (prec > 20 && n > 0.4 * pow(prec, 0.8))
else if (prec > 20 && n > 6 && n > 0.4 * pow(prec, 0.8))
{
fmprb_zeta_ui_euler_product(x, n, prec);
}
@ -115,7 +115,7 @@ fmprb_zeta_ui_vec(fmprb_struct * x, ulong start, long num, long prec)
start_odd = start % 2;
fmprb_zeta_ui_vec_even(x, start + start_odd, num_even, prec);
fmprb_zeta_ui_vec_odd(x + num_even, start + !start_odd, num_even, prec);
fmprb_zeta_ui_vec_odd(x + num_even, start + !start_odd, num_odd, prec);
/* interleave */
tmp = flint_malloc(sizeof(fmprb_struct) * num);