mirror of
https://github.com/vale981/arb
synced 2025-03-05 09:21:38 -05:00
Merge pull request #334 from p15-git-acc/nan-to-inf
Treat infinite inputs more carefully.
This commit is contained in:
Fredrik Johansson
GitHub
commit
a976078288
9 changed files with 434 additions and 12 deletions
15
arb/cot_pi.c
15
arb/cot_pi.c
|
|
@ -13,6 +13,20 @@
|
|||
|
||||
void
|
||||
arb_cot_pi(arb_t y, const arb_t x, slong prec)
|
||||
{
|
||||
|
||||
if (!arb_is_finite(x))
|
||||
{
|
||||
arb_indeterminate(y);
|
||||
}
|
||||
else if (arb_is_int_2exp_si(x, -1))
|
||||
{
|
||||
if (arb_is_int(x))
|
||||
arb_indeterminate(y);
|
||||
else
|
||||
arb_zero(y);
|
||||
}
|
||||
else
|
||||
{
|
||||
arb_t u;
|
||||
arb_init(u);
|
||||
|
|
@ -20,4 +34,5 @@ arb_cot_pi(arb_t y, const arb_t x, slong prec)
|
|||
arb_div(y, u, y, prec);
|
||||
arb_clear(u);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
|
|
@ -67,12 +67,24 @@ arb_log_hypot(arb_t res, const arb_t a, const arb_t b, slong prec)
|
|||
}
|
||||
|
||||
if (!arb_is_finite(a) || !arb_is_finite(b))
|
||||
{
|
||||
if (arf_is_nan(arb_midref(a)) || arf_is_nan(arb_midref(b)))
|
||||
{
|
||||
arb_indeterminate(res);
|
||||
}
|
||||
else if ((!arb_is_finite(a) && !arb_contains_zero(a)) ||
|
||||
(!arb_is_finite(b) && !arb_contains_zero(b)))
|
||||
{
|
||||
arb_pos_inf(res);
|
||||
}
|
||||
else
|
||||
{
|
||||
arb_indeterminate(res);
|
||||
}
|
||||
return;
|
||||
}
|
||||
|
||||
/* a close to 1 -- for accurate acb_log1p */
|
||||
/* a close to 1 -- for accurate arb_log1p */
|
||||
if (mag_cmp_2exp_si(arb_radref(a), -3) < 0 &&
|
||||
mag_cmp_2exp_si(arb_radref(b), -3) < 0 &&
|
||||
arf_cmpabs_2exp_si(arb_midref(b), -3) < 0 &&
|
||||
|
|
|
|||
|
|
@ -17,6 +17,17 @@ arb_sinc_pi(arb_t res, const arb_t x, slong prec)
|
|||
mag_t m;
|
||||
arb_t t;
|
||||
|
||||
if (!arb_is_finite(x))
|
||||
{
|
||||
if (arf_is_nan(arb_midref(x)))
|
||||
arb_indeterminate(res);
|
||||
else if (arb_contains_zero(x))
|
||||
arb_zero_pm_one(res);
|
||||
else
|
||||
arb_zero(res);
|
||||
return;
|
||||
}
|
||||
|
||||
if (arb_is_zero(x))
|
||||
{
|
||||
arb_one(res);
|
||||
|
|
|
|||
14
arb/tan_pi.c
14
arb/tan_pi.c
|
|
@ -13,6 +13,19 @@
|
|||
|
||||
void
|
||||
arb_tan_pi(arb_t y, const arb_t x, slong prec)
|
||||
{
|
||||
if (!arb_is_finite(x))
|
||||
{
|
||||
arb_indeterminate(y);
|
||||
}
|
||||
else if (arb_is_int_2exp_si(x, -1))
|
||||
{
|
||||
if (arb_is_int(x))
|
||||
arb_zero(y);
|
||||
else
|
||||
arb_indeterminate(y);
|
||||
}
|
||||
else
|
||||
{
|
||||
arb_t u;
|
||||
arb_init(u);
|
||||
|
|
@ -20,4 +33,5 @@ arb_tan_pi(arb_t y, const arb_t x, slong prec)
|
|||
arb_div(y, y, u, prec);
|
||||
arb_clear(u);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
131
arb/test/t-cot_pi.c
Normal file
131
arb/test/t-cot_pi.c
Normal file
|
|
@ -0,0 +1,131 @@
|
|||
/*
|
||||
Copyright (C) 2013 Fredrik Johansson
|
||||
|
||||
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 "arb.h"
|
||||
|
||||
int main()
|
||||
{
|
||||
slong iter;
|
||||
flint_rand_t state;
|
||||
|
||||
flint_printf("cot_pi....");
|
||||
fflush(stdout);
|
||||
|
||||
flint_randinit(state);
|
||||
|
||||
for (iter = 0; iter < 10000 * arb_test_multiplier(); iter++)
|
||||
{
|
||||
arb_t x, y, a, b, c;
|
||||
slong prec1, prec2;
|
||||
|
||||
prec1 = 2 + n_randint(state, 1000);
|
||||
prec2 = prec1 + 30;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
arb_init(a);
|
||||
arb_init(b);
|
||||
arb_init(c);
|
||||
|
||||
arb_randtest(x, state, 1 + n_randint(state, 1000), 2 + n_randint(state, 100));
|
||||
|
||||
arb_cot_pi(a, x, prec1);
|
||||
arb_cot_pi(b, x, prec2);
|
||||
|
||||
/* check consistency */
|
||||
if (!arb_overlaps(a, b))
|
||||
{
|
||||
flint_printf("FAIL: overlap\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("b = "); arb_print(b); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
/* compare with cot */
|
||||
arb_const_pi(c, prec1);
|
||||
arb_mul(y, x, c, prec1);
|
||||
arb_cot(c, y, prec1);
|
||||
|
||||
if (!arb_overlaps(a, c))
|
||||
{
|
||||
flint_printf("FAIL: functional equation\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_print(y); flint_printf("\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("c = "); arb_print(c); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_cot_pi(x, x, prec1);
|
||||
|
||||
if (!arb_overlaps(a, x))
|
||||
{
|
||||
flint_printf("FAIL: aliasing\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_clear(x);
|
||||
arb_clear(y);
|
||||
arb_clear(a);
|
||||
arb_clear(b);
|
||||
arb_clear(c);
|
||||
}
|
||||
|
||||
/* Check a few special points. */
|
||||
{
|
||||
arb_t x, y;
|
||||
slong i, prec;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
|
||||
prec = 32;
|
||||
|
||||
for (i = -1; i <= 1; i++)
|
||||
{
|
||||
int arbitrary_integer = 7;
|
||||
|
||||
/* integer */
|
||||
arb_set_d(x, arbitrary_integer*i);
|
||||
arb_cot_pi(y, x, prec);
|
||||
if (arb_is_finite(y))
|
||||
{
|
||||
flint_printf("FAIL: (integer)\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
/* integer + 1/2 */
|
||||
arb_set_d(x, arbitrary_integer*i + 0.5);
|
||||
arb_cot_pi(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: (integer + 1/2)\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
}
|
||||
|
||||
arb_clear(x);
|
||||
arb_clear(y);
|
||||
}
|
||||
|
||||
flint_randclear(state);
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
||||
|
|
@ -94,6 +94,25 @@ int main()
|
|||
flint_abort();
|
||||
}
|
||||
|
||||
if (!arf_is_nan(arb_midref(a)) && !arf_is_nan(arb_midref(b)))
|
||||
{
|
||||
if ((!arb_is_finite(a) && !arb_contains_zero(a)) ||
|
||||
(!arb_is_finite(b) && !arb_contains_zero(b)))
|
||||
{
|
||||
if (arf_is_nan(arb_midref(r)) || arb_is_finite(r) ||
|
||||
!arb_is_positive(r))
|
||||
{
|
||||
flint_printf("FAIL: infinite\n\n");
|
||||
flint_printf("prec1 = %wd, acc1 = %wd, acc2 = %wd\n\n", prec1, acc1, acc2);
|
||||
flint_printf("a = "); arb_printn(a, 50, 0); flint_printf("\n\n");
|
||||
flint_printf("b = "); arb_printn(b, 50, 0); flint_printf("\n\n");
|
||||
flint_printf("r = "); arb_printn(r, 50, 0); flint_printf("\n\n");
|
||||
flint_printf("s = "); arb_printn(s, 50, 0); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
arb_clear(a);
|
||||
arb_clear(b);
|
||||
arb_clear(r);
|
||||
|
|
|
|||
|
|
@ -95,6 +95,53 @@ int main()
|
|||
arb_clear(a);
|
||||
}
|
||||
|
||||
/* Check a few special intervals. */
|
||||
{
|
||||
arb_t x, y;
|
||||
slong prec;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
|
||||
prec = 32;
|
||||
|
||||
arb_neg_inf(x);
|
||||
arb_sinc(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: negative infinity\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_pos_inf(x);
|
||||
arb_sinc(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: positive infinity\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_zero_pm_inf(x);
|
||||
arb_sinc(y, x, prec);
|
||||
if (!(arb_is_finite(y) &&
|
||||
arb_contains_negative(y) &&
|
||||
arb_contains_positive(y) &&
|
||||
arb_contains_zero(y)))
|
||||
{
|
||||
flint_printf("FAIL: the whole extended real line\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_clear(x);
|
||||
arb_clear(y);
|
||||
}
|
||||
|
||||
flint_randclear(state);
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
|
|
|
|||
|
|
@ -61,6 +61,53 @@ int main()
|
|||
arb_clear(z);
|
||||
}
|
||||
|
||||
/* Check a few special intervals. */
|
||||
{
|
||||
arb_t x, y;
|
||||
slong prec;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
|
||||
prec = 32;
|
||||
|
||||
arb_neg_inf(x);
|
||||
arb_sinc_pi(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: negative infinity\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_pos_inf(x);
|
||||
arb_sinc_pi(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: positive infinity\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_zero_pm_inf(x);
|
||||
arb_sinc_pi(y, x, prec);
|
||||
if (!(arb_is_finite(y) &&
|
||||
arb_contains_negative(y) &&
|
||||
arb_contains_positive(y) &&
|
||||
arb_contains_zero(y)))
|
||||
{
|
||||
flint_printf("FAIL: the whole extended real line\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_clear(x);
|
||||
arb_clear(y);
|
||||
}
|
||||
|
||||
flint_randclear(state);
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
|
|
|
|||
126
arb/test/t-tan_pi.c
Normal file
126
arb/test/t-tan_pi.c
Normal file
|
|
@ -0,0 +1,126 @@
|
|||
/*
|
||||
Copyright (C) 2013 Fredrik Johansson
|
||||
|
||||
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 "arb.h"
|
||||
|
||||
int main()
|
||||
{
|
||||
slong iter;
|
||||
flint_rand_t state;
|
||||
|
||||
flint_printf("tan_pi....");
|
||||
fflush(stdout);
|
||||
|
||||
flint_randinit(state);
|
||||
|
||||
for (iter = 0; iter < 10000 * arb_test_multiplier(); iter++)
|
||||
{
|
||||
arb_t x, y, a, b, c;
|
||||
slong prec1, prec2;
|
||||
|
||||
prec1 = 2 + n_randint(state, 1000);
|
||||
prec2 = prec1 + 30;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
arb_init(a);
|
||||
arb_init(b);
|
||||
arb_init(c);
|
||||
|
||||
arb_randtest(x, state, 1 + n_randint(state, 1000), 2 + n_randint(state, 100));
|
||||
|
||||
arb_tan_pi(a, x, prec1);
|
||||
arb_tan_pi(b, x, prec2);
|
||||
|
||||
/* check consistency */
|
||||
if (!arb_overlaps(a, b))
|
||||
{
|
||||
flint_printf("FAIL: overlap\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("b = "); arb_print(b); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
/* compare with tan */
|
||||
arb_const_pi(c, prec1);
|
||||
arb_mul(y, x, c, prec1);
|
||||
arb_tan(c, y, prec1);
|
||||
|
||||
if (!arb_overlaps(a, c))
|
||||
{
|
||||
flint_printf("FAIL: functional equation\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_print(y); flint_printf("\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("c = "); arb_print(c); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_tan_pi(x, x, prec1);
|
||||
|
||||
if (!arb_overlaps(a, x))
|
||||
{
|
||||
flint_printf("FAIL: aliasing\n\n");
|
||||
flint_printf("a = "); arb_print(a); flint_printf("\n\n");
|
||||
flint_printf("x = "); arb_print(x); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_clear(x);
|
||||
arb_clear(y);
|
||||
arb_clear(a);
|
||||
arb_clear(b);
|
||||
arb_clear(c);
|
||||
}
|
||||
|
||||
/* Check a few special points. */
|
||||
{
|
||||
arb_t x, y;
|
||||
slong i, prec;
|
||||
|
||||
arb_init(x);
|
||||
arb_init(y);
|
||||
|
||||
prec = 32;
|
||||
|
||||
for (i = -1; i <= 1; i++)
|
||||
{
|
||||
int arbitrary_integer = 7;
|
||||
|
||||
arb_set_d(x, arbitrary_integer*i);
|
||||
arb_tan_pi(y, x, prec);
|
||||
if (!arb_is_zero(y))
|
||||
{
|
||||
flint_printf("FAIL: (integer)\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
|
||||
arb_set_d(x, arbitrary_integer*i + 0.5);
|
||||
arb_tan_pi(y, x, prec);
|
||||
if (arb_is_finite(y))
|
||||
{
|
||||
flint_printf("FAIL: (integer + 1/2)\n");
|
||||
flint_printf("x = "); arb_printd(x, 30); flint_printf("\n\n");
|
||||
flint_printf("y = "); arb_printd(y, 30); flint_printf("\n\n");
|
||||
flint_abort();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
flint_randclear(state);
|
||||
flint_cleanup();
|
||||
flint_printf("PASS\n");
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
||||
Loading…
Add table
Reference in a new issue