Changeset a5b6997 in sasmodels for sasmodels/models/lib

Timestamp:

Oct 14, 2016 6:07:28 PM (8 years ago)

Author:

Paul Kienzle <pkienzle@…>

Branches:

master, core_shell_microgels, costrafo411, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests

Children:

4962519

Parents:

3a48772

Message:

Si: Use Horner's method when computing the Sine integral

File:

• sasmodels/models/lib/Si.c (modified) (1 diff)

sasmodels/models/lib/Si.c

@@ -3 +3 @@
 double Si(double x)
 {
-    double out;
+    if (x >= M_PI*6.2/4.0){
+        const double z = 1./(x*x);
+        // Explicitly writing factorial values triples the speed of the calculation
+        const double out_cos = (((-720.*z + 24.)*z - 2.)*z + 1.)/x;
+        const double out_sin = (((-5040.*z + 120.)*z - 6.)*z + 1)*z;
 
-    if (x >= M_PI*6.2/4.0){
-        double out_sin = 0.0;
-        double out_cos = 0.0;
-        out = M_PI/2.0;
-
+        double cos_x, sin_x;
+        SINCOS(x, cos_x, sin_x);
+        return M_PI_2 - cos_x*out_cos - sin_x*out_sin;
+    } else {
+        const double z = x*x;
         // Explicitly writing factorial values triples the speed of the calculation
-        out_cos = 1./x - 2./pow(x,3) + 24./pow(x,5) - 720./pow(x,7);
-        out_sin = 1./pow(x,2) - 6./pow(x,4) + 120./pow(x,6) - 5040./pow(x,8);
-
-        out -= cos(x) * out_cos;
-        out -= sin(x) * out_sin;
-        return out;
+        return (((((-1./439084800.*z
+            + 1./3265920.)*z
+            - 1./35280.)*z
+            + 1./600.)*z
+            - 1./18.)*z
+            + 1.)*x;
     }
-
-    // Explicitly writing factorial values triples the speed of the calculation
-    out = x - pow(x, 3)/18. + pow(x,5)/600. - pow(x,7)/35280. + pow(x,9)/3265920. - pow(x,11)/439084800.;
-
-    return out;
 }