[gnumeric] SINPI, COSPI: New functions.



commit 38b63891d1acd0ba3c7305c1c51d2224412991d2
Author: Morten Welinder <terra gnome org>
Date:   Thu Nov 14 20:16:11 2013 -0500

    SINPI, COSPI: New functions.
    
    These functions are quite handy and don't suffer argument reduction
    losses.  They are recommended in the 2008 version of IEEE 754.

 ChangeLog                     |    5 ++
 NEWS                          |    2 +
 doc/C/func.defs               |   14 +++++
 doc/C/functions.xml           |   54 ++++++++++++++++++++
 plugins/fn-math/ChangeLog     |    4 ++
 plugins/fn-math/functions.c   |   40 +++++++++++++++
 plugins/fn-math/plugin.xml.in |    2 +
 src/mathfunc.c                |  108 +++++++++++++++++++++++++++++++++++------
 src/mathfunc.h                |    3 +
 9 files changed, 217 insertions(+), 15 deletions(-)
---
diff --git a/ChangeLog b/ChangeLog
index f6cd360..0641339 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,8 @@
+2013-11-14  Morten Welinder  <terra gnome org>
+
+       * src/mathfunc.c (gnm_sinpi, gnm_cospi): New functions.
+       (bessel_i, etc, qfactf, lgamma_r): Use gnm_sinpi and gnm_cospi.
+
 2013-11-13  Morten Welinder  <terra gnome org>
 
        * src/mathfunc.c (lgamma_r): Fix fallback.  Didn't work for x<10.
diff --git a/NEWS b/NEWS
index 57395a7..9076c6a 100644
--- a/NEWS
+++ b/NEWS
@@ -17,7 +17,9 @@ Morten:
        * Restore sheet reordering by drag.
        * Fix BETA on win32.
        * Fix win32 registry initialization.
+       * Fix win32 gdk-pixbuf install heisen-crash.
        * Incorporate new tests from crlibm.
