X-Git-Url: http://review.tizen.org/git/?a=blobdiff_plain;f=src%2Ftan.c;h=2154471e86dc859a82d810ec6c329884da664259;hb=f34506723f150a03d5e6e389820bd369b89e30b0;hp=24cd92b7ce5a3e8383cf6943a347d8f9c34f996e;hpb=c9266275a7240e07da227b662cc6afcb0000d7a6;p=platform%2Fupstream%2Fmpc.git diff --git a/src/tan.c b/src/tan.c index 24cd92b..2154471 100644 --- a/src/tan.c +++ b/src/tan.c @@ -1,62 +1,61 @@ /* mpc_tan -- tangent of a complex number. -Copyright (C) 2008, 2009, 2010, 2011 INRIA +Copyright (C) 2008, 2009 Philippe Th\'eveny, Andreas Enge -This file is part of GNU MPC. +This file is part of the MPC Library. -GNU MPC is free software; you can redistribute it and/or modify it under -the terms of the GNU Lesser General Public License as published by the -Free Software Foundation; either version 3 of the License, or (at your +The MPC Library is free software; you can redistribute it and/or modify +it under the terms of the GNU Lesser General Public License as published by +the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version. -GNU MPC 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 Lesser General Public License for -more details. +The MPC Library 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 Lesser General Public +License for more details. You should have received a copy of the GNU Lesser General Public License -along with this program. If not, see http://www.gnu.org/licenses/ . -*/ +along with the MPC Library; see the file COPYING.LIB. If not, write to +the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, +MA 02111-1307, USA. */ -#include /* for MPC_ASSERT */ -#include #include "mpc-impl.h" int mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) { mpc_t x, y; - mpfr_prec_t prec; - mpfr_exp_t err; + mp_prec_t prec; + mp_exp_t err; int ok = 0; int inex; /* special values */ - if (!mpc_fin_p (op)) + if (!mpfr_number_p (MPC_RE (op)) || !mpfr_number_p (MPC_IM (op))) { - if (mpfr_nan_p (mpc_realref (op))) + if (mpfr_nan_p (MPC_RE (op))) { - if (mpfr_inf_p (mpc_imagref (op))) + if (mpfr_inf_p (MPC_IM (op))) /* tan(NaN -i*Inf) = +/-0 -i */ /* tan(NaN +i*Inf) = +/-0 +i */ { /* exact unless 1 is not in exponent range */ inex = mpc_set_si_si (rop, 0, - (MPFR_SIGN (mpc_imagref (op)) < 0) ? -1 : +1, + (MPFR_SIGN (MPC_IM (op)) < 0) ? -1 : +1, rnd); } else /* tan(NaN +i*y) = NaN +i*NaN, when y is finite */ /* tan(NaN +i*NaN) = NaN +i*NaN */ { - mpfr_set_nan (mpc_realref (rop)); - mpfr_set_nan (mpc_imagref (rop)); + mpfr_set_nan (MPC_RE (rop)); + mpfr_set_nan (MPC_IM (rop)); inex = MPC_INEX (0, 0); /* always exact */ } } - else if (mpfr_nan_p (mpc_imagref (op))) + else if (mpfr_nan_p (MPC_IM (op))) { - if (mpfr_cmp_ui (mpc_realref (op), 0) == 0) + if (mpfr_cmp_ui (MPC_RE (op), 0) == 0) /* tan(-0 +i*NaN) = -0 +i*NaN */ /* tan(+0 +i*NaN) = +0 +i*NaN */ { @@ -66,28 +65,28 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) else /* tan(x +i*NaN) = NaN +i*NaN, when x != 0 */ { - mpfr_set_nan (mpc_realref (rop)); - mpfr_set_nan (mpc_imagref (rop)); + mpfr_set_nan (MPC_RE (rop)); + mpfr_set_nan (MPC_IM (rop)); inex = MPC_INEX (0, 0); /* always exact */ } } - else if (mpfr_inf_p (mpc_realref (op))) + else if (mpfr_inf_p (MPC_RE (op))) { - if (mpfr_inf_p (mpc_imagref (op))) + if (mpfr_inf_p (MPC_IM (op))) /* tan(-Inf -i*Inf) = -/+0 -i */ /* tan(-Inf +i*Inf) = -/+0 +i */ /* tan(+Inf -i*Inf) = +/-0 -i */ /* tan(+Inf +i*Inf) = +/-0 +i */ { - const int sign_re = mpfr_signbit (mpc_realref (op)); + const int sign_re = mpfr_signbit (MPC_RE (op)); int inex_im; - mpfr_set_ui (mpc_realref (rop), 0, MPC_RND_RE (rnd)); - mpfr_setsign (mpc_realref (rop), mpc_realref (rop), sign_re, GMP_RNDN); + mpfr_set_ui (MPC_RE (rop), 0, MPC_RND_RE (rnd)); + mpfr_setsign (MPC_RE (rop), MPC_RE (rop), sign_re, GMP_RNDN); /* exact, unless 1 is not in exponent range */ - inex_im = mpfr_set_si (mpc_imagref (rop), - mpfr_signbit (mpc_imagref (op)) ? -1 : +1, + inex_im = mpfr_set_si (MPC_IM (rop), + mpfr_signbit (MPC_IM (op)) ? -1 : +1, MPC_RND_IM (rnd)); inex = MPC_INEX (0, inex_im); } @@ -95,8 +94,8 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) /* tan(-Inf +i*y) = tan(+Inf +i*y) = NaN +i*NaN, when y is finite */ { - mpfr_set_nan (mpc_realref (rop)); - mpfr_set_nan (mpc_imagref (rop)); + mpfr_set_nan (MPC_RE (rop)); + mpfr_set_nan (MPC_IM (rop)); inex = MPC_INEX (0, 0); /* always exact */ } } @@ -111,13 +110,13 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) mpfr_init (c); mpfr_init (s); - mpfr_sin_cos (s, c, mpc_realref (op), GMP_RNDN); - mpfr_set_ui (mpc_realref (rop), 0, MPC_RND_RE (rnd)); - mpfr_setsign (mpc_realref (rop), mpc_realref (rop), + mpfr_sin_cos (s, c, MPC_RE (op), GMP_RNDN); + mpfr_set_ui (MPC_RE (rop), 0, MPC_RND_RE (rnd)); + mpfr_setsign (MPC_RE (rop), MPC_RE (rop), mpfr_signbit (c) != mpfr_signbit (s), GMP_RNDN); /* exact, unless 1 is not in exponent range */ - inex_im = mpfr_set_si (mpc_imagref (rop), - (mpfr_signbit (mpc_imagref (op)) ? -1 : +1), + inex_im = mpfr_set_si (MPC_IM (rop), + (mpfr_signbit (MPC_IM (op)) ? -1 : +1), MPC_RND_IM (rnd)); inex = MPC_INEX (0, inex_im); @@ -128,26 +127,26 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) return inex; } - if (mpfr_zero_p (mpc_realref (op))) + if (mpfr_zero_p (MPC_RE (op))) /* tan(-0 -i*y) = -0 +i*tanh(y), when y is finite. */ /* tan(+0 +i*y) = +0 +i*tanh(y), when y is finite. */ { int inex_im; - mpfr_set (mpc_realref (rop), mpc_realref (op), MPC_RND_RE (rnd)); - inex_im = mpfr_tanh (mpc_imagref (rop), mpc_imagref (op), MPC_RND_IM (rnd)); + mpfr_set (MPC_RE (rop), MPC_RE (op), MPC_RND_RE (rnd)); + inex_im = mpfr_tanh (MPC_IM (rop), MPC_IM (op), MPC_RND_IM (rnd)); return MPC_INEX (0, inex_im); } - if (mpfr_zero_p (mpc_imagref (op))) + if (mpfr_zero_p (MPC_IM (op))) /* tan(x -i*0) = tan(x) -i*0, when x is finite. */ /* tan(x +i*0) = tan(x) +i*0, when x is finite. */ { int inex_re; - inex_re = mpfr_tan (mpc_realref (rop), mpc_realref (op), MPC_RND_RE (rnd)); - mpfr_set (mpc_imagref (rop), mpc_imagref (op), MPC_RND_IM (rnd)); + inex_re = mpfr_tan (MPC_RE (rop), MPC_RE (op), MPC_RND_RE (rnd)); + mpfr_set (MPC_IM (rop), MPC_IM (op), MPC_RND_IM (rnd)); return MPC_INEX (inex_re, 0); } @@ -181,7 +180,8 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) do { - mpfr_exp_t k, exr, eyr, eyi, ezr; + mp_exp_t k; + mp_exp_t exr, eyr, eyi, ezr; ok = 0; @@ -193,45 +193,19 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) /* rounding away from zero: except in the cases x=0 or y=0 (processed above), sin x and cos y are never exact, so rounding away from 0 is rounding towards 0 and adding one ulp to the absolute value */ - mpc_sin_cos (x, y, op, MPC_RNDZZ, MPC_RNDZZ); - MPFR_ADD_ONE_ULP (mpc_realref (x)); - MPFR_ADD_ONE_ULP (mpc_imagref (x)); - MPFR_ADD_ONE_ULP (mpc_realref (y)); - MPFR_ADD_ONE_ULP (mpc_imagref (y)); - MPC_ASSERT (mpfr_zero_p (mpc_realref (x)) == 0); - - if ( mpfr_inf_p (mpc_realref (x)) || mpfr_inf_p (mpc_imagref (x)) - || mpfr_inf_p (mpc_realref (y)) || mpfr_inf_p (mpc_imagref (y))) { - /* If the real or imaginary part of x is infinite, it means that - Im(op) was large, in which case the result is - sign(tan(Re(op)))*0 + sign(Im(op))*I, - where sign(tan(Re(op))) = sign(Re(x))*sign(Re(y)). */ - int inex_re, inex_im; - mpfr_set_ui (mpc_realref (rop), 0, GMP_RNDN); - if (mpfr_sgn (mpc_realref (x)) * mpfr_sgn (mpc_realref (y)) < 0) - { - mpfr_neg (mpc_realref (rop), mpc_realref (rop), GMP_RNDN); - inex_re = 1; - } - else - inex_re = -1; /* +0 is rounded down */ - if (mpfr_sgn (mpc_imagref (op)) > 0) - { - mpfr_set_ui (mpc_imagref (rop), 1, GMP_RNDN); - inex_im = 1; - } - else - { - mpfr_set_si (mpc_imagref (rop), -1, GMP_RNDN); - inex_im = -1; - } - inex = MPC_INEX(inex_re, inex_im); - goto end; - } - - exr = mpfr_get_exp (mpc_realref (x)); - eyr = mpfr_get_exp (mpc_realref (y)); - eyi = mpfr_get_exp (mpc_imagref (y)); + mpc_sin (x, op, MPC_RNDZZ); + mpfr_signbit (MPC_RE (x)) ? + mpfr_nextbelow (MPC_RE (x)) : mpfr_nextabove (MPC_RE (x)); + mpfr_signbit (MPC_IM (x)) ? + mpfr_nextbelow (MPC_IM (x)) : mpfr_nextabove (MPC_IM (x)); + exr = mpfr_get_exp (MPC_RE (x)); + mpc_cos (y, op, MPC_RNDZZ); + mpfr_signbit (MPC_RE (y)) ? + mpfr_nextbelow (MPC_RE (y)) : mpfr_nextabove (MPC_RE (y)); + mpfr_signbit (MPC_IM (y)) ? + mpfr_nextbelow (MPC_IM (y)) : mpfr_nextabove (MPC_IM (y)); + eyr = mpfr_get_exp (MPC_RE (y)); + eyi = mpfr_get_exp (MPC_IM (y)); /* some parts of the quotient may be exact */ inex = mpc_div (x, x, y, MPC_RNDZZ); @@ -241,17 +215,18 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) tan(1+14*I) = 1.26e-10 + 1.00*I. For small precision sin(op) and cos(op) differ only by a factor I, thus after mpc_div x = I and its real part is zero. */ - if (mpfr_zero_p (mpc_realref (x)) || mpfr_zero_p (mpc_imagref (x))) + if (mpfr_zero_p (MPC_RE (x)) || mpfr_zero_p (MPC_IM (x))) { err = prec; /* double precision */ continue; } if (MPC_INEX_RE (inex)) - MPFR_ADD_ONE_ULP (mpc_realref (x)); + mpfr_signbit (MPC_RE (x)) ? + mpfr_nextbelow (MPC_RE (x)) : mpfr_nextabove (MPC_RE (x)); if (MPC_INEX_IM (inex)) - MPFR_ADD_ONE_ULP (mpc_imagref (x)); - MPC_ASSERT (mpfr_zero_p (mpc_realref (x)) == 0); - ezr = mpfr_get_exp (mpc_realref (x)); + mpfr_signbit (MPC_IM (x)) ? + mpfr_nextbelow (MPC_IM (x)) : mpfr_nextabove (MPC_IM (x)); + ezr = mpfr_get_exp (MPC_RE (x)); /* FIXME: compute k = Exp(Re(x))+Exp(Re(y))-2min{Exp(Re(y)), Exp(Im(y))}-Exp(Re(x/y)) @@ -259,24 +234,23 @@ mpc_tan (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) k = exr - ezr + MPC_MAX(-eyr, eyr - 2 * eyi); err = k < 2 ? 7 : (k == 2 ? 8 : (5 + k)); - /* Can the real part be rounded? */ - ok = (!mpfr_number_p (mpc_realref (x))) - || mpfr_can_round (mpc_realref(x), prec - err, GMP_RNDN, GMP_RNDZ, - MPC_PREC_RE(rop) + (MPC_RND_RE(rnd) == GMP_RNDN)); + /* Can the real part be rounded ? */ + ok = mpfr_inf_p (MPC_RE (x)) + || mpfr_can_round (MPC_RE(x), prec - err, GMP_RNDN, GMP_RNDZ, + MPFR_PREC(MPC_RE(rop)) + (MPC_RND_RE(rnd) == GMP_RNDN)); if (ok) { - /* Can the imaginary part be rounded? */ - ok = (!mpfr_number_p (mpc_imagref (x))) - || mpfr_can_round (mpc_imagref(x), prec - 6, GMP_RNDN, GMP_RNDZ, - MPC_PREC_IM(rop) + (MPC_RND_IM(rnd) == GMP_RNDN)); + /* Can the imaginary part be rounded ? */ + ok = mpfr_inf_p (MPC_IM (x)) + || mpfr_can_round (MPC_IM(x), prec - 6, GMP_RNDN, GMP_RNDZ, + MPFR_PREC(MPC_IM(rop)) + (MPC_RND_IM(rnd) == GMP_RNDN)); } } while (ok == 0); inex = mpc_set (rop, x, rnd); - end: mpc_clear (x); mpc_clear (y);