diff options
Diffstat (limited to 'gsl-1.9/specfunc/trig.c')
-rw-r--r-- | gsl-1.9/specfunc/trig.c | 771 |
1 files changed, 771 insertions, 0 deletions
diff --git a/gsl-1.9/specfunc/trig.c b/gsl-1.9/specfunc/trig.c new file mode 100644 index 0000000..a604cd8 --- /dev/null +++ b/gsl-1.9/specfunc/trig.c @@ -0,0 +1,771 @@ +/* specfunc/trig.c + * + * Copyright (C) 1996, 1997, 1998, 1999, 2000 Gerard Jungman + * + * This program 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. + * + * This program 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 this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + */ + +/* Author: G. Jungman */ + +#include <config.h> +#include <gsl/gsl_math.h> +#include <gsl/gsl_errno.h> +#include <gsl/gsl_sf_log.h> +#include <gsl/gsl_sf_trig.h> + +#include "error.h" + +#include "chebyshev.h" +#include "cheb_eval.c" + +/* sinh(x) series + * double-precision for |x| < 1.0 + */ +inline +static +int +sinh_series(const double x, double * result) +{ + const double y = x*x; + const double c0 = 1.0/6.0; + const double c1 = 1.0/120.0; + const double c2 = 1.0/5040.0; + const double c3 = 1.0/362880.0; + const double c4 = 1.0/39916800.0; + const double c5 = 1.0/6227020800.0; + const double c6 = 1.0/1307674368000.0; + const double c7 = 1.0/355687428096000.0; + *result = x*(1.0 + y*(c0+y*(c1+y*(c2+y*(c3+y*(c4+y*(c5+y*(c6+y*c7)))))))); + return GSL_SUCCESS; +} + + +/* cosh(x)-1 series + * double-precision for |x| < 1.0 + */ +inline +static +int +cosh_m1_series(const double x, double * result) +{ + const double y = x*x; + const double c0 = 0.5; + const double c1 = 1.0/24.0; + const double c2 = 1.0/720.0; + const double c3 = 1.0/40320.0; + const double c4 = 1.0/3628800.0; + const double c5 = 1.0/479001600.0; + const double c6 = 1.0/87178291200.0; + const double c7 = 1.0/20922789888000.0; + const double c8 = 1.0/6402373705728000.0; + *result = y*(c0+y*(c1+y*(c2+y*(c3+y*(c4+y*(c5+y*(c6+y*(c7+y*c8)))))))); + return GSL_SUCCESS; +} + + +/* Chebyshev expansion for f(t) = sinc((t+1)/2), -1 < t < 1 + */ +static double sinc_data[17] = { + 1.133648177811747875422, + -0.532677564732557348781, + -0.068293048346633177859, + 0.033403684226353715020, + 0.001485679893925747818, + -0.000734421305768455295, + -0.000016837282388837229, + 0.000008359950146618018, + 0.000000117382095601192, + -0.000000058413665922724, + -0.000000000554763755743, + 0.000000000276434190426, + 0.000000000001895374892, + -0.000000000000945237101, + -0.000000000000004900690, + 0.000000000000002445383, + 0.000000000000000009925 +}; +static cheb_series sinc_cs = { + sinc_data, + 16, + -1, 1, + 10 +}; + + +/* Chebyshev expansion for f(t) = g((t+1)Pi/8), -1<t<1 + * g(x) = (sin(x)/x - 1)/(x*x) + */ +static double sin_data[12] = { + -0.3295190160663511504173, + 0.0025374284671667991990, + 0.0006261928782647355874, + -4.6495547521854042157541e-06, + -5.6917531549379706526677e-07, + 3.7283335140973803627866e-09, + 3.0267376484747473727186e-10, + -1.7400875016436622322022e-12, + -1.0554678305790849834462e-13, + 5.3701981409132410797062e-16, + 2.5984137983099020336115e-17, + -1.1821555255364833468288e-19 +}; +static cheb_series sin_cs = { + sin_data, + 11, + -1, 1, + 11 +}; + +/* Chebyshev expansion for f(t) = g((t+1)Pi/8), -1<t<1 + * g(x) = (2(cos(x) - 1)/(x^2) + 1) / x^2 + */ +static double cos_data[11] = { + 0.165391825637921473505668118136, + -0.00084852883845000173671196530195, + -0.000210086507222940730213625768083, + 1.16582269619760204299639757584e-6, + 1.43319375856259870334412701165e-7, + -7.4770883429007141617951330184e-10, + -6.0969994944584252706997438007e-11, + 2.90748249201909353949854872638e-13, + 1.77126739876261435667156490461e-14, + -7.6896421502815579078577263149e-17, + -3.7363121133079412079201377318e-18 +}; +static cheb_series cos_cs = { + cos_data, + 10, + -1, 1, + 10 +}; + + +/*-*-*-*-*-*-*-*-*-*-*-* Functions with Error Codes *-*-*-*-*-*-*-*-*-*-*-*/ + +/* I would have prefered just using the library sin() function. + * But after some experimentation I decided that there was + * no good way to understand the error; library sin() is just a black box. + * So we have to roll our own. + */ +int +gsl_sf_sin_e(double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + { + const double P1 = 7.85398125648498535156e-1; + const double P2 = 3.77489470793079817668e-8; + const double P3 = 2.69515142907905952645e-15; + + const double sgn_x = GSL_SIGN(x); + const double abs_x = fabs(x); + + if(abs_x < GSL_ROOT4_DBL_EPSILON) { + const double x2 = x*x; + result->val = x * (1.0 - x2/6.0); + result->err = fabs(x*x2*x2 / 100.0); + return GSL_SUCCESS; + } + else { + double sgn_result = sgn_x; + double y = floor(abs_x/(0.25*M_PI)); + int octant = y - ldexp(floor(ldexp(y,-3)),3); + int stat_cs; + double z; + + if(GSL_IS_ODD(octant)) { + octant += 1; + octant &= 07; + y += 1.0; + } + + if(octant > 3) { + octant -= 4; + sgn_result = -sgn_result; + } + + z = ((abs_x - y * P1) - y * P2) - y * P3; + + if(octant == 0) { + gsl_sf_result sin_cs_result; + const double t = 8.0*fabs(z)/M_PI - 1.0; + stat_cs = cheb_eval_e(&sin_cs, t, &sin_cs_result); + result->val = z * (1.0 + z*z * sin_cs_result.val); + } + else { /* octant == 2 */ + gsl_sf_result cos_cs_result; + const double t = 8.0*fabs(z)/M_PI - 1.0; + stat_cs = cheb_eval_e(&cos_cs, t, &cos_cs_result); + result->val = 1.0 - 0.5*z*z * (1.0 - z*z * cos_cs_result.val); + } + + result->val *= sgn_result; + + if(abs_x > 1.0/GSL_DBL_EPSILON) { + result->err = fabs(result->val); + } + else if(abs_x > 100.0/GSL_SQRT_DBL_EPSILON) { + result->err = 2.0 * abs_x * GSL_DBL_EPSILON * fabs(result->val); + } + else if(abs_x > 0.1/GSL_SQRT_DBL_EPSILON) { + result->err = 2.0 * GSL_SQRT_DBL_EPSILON * fabs(result->val); + } + else { + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + } + + return stat_cs; + } + } +} + + +int +gsl_sf_cos_e(double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + { + const double P1 = 7.85398125648498535156e-1; + const double P2 = 3.77489470793079817668e-8; + const double P3 = 2.69515142907905952645e-15; + + const double abs_x = fabs(x); + + if(abs_x < GSL_ROOT4_DBL_EPSILON) { + const double x2 = x*x; + result->val = 1.0 - 0.5*x2; + result->err = fabs(x2*x2/12.0); + return GSL_SUCCESS; + } + else { + double sgn_result = 1.0; + double y = floor(abs_x/(0.25*M_PI)); + int octant = y - ldexp(floor(ldexp(y,-3)),3); + int stat_cs; + double z; + + if(GSL_IS_ODD(octant)) { + octant += 1; + octant &= 07; + y += 1.0; + } + + if(octant > 3) { + octant -= 4; + sgn_result = -sgn_result; + } + + if(octant > 1) { + sgn_result = -sgn_result; + } + + z = ((abs_x - y * P1) - y * P2) - y * P3; + + if(octant == 0) { + gsl_sf_result cos_cs_result; + const double t = 8.0*fabs(z)/M_PI - 1.0; + stat_cs = cheb_eval_e(&cos_cs, t, &cos_cs_result); + result->val = 1.0 - 0.5*z*z * (1.0 - z*z * cos_cs_result.val); + } + else { /* octant == 2 */ + gsl_sf_result sin_cs_result; + const double t = 8.0*fabs(z)/M_PI - 1.0; + stat_cs = cheb_eval_e(&sin_cs, t, &sin_cs_result); + result->val = z * (1.0 + z*z * sin_cs_result.val); + } + + result->val *= sgn_result; + + if(abs_x > 1.0/GSL_DBL_EPSILON) { + result->err = fabs(result->val); + } + else if(abs_x > 100.0/GSL_SQRT_DBL_EPSILON) { + result->err = 2.0 * abs_x * GSL_DBL_EPSILON * fabs(result->val); + } + else if(abs_x > 0.1/GSL_SQRT_DBL_EPSILON) { + result->err = 2.0 * GSL_SQRT_DBL_EPSILON * fabs(result->val); + } + else { + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + } + + return stat_cs; + } + } +} + + +int +gsl_sf_hypot_e(const double x, const double y, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + if(x == 0.0 && y == 0.0) { + result->val = 0.0; + result->err = 0.0; + return GSL_SUCCESS; + } + else { + const double a = fabs(x); + const double b = fabs(y); + const double min = GSL_MIN_DBL(a,b); + const double max = GSL_MAX_DBL(a,b); + const double rat = min/max; + const double root_term = sqrt(1.0 + rat*rat); + + if(max < GSL_DBL_MAX/root_term) { + result->val = max * root_term; + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else { + OVERFLOW_ERROR(result); + } + } +} + + +int +gsl_sf_complex_sin_e(const double zr, const double zi, + gsl_sf_result * szr, gsl_sf_result * szi) +{ + /* CHECK_POINTER(szr) */ + /* CHECK_POINTER(szi) */ + + if(fabs(zi) < 1.0) { + double ch_m1, sh; + sinh_series(zi, &sh); + cosh_m1_series(zi, &ch_m1); + szr->val = sin(zr)*(ch_m1 + 1.0); + szi->val = cos(zr)*sh; + szr->err = 2.0 * GSL_DBL_EPSILON * fabs(szr->val); + szi->err = 2.0 * GSL_DBL_EPSILON * fabs(szi->val); + return GSL_SUCCESS; + } + else if(fabs(zi) < GSL_LOG_DBL_MAX) { + double ex = exp(zi); + double ch = 0.5*(ex+1.0/ex); + double sh = 0.5*(ex-1.0/ex); + szr->val = sin(zr)*ch; + szi->val = cos(zr)*sh; + szr->err = 2.0 * GSL_DBL_EPSILON * fabs(szr->val); + szi->err = 2.0 * GSL_DBL_EPSILON * fabs(szi->val); + return GSL_SUCCESS; + } + else { + OVERFLOW_ERROR_2(szr, szi); + } +} + + +int +gsl_sf_complex_cos_e(const double zr, const double zi, + gsl_sf_result * czr, gsl_sf_result * czi) +{ + /* CHECK_POINTER(czr) */ + /* CHECK_POINTER(czi) */ + + if(fabs(zi) < 1.0) { + double ch_m1, sh; + sinh_series(zi, &sh); + cosh_m1_series(zi, &ch_m1); + czr->val = cos(zr)*(ch_m1 + 1.0); + czi->val = -sin(zr)*sh; + czr->err = 2.0 * GSL_DBL_EPSILON * fabs(czr->val); + czi->err = 2.0 * GSL_DBL_EPSILON * fabs(czi->val); + return GSL_SUCCESS; + } + else if(fabs(zi) < GSL_LOG_DBL_MAX) { + double ex = exp(zi); + double ch = 0.5*(ex+1.0/ex); + double sh = 0.5*(ex-1.0/ex); + czr->val = cos(zr)*ch; + czi->val = -sin(zr)*sh; + czr->err = 2.0 * GSL_DBL_EPSILON * fabs(czr->val); + czi->err = 2.0 * GSL_DBL_EPSILON * fabs(czi->val); + return GSL_SUCCESS; + } + else { + OVERFLOW_ERROR_2(czr,czi); + } +} + + +int +gsl_sf_complex_logsin_e(const double zr, const double zi, + gsl_sf_result * lszr, gsl_sf_result * lszi) +{ + /* CHECK_POINTER(lszr) */ + /* CHECK_POINTER(lszi) */ + + if(zi > 60.0) { + lszr->val = -M_LN2 + zi; + lszi->val = 0.5*M_PI - zr; + lszr->err = 2.0 * GSL_DBL_EPSILON * fabs(lszr->val); + lszi->err = 2.0 * GSL_DBL_EPSILON * fabs(lszi->val); + } + else if(zi < -60.0) { + lszr->val = -M_LN2 - zi; + lszi->val = -0.5*M_PI + zr; + lszr->err = 2.0 * GSL_DBL_EPSILON * fabs(lszr->val); + lszi->err = 2.0 * GSL_DBL_EPSILON * fabs(lszi->val); + } + else { + gsl_sf_result sin_r, sin_i; + int status; + gsl_sf_complex_sin_e(zr, zi, &sin_r, &sin_i); /* ok by construction */ + status = gsl_sf_complex_log_e(sin_r.val, sin_i.val, lszr, lszi); + if(status == GSL_EDOM) { + DOMAIN_ERROR_2(lszr, lszi); + } + } + return gsl_sf_angle_restrict_symm_e(&(lszi->val)); +} + + +int +gsl_sf_lnsinh_e(const double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + if(x <= 0.0) { + DOMAIN_ERROR(result); + } + else if(fabs(x) < 1.0) { + double eps; + sinh_series(x, &eps); + result->val = log(eps); + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else if(x < -0.5*GSL_LOG_DBL_EPSILON) { + result->val = x + log(0.5*(1.0 - exp(-2.0*x))); + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else { + result->val = -M_LN2 + x; + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } +} + + +int gsl_sf_lncosh_e(const double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + if(fabs(x) < 1.0) { + double eps; + cosh_m1_series(x, &eps); + return gsl_sf_log_1plusx_e(eps, result); + } + else if(x < -0.5*GSL_LOG_DBL_EPSILON) { + result->val = x + log(0.5*(1.0 + exp(-2.0*x))); + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else { + result->val = -M_LN2 + x; + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } +} + + +/* +inline int gsl_sf_sincos_e(const double theta, double * s, double * c) +{ + double tan_half = tan(0.5 * theta); + double den = 1. + tan_half*tan_half; + double cos_theta = (1.0 - tan_half*tan_half) / den; + double sin_theta = 2.0 * tan_half / den; +} +*/ + +int +gsl_sf_polar_to_rect(const double r, const double theta, + gsl_sf_result * x, gsl_sf_result * y) +{ + double t = theta; + int status = gsl_sf_angle_restrict_symm_e(&t); + double c = cos(t); + double s = sin(t); + x->val = r * cos(t); + y->val = r * sin(t); + x->err = r * fabs(s * GSL_DBL_EPSILON * t); + x->err += 2.0 * GSL_DBL_EPSILON * fabs(x->val); + y->err = r * fabs(c * GSL_DBL_EPSILON * t); + y->err += 2.0 * GSL_DBL_EPSILON * fabs(y->val); + return status; +} + + +int +gsl_sf_rect_to_polar(const double x, const double y, + gsl_sf_result * r, gsl_sf_result * theta) +{ + int stat_h = gsl_sf_hypot_e(x, y, r); + if(r->val > 0.0) { + theta->val = atan2(y, x); + theta->err = 2.0 * GSL_DBL_EPSILON * fabs(theta->val); + return stat_h; + } + else { + DOMAIN_ERROR(theta); + } +} + + +int gsl_sf_angle_restrict_symm_err_e(const double theta, gsl_sf_result * result) +{ + /* synthetic extended precision constants */ + const double P1 = 4 * 7.8539812564849853515625e-01; + const double P2 = 4 * 3.7748947079307981766760e-08; + const double P3 = 4 * 2.6951514290790594840552e-15; + const double TwoPi = 2*(P1 + P2 + P3); + + const double y = GSL_SIGN(theta) * 2 * floor(fabs(theta)/TwoPi); + double r = ((theta - y*P1) - y*P2) - y*P3; + + if(r > M_PI) { r = (((r-2*P1)-2*P2)-2*P3); } /* r-TwoPi */ + else if (r < -M_PI) r = (((r+2*P1)+2*P2)+2*P3); /* r+TwoPi */ + + result->val = r; + + if(fabs(theta) > 0.0625/GSL_DBL_EPSILON) { + result->val = GSL_NAN; + result->err = GSL_NAN; + GSL_ERROR ("error", GSL_ELOSS); + } + else if(fabs(theta) > 0.0625/GSL_SQRT_DBL_EPSILON) { + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val - theta); + return GSL_SUCCESS; + } + else { + double delta = fabs(result->val - theta); + result->err = 2.0 * GSL_DBL_EPSILON * ((delta < M_PI) ? delta : M_PI); + return GSL_SUCCESS; + } +} + + +int gsl_sf_angle_restrict_pos_err_e(const double theta, gsl_sf_result * result) +{ + /* synthetic extended precision constants */ + const double P1 = 4 * 7.85398125648498535156e-01; + const double P2 = 4 * 3.77489470793079817668e-08; + const double P3 = 4 * 2.69515142907905952645e-15; + const double TwoPi = 2*(P1 + P2 + P3); + + const double y = 2*floor(theta/TwoPi); + + double r = ((theta - y*P1) - y*P2) - y*P3; + + if(r > TwoPi) {r = (((r-2*P1)-2*P2)-2*P3); } /* r-TwoPi */ + else if (r < 0) { /* may happen due to FP rounding */ + r = (((r+2*P1)+2*P2)+2*P3); /* r+TwoPi */ + } + + result->val = r; + + if(fabs(theta) > 0.0625/GSL_DBL_EPSILON) { + result->val = GSL_NAN; + result->err = fabs(result->val); + GSL_ERROR ("error", GSL_ELOSS); + } + else if(fabs(theta) > 0.0625/GSL_SQRT_DBL_EPSILON) { + result->err = GSL_DBL_EPSILON * fabs(result->val - theta); + return GSL_SUCCESS; + } + else { + double delta = fabs(result->val - theta); + result->err = 2.0 * GSL_DBL_EPSILON * ((delta < M_PI) ? delta : M_PI); + return GSL_SUCCESS; + } +} + + +int gsl_sf_angle_restrict_symm_e(double * theta) +{ + gsl_sf_result r; + int stat = gsl_sf_angle_restrict_symm_err_e(*theta, &r); + *theta = r.val; + return stat; +} + + +int gsl_sf_angle_restrict_pos_e(double * theta) +{ + gsl_sf_result r; + int stat = gsl_sf_angle_restrict_pos_err_e(*theta, &r); + *theta = r.val; + return stat; +} + + +int gsl_sf_sin_err_e(const double x, const double dx, gsl_sf_result * result) +{ + int stat_s = gsl_sf_sin_e(x, result); + result->err += fabs(cos(x) * dx); + result->err += GSL_DBL_EPSILON * fabs(result->val); + return stat_s; +} + + +int gsl_sf_cos_err_e(const double x, const double dx, gsl_sf_result * result) +{ + int stat_c = gsl_sf_cos_e(x, result); + result->err += fabs(sin(x) * dx); + result->err += GSL_DBL_EPSILON * fabs(result->val); + return stat_c; +} + + +#if 0 +int +gsl_sf_sin_pi_x_e(const double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + if(-100.0 < x && x < 100.0) { + result->val = sin(M_PI * x) / (M_PI * x); + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else { + const double N = floor(x + 0.5); + const double f = x - N; + + if(N < INT_MAX && N > INT_MIN) { + /* Make it an integer if we can. Saves another + * call to floor(). + */ + const int intN = (int)N; + const double sign = ( GSL_IS_ODD(intN) ? -1.0 : 1.0 ); + result->val = sign * sin(M_PI * f); + result->err = GSL_DBL_EPSILON * fabs(result->val); + } + else if(N > 2.0/GSL_DBL_EPSILON || N < -2.0/GSL_DBL_EPSILON) { + /* All integer-valued floating point numbers + * bigger than 2/eps=2^53 are actually even. + */ + result->val = 0.0; + result->err = 0.0; + } + else { + const double resN = N - 2.0*floor(0.5*N); /* 0 for even N, 1 for odd N */ + const double sign = ( fabs(resN) > 0.5 ? -1.0 : 1.0 ); + result->val = sign * sin(M_PI*f); + result->err = GSL_DBL_EPSILON * fabs(result->val); + } + + return GSL_SUCCESS; + } +} +#endif + + +int gsl_sf_sinc_e(double x, gsl_sf_result * result) +{ + /* CHECK_POINTER(result) */ + + { + const double ax = fabs(x); + + if(ax < 0.8) { + /* Do not go to the limit of the fit since + * there is a zero there and the Chebyshev + * accuracy will go to zero. + */ + return cheb_eval_e(&sinc_cs, 2.0*ax-1.0, result); + } + else if(ax < 100.0) { + /* Small arguments are no problem. + * We trust the library sin() to + * roughly machine precision. + */ + result->val = sin(M_PI * ax)/(M_PI * ax); + result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return GSL_SUCCESS; + } + else { + /* Large arguments must be handled separately. + */ + const double r = M_PI*ax; + gsl_sf_result s; + int stat_s = gsl_sf_sin_e(r, &s); + result->val = s.val/r; + result->err = s.err/r + 2.0 * GSL_DBL_EPSILON * fabs(result->val); + return stat_s; + } + } +} + + + +/*-*-*-*-*-*-*-*-*-* Functions w/ Natural Prototypes *-*-*-*-*-*-*-*-*-*-*/ + +#include "eval.h" + +double gsl_sf_sin(const double x) +{ + EVAL_RESULT(gsl_sf_sin_e(x, &result)); +} + +double gsl_sf_cos(const double x) +{ + EVAL_RESULT(gsl_sf_cos_e(x, &result)); +} + +double gsl_sf_hypot(const double x, const double y) +{ + EVAL_RESULT(gsl_sf_hypot_e(x, y, &result)); +} + +double gsl_sf_lnsinh(const double x) +{ + EVAL_RESULT(gsl_sf_lnsinh_e(x, &result)); +} + +double gsl_sf_lncosh(const double x) +{ + EVAL_RESULT(gsl_sf_lncosh_e(x, &result)); +} + +double gsl_sf_angle_restrict_symm(const double theta) +{ + double result = theta; + EVAL_DOUBLE(gsl_sf_angle_restrict_symm_e(&result)); +} + +double gsl_sf_angle_restrict_pos(const double theta) +{ + double result = theta; + EVAL_DOUBLE(gsl_sf_angle_restrict_pos_e(&result)); +} + +#if 0 +double gsl_sf_sin_pi_x(const double x) +{ + EVAL_RESULT(gsl_sf_sin_pi_x_e(x, &result)); +} +#endif + +double gsl_sf_sinc(const double x) +{ + EVAL_RESULT(gsl_sf_sinc_e(x, &result)); +} |