+       * New functions SINPI and COSPI.
 
 Xabier Rodríguez Calvar:
        * Fix dialog button order. [#710378]
diff --git a/doc/C/func.defs b/doc/C/func.defs
index 6226c46..d864533 100644
--- a/doc/C/func.defs
+++ b/doc/C/func.defs
@@ -3232,6 +3232,13 @@ The depreciation coefficient used is:
 @SEEALSO=SIN,TAN,SINH,COSH,TANH
 
 @CATEGORY=Mathematics
+ FUNCTION=COSPI
+ SHORTDESC=the cosine of Pi* {x}
+ SYNTAX=COSPI(x)
+ ARGUMENTDESCRIPTION=@{x}: number of half turns
+ SEEALSO=COS
+
+ CATEGORY=Mathematics
 @FUNCTION=COT
 @SHORTDESC=the cotangent of @{x}
 @SYNTAX=COT(x)
@@ -3722,6 +3729,13 @@ If @{d} is less than zero, @{x} is rounded away from 0 to the left of the decima
 @SEEALSO=SIN,COSH,ASINH
 
 @CATEGORY=Mathematics
+ FUNCTION=SINPI
+ SHORTDESC=the sine of Pi* {x}
+ SYNTAX=SINPI(x)
+ ARGUMENTDESCRIPTION=@{x}: number of half turns
+ SEEALSO=SIN
+
+ CATEGORY=Mathematics
 @FUNCTION=SQRT
 @SHORTDESC=square root of @{x}
 @SYNTAX=SQRT(x)
diff --git a/doc/C/functions.xml b/doc/C/functions.xml
index be4f783..c71b1fa 100644
--- a/doc/C/functions.xml
+++ b/doc/C/functions.xml
@@ -10784,6 +10784,33 @@
       </para>
       </refsect1>
     </refentry>
+    <refentry id="gnumeric-function-COSPI">
+      <refmeta>
+        <refentrytitle>
+          <function>COSPI</function>
+        </refentrytitle>
+      </refmeta>
+      <refnamediv>
+        <refname>
+          <function>COSPI</function>
+        </refname>
+        <refpurpose>
+        the cosine of Pi*<parameter>x</parameter>
+      </refpurpose>
+      </refnamediv>
+      <refsynopsisdiv>
+        <synopsis><function>COSPI</function>(<parameter>x</parameter>)</synopsis>
+      </refsynopsisdiv>
+      <refsect1>
+        <title>Arguments</title>
+        <para><parameter>x</parameter>: number of half turns</para>
+      </refsect1>
+      <refsect1>
+        <title>See also</title>
+        <para><link linkend="gnumeric-function-COS"><function>COS</function></link>.
+      </para>
+      </refsect1>
+    </refentry>
     <refentry id="gnumeric-function-COT">
       <refmeta>
         <refentrytitle>
@@ -12668,6 +12695,33 @@
       </para>
       </refsect1>
     </refentry>
+    <refentry id="gnumeric-function-SINPI">
+      <refmeta>
+        <refentrytitle>
+          <function>SINPI</function>
+        </refentrytitle>
+      </refmeta>
+      <refnamediv>
+        <refname>
+          <function>SINPI</function>
+        </refname>
+        <refpurpose>
+        the sine of Pi*<parameter>x</parameter>
+      </refpurpose>
+      </refnamediv>
+      <refsynopsisdiv>
+        <synopsis><function>SINPI</function>(<parameter>x</parameter>)</synopsis>
+      </refsynopsisdiv>
+      <refsect1>
+        <title>Arguments</title>
+        <para><parameter>x</parameter>: number of half turns</para>
+      </refsect1>
+      <refsect1>
+        <title>See also</title>
+        <para><link linkend="gnumeric-function-SIN"><function>SIN</function></link>.
+      </para>
+      </refsect1>
+    </refentry>
     <refentry id="gnumeric-function-SQRT">
       <refmeta>
         <refentrytitle>
diff --git a/plugins/fn-math/ChangeLog b/plugins/fn-math/ChangeLog
index eb696da..482255d 100644
--- a/plugins/fn-math/ChangeLog
+++ b/plugins/fn-math/ChangeLog
@@ -1,3 +1,7 @@
+2013-11-14  Morten Welinder  <terra gnome org>
+
+       * functions.c (gnumeric_cospi, gnumeric_sinpi): New functions.
+
 2013-10-07  Morten Welinder <terra gnome org>
 
        * Release 1.12.8
diff --git a/plugins/fn-math/functions.c b/plugins/fn-math/functions.c
index 993a7f1..afd4537 100644
--- a/plugins/fn-math/functions.c
+++ b/plugins/fn-math/functions.c
@@ -780,6 +780,23 @@ gnumeric_cos (GnmFuncEvalInfo *ei, GnmValue const * const *argv)
 
 /***************************************************************************/
 
+static GnmFuncHelp const help_cospi[] = {
+       { GNM_FUNC_HELP_NAME, F_("COSPI:the cosine of Pi* {x}")},
+       { GNM_FUNC_HELP_ARG, F_("x:number of half turns")},
+       { GNM_FUNC_HELP_EXAMPLES, "=COSPI(0.5)" },
+       { GNM_FUNC_HELP_EXAMPLES, "=COSPI(1)" },
+       { GNM_FUNC_HELP_SEEALSO, "COS" },
+       { GNM_FUNC_HELP_END }
+};
+
+static GnmValue *
+gnumeric_cospi (GnmFuncEvalInfo *ei, GnmValue const * const *argv)
+{
+       return value_new_float (gnm_cospi (value_get_as_float (argv[0])));
+}
+
+/***************************************************************************/
+
 static GnmFuncHelp const help_cosh[] = {
        { GNM_FUNC_HELP_NAME, F_("COSH:the hyperbolic cosine of @{x}")},
        { GNM_FUNC_HELP_ARG, F_("x:number")},
@@ -1384,6 +1401,23 @@ gnumeric_sin (GnmFuncEvalInfo *ei, GnmValue const * const *argv)
 
 /***************************************************************************/
 
+static GnmFuncHelp const help_sinpi[] = {
+       { GNM_FUNC_HELP_NAME, F_("SINPI:the sine of Pi* {x}")},
+       { GNM_FUNC_HELP_ARG, F_("x:number of half turns")},
+       { GNM_FUNC_HELP_EXAMPLES, "=SINPI(0.5)" },
+       { GNM_FUNC_HELP_EXAMPLES, "=SINPI(1)" },
+       { GNM_FUNC_HELP_SEEALSO, "SIN" },
+       { GNM_FUNC_HELP_END }
+};
+
+static GnmValue *
+gnumeric_sinpi (GnmFuncEvalInfo *ei, GnmValue const * const *argv)
+{
+       return value_new_float (gnm_sinpi (value_get_as_float (argv[0])));
+}
+
+/***************************************************************************/
+
 static GnmFuncHelp const help_csc[] = {
        { GNM_FUNC_HELP_NAME, F_("CSC:the cosecant of @{x}")},
        { GNM_FUNC_HELP_ARG, F_("x:angle in radians")},
@@ -3215,6 +3249,9 @@ GnmFuncDescriptor const math_functions[] = {
        { "cosh",    "f",     help_cosh,
          gnumeric_cosh, NULL, NULL, NULL,
          GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_COMPLETE, GNM_FUNC_TEST_STATUS_EXHAUSTIVE },
+       { "cospi",   "f",     help_cospi,
+         gnumeric_cospi, NULL, NULL, NULL,
+         GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_UNIQUE_TO_GNUMERIC, GNM_FUNC_TEST_STATUS_EXHAUSTIVE },
        { "cot",     "f",     help_cot,
          gnumeric_cot, NULL, NULL, NULL,
          GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_COMPLETE, GNM_FUNC_TEST_STATUS_NO_TESTSUITE },
@@ -3374,6 +3411,9 @@ GnmFuncDescriptor const math_functions[] = {
        { "sinh",    "f",     help_sinh,
          gnumeric_sinh, NULL, NULL, NULL,
          GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_COMPLETE, GNM_FUNC_TEST_STATUS_EXHAUSTIVE },
+       { "sinpi",   "f",     help_sinpi,
+         gnumeric_sinpi, NULL, NULL, NULL,
+         GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_UNIQUE_TO_GNUMERIC, GNM_FUNC_TEST_STATUS_EXHAUSTIVE },
        { "sqrt",    "f",     help_sqrt,
          gnumeric_sqrt, NULL, NULL, NULL,
          GNM_FUNC_SIMPLE, GNM_FUNC_IMPL_STATUS_COMPLETE, GNM_FUNC_TEST_STATUS_BASIC },
diff --git a/plugins/fn-math/plugin.xml.in b/plugins/fn-math/plugin.xml.in
index 2559246..d0f3f66 100644
--- a/plugins/fn-math/plugin.xml.in
+++ b/plugins/fn-math/plugin.xml.in
@@ -32,6 +32,7 @@
                                <function name="combina"/>
                                <function name="cos"/>
                                <function name="cosh"/>
+                               <function name="cospi"/>
                                <function name="cot"/>
                                <function name="coth"/>
                                <function name="countif"/>
@@ -85,6 +86,7 @@
                                <function name="sign"/>
                                <function name="sin"/>
                                <function name="sinh"/>
+                               <function name="sinpi"/>
                                <function name="sqrt"/>
                                <function name="sqrtpi"/>
                                <function name="suma"/>
diff --git a/src/mathfunc.c b/src/mathfunc.c
index 3c11b5a..6c03c98 100644
--- a/src/mathfunc.c
+++ b/src/mathfunc.c
@@ -165,6 +165,86 @@ gnm_acoth (gnm_float x)
                : gnm_log ((x - 1) / (x + 1)) / -2;
 }
 
+/**
+ * gnm_sinpi:
+ * @x: a number
+ *
+ * Returns: the sine of Pi times @x, but with less error than doing the
+ * multiplication outright.
+ */
+gnm_float
+gnm_sinpi (gnm_float x)
+{
+       int k;
+
+       if (gnm_isnan (x))
+               return x;
+
+       if (!gnm_finite (x))
+               return gnm_nan;
+
+       k = (x < 0) ? 2 : 0;
+       x = gnm_fmod (gnm_abs (x), 2);
+       if (x >= 1) {
+               x -= 1;
+               k ^= 2;
+       }
+       if (x >= 0.5) {
+               x -= 0.5;
+               k += 1;
+       }
+       if (x == 0) {
+               static const gnm_float ys[4] = { 0, 1, -0, -1 };
+               return ys[k];
+       } else {
+               switch (k) {
+               default:
+               case 0: return +gnm_sin (M_PIgnum * x);
+               case 1: return +gnm_cos (M_PIgnum * x);
+               case 2: return -gnm_sin (M_PIgnum * x);
+               case 3: return -gnm_cos (M_PIgnum * x);
+               }
+       }
+}
+
+/**
+ * gnm_cospi:
+ * @x: a number
+ *
+ * Returns: the cosine of Pi times @x, but with less error than doing the
+ * multiplication outright.
+ */
+gnm_float
+gnm_cospi (gnm_float x)
+{
+       int k = 0;
+
+       if (!gnm_finite (x))
+               return gnm_nan;
+
+       x = gnm_fmod (gnm_abs (x), 2);
+       if (x >= 1) {
+               x -= 1;
+               k ^= 2;
+       }
+       if (x >= 0.5) {
+               x -= 0.5;
+               k += 1;
+       }
+       if (x == 0) {
+               static const gnm_float ys[4] = { 1, 0, -1, -0 };
+               return ys[k];
+       } else {
+               switch (k) {
+               default:
+               case 0: return +gnm_cos (M_PIgnum * x);
+               case 1: return -gnm_sin (M_PIgnum * x);
+               case 2: return -gnm_cos (M_PIgnum * x);
+               case 3: return +gnm_sin (M_PIgnum * x);
+               }
+       }
+}
+
 /* ------------------------------------------------------------------------- */
 /* --- BEGIN MAGIC R SOURCE MARKER --- */
 
@@ -3924,7 +4004,7 @@ static gnm_float bessel_i(gnm_float x, gnm_float alpha, gnm_float expo)
         * this may not be quite optimal (CPU and accuracy wise) */
        return(bessel_i(x, -alpha, expo) +
               bessel_k(x, -alpha, expo) * ((ize == 1)? 2. : 2.*gnm_exp(-x))/M_PIgnum
-              * gnm_sin(-M_PIgnum * gnm_fmod (alpha, 2)));
+              * gnm_sinpi(-alpha));
     }
     nb = 1+ (long)gnm_floor(alpha);/* nb-1 <= alpha < nb */
     alpha -= (nb-1);
@@ -4399,9 +4479,9 @@ static gnm_float bessel_j(gnm_float x, gnm_float alpha)
     if (alpha < 0) {
        /* Using Abramowitz & Stegun  9.1.2
         * this may not be quite optimal (CPU and accuracy wise) */
-       return(bessel_j(x, -alpha) * gnm_cos(M_PIgnum * gnm_fmod (alpha, 2)) +
+       return(bessel_j(x, -alpha) * gnm_cospi(alpha) +
               ((alpha == na) ? 0 :
-              bessel_y(x, -alpha) * gnm_sin(M_PIgnum * gnm_fmod (alpha, 2))));
+              bessel_y(x, -alpha) * gnm_sinpi(alpha)));
     }
     nb = 1 + (long)na; /* nb-1 <= alpha < nb */
     alpha -= (gnm_float)(nb-1);
@@ -4449,9 +4529,9 @@ static gnm_float bessel_j_ex(gnm_float x, gnm_float alpha, gnm_float *bj)
     if (alpha < 0) {
        /* Using Abramowitz & Stegun  9.1.2
         * this may not be quite optimal (CPU and accuracy wise) */
-       return(bessel_j_ex(x, -alpha, bj) * gnm_cos(M_PIgnum * gnm_fmod (alpha, 2)) +
+       return(bessel_j_ex(x, -alpha, bj) * gnm_cospi(alpha) +
               ((alpha == na) ? 0 :
-               bessel_y_ex(x, -alpha, bj) * gnm_sin(M_PIgnum * gnm_fmod (alpha, 2))));
+               bessel_y_ex(x, -alpha, bj) * gnm_sinpi(alpha)));
     }
     nb = 1 + (long)na; /* nb-1 <= alpha < nb */
     alpha -= (gnm_float)(nb-1);
@@ -5499,9 +5579,9 @@ static gnm_float bessel_y(gnm_float x, gnm_float alpha)
     if (alpha < 0) {
        /* Using Abramowitz & Stegun  9.1.2
         * this may not be quite optimal (CPU and accuracy wise) */
-       return(bessel_y(x, -alpha) * gnm_cos(M_PIgnum * gnm_fmod (alpha, 2)) -
+       return(bessel_y(x, -alpha) * gnm_cospi(alpha) -
               ((alpha == na) ? 0 :
-               bessel_j(x, -alpha) * gnm_sin(M_PIgnum * gnm_fmod (alpha, 2))));
+               bessel_j(x, -alpha) * gnm_sinpi(alpha)));
     }
     nb = 1+ (long)na;/* nb-1 <= alpha < nb */
     alpha -= (gnm_float)(nb-1);
@@ -5551,9 +5631,9 @@ static gnm_float bessel_y_ex(gnm_float x, gnm_float alpha, gnm_float *by)
     if (alpha < 0) {
        /* Using Abramowitz & Stegun  9.1.2
         * this may not be quite optimal (CPU and accuracy wise) */
-       return(bessel_y_ex(x, -alpha, by) * gnm_cos(M_PIgnum * gnm_fmod (alpha, 2)) -
+       return(bessel_y_ex(x, -alpha, by) * gnm_cospi(alpha) -
               ((alpha == na) ? 0 :
-               bessel_j_ex(x, -alpha, by) * gnm_sin(M_PIgnum * gnm_fmod (alpha, 2))));
+               bessel_j_ex(x, -alpha, by) * gnm_sinpi(alpha)));
     }
     nb = 1+ (long)na;/* nb-1 <= alpha < nb */
     alpha -= (gnm_float)(nb-1);
@@ -8158,13 +8238,11 @@ qfactf (gnm_float x, GnmQuad *mant, int *exp2)
                if (qfactf (-x - 1, mant, exp2))
                        res = 1;
                else {
-                       GnmQuad a, b;
-                       gnm_float xm2 = gnm_fmod (x, 2);
+                       GnmQuad b;
 
-                       gnm_quad_init (&a, -M_PIgnum); /* FIXME: Do better */
-                       gnm_quad_init (&b, gnm_sin (xm2 * M_PIgnum)); /* ? */
+                       gnm_quad_init (&b, -gnm_sinpi (x)); /* ? */
                        gnm_quad_mul (&b, &b, mant);
-                       gnm_quad_div (mant, &a, &b);
+                       gnm_quad_div (mant, &gnm_quad_pi, &b);
                        *exp2 = -*exp2;
                }
        } else if (x >= QFACTI_LIMIT - 0.5) {
@@ -8698,7 +8776,7 @@ lgamma_r (double x, int *signp)
                        x + lgammacor(x)) - gnm_log (f);
        } else {
                gnm_float axm2 = gnm_fmod (-x, 2.0);
-               gnm_float y = gnm_sin (M_PIgnum * axm2) / M_PIgnum;
+               gnm_float y = gnm_sinpi (axm2) / M_PIgnum;
 
                *signp = axm2 > 1.0 ? +1 : -1;
 
diff --git a/src/mathfunc.h b/src/mathfunc.h
index 031e526..0c2b019 100644
--- a/src/mathfunc.h
+++ b/src/mathfunc.h
@@ -44,6 +44,9 @@ gnm_float gnm_acot (gnm_float x);
 gnm_float gnm_coth (gnm_float x);
 gnm_float gnm_acoth (gnm_float x);
 
+gnm_float gnm_sinpi (gnm_float x);
+gnm_float gnm_cospi (gnm_float x);
+
 gnm_float beta (gnm_float a, gnm_float b);
 gnm_float lbeta3 (gnm_float a, gnm_float b, int *sign);
 


[Date Prev][Date Next]   [Thread Prev][Thread Next]   [Thread Index] [Date Index] [Author Index]