Changeset 3371d81 in sasview for src/sas/sascalc/calculator/c_extensions

Timestamp:

May 9, 2016 6:59:33 AM (8 years ago)

Author:

wojciech

Branches:

master, ESS_GUI, ESS_GUI_Docs, ESS_GUI_batch_fitting, ESS_GUI_bumps_abstraction, ESS_GUI_iss1116, ESS_GUI_iss879, ESS_GUI_iss959, ESS_GUI_opencl, ESS_GUI_ordering, ESS_GUI_sync_sascalc, costrafo411, magnetic_scatt, release-4.1.1, release-4.1.2, release-4.2.2, release_4.0.1, ticket-1009, ticket-1094-headless, ticket-1242-2d-resolution, ticket-1243, ticket-1249, ticket885, unittest-saveload

Children:

aa4a549

Parents:

aff7904 (diff), 1d5f5c2 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.

Message:

setup.py fixed. Merged with Master

Location:

src/sas/sascalc/calculator/c_extensions

Files:

• librefl.c (modified) (1 diff)
• sld2i.cpp (modified) (2 diffs)
• sld2i.hh (modified) (1 diff)
• sld2i_module.cpp (modified) (1 diff)
• winFuncs.c (deleted)
• winFuncs.h (deleted)

src/sas/sascalc/calculator/c_extensions/librefl.c

@@ -8 +8 @@
 #include <stdlib.h>
 #if defined(_MSC_VER)
-#include "winFuncs.h"
+#define NEED_ERF
 #endif
+
+
+
+#if defined(NEED_ERF)
+/* erf.c  - public domain implementation of error function erf(3m)
+
+reference - Haruhiko Okumura: C-gengo niyoru saishin algorithm jiten
+            (New Algorithm handbook in C language) (Gijyutsu hyouron
+            sha, Tokyo, 1991) p.227 [in Japanese]                 */
+
+
+#ifdef _WIN32
+# include <float.h>
+# if !defined __MINGW32__ || defined __NO_ISOCEXT
+#  ifndef isnan
+#   define isnan(x) _isnan(x)
+#  endif
+#  ifndef isinf
+#   define isinf(x) (!_finite(x) && !_isnan(x))
+#  endif
+#  ifndef finite
+#   define finite(x) _finite(x)
+#  endif
+# endif
+#endif
+
+static double q_gamma(double, double, double);
+
+/* Incomplete gamma function
+   1 / Gamma(a) * Int_0^x exp(-t) t^(a-1) dt  */
+static double p_gamma(double a, double x, double loggamma_a)
+{
+    int k;
+    double result, term, previous;
+
+    if (x >= 1 + a) return 1 - q_gamma(a, x, loggamma_a);
+    if (x == 0)     return 0;
+    result = term = exp(a * log(x) - x - loggamma_a) / a;
+    for (k = 1; k < 1000; k++) {
+        term *= x / (a + k);
+        previous = result;  result += term;
+        if (result == previous) return result;
+    }
+    fprintf(stderr, "erf.c:%d:p_gamma() could not converge.", __LINE__);
+    return result;
+}
+
+/* Incomplete gamma function
+   1 / Gamma(a) * Int_x^inf exp(-t) t^(a-1) dt  */
+static double q_gamma(double a, double x, double loggamma_a)
+{
+    int k;
+    double result, w, temp, previous;
+    double la = 1, lb = 1 + x - a;  /* Laguerre polynomial */
+
+    if (x < 1 + a) return 1 - p_gamma(a, x, loggamma_a);
+    w = exp(a * log(x) - x - loggamma_a);
+    result = w / lb;
+    for (k = 2; k < 1000; k++) {
+        temp = ((k - 1 - a) * (lb - la) + (k + x) * lb) / k;
+        la = lb;  lb = temp;
+        w *= (k - 1 - a) / k;
+        temp = w / (la * lb);
+        previous = result;  result += temp;
+        if (result == previous) return result;
+    }
+    fprintf(stderr, "erf.c:%d:q_gamma() could not converge.", __LINE__);
+    return result;
+}
+
+#define LOG_PI_OVER_2 0.572364942924700087071713675675 /* log_e(PI)/2 */
+
+double erf(double x)
+{
+    if (!finite(x)) {
+        if (isnan(x)) return x;      /* erf(NaN)   = NaN   */
+        return (x>0 ? 1.0 : -1.0);   /* erf(+-inf) = +-1.0 */
+    }
+    if (x >= 0) return   p_gamma(0.5, x * x, LOG_PI_OVER_2);
+    else        return - p_gamma(0.5, x * x, LOG_PI_OVER_2);
+}
+
+double erfc(double x)
+{
+    if (!finite(x)) {
+        if (isnan(x)) return x;      /* erfc(NaN)   = NaN      */
+        return (x>0 ? 0.0 : 2.0);    /* erfc(+-inf) = 0.0, 2.0 */
+    }
+    if (x >= 0) return  q_gamma(0.5, x * x, LOG_PI_OVER_2);
+    else        return  1 + p_gamma(0.5, x * x, LOG_PI_OVER_2);
+}
+#endif // NEED_ERF
 
 complex cassign(real, imag)

src/sas/sascalc/calculator/c_extensions/sld2i.cpp

@@ -164 +164 @@
         //Assume that pixel volumes are given in vol_pix in A^3 unit
         // Loop over q-values and multiply apply matrix
-        for(size_t i=0; i<npoints; i++){
+        for(int i=0; i<npoints; i++){
                 sumj =0.0;              
-                for(size_t j=0; j<n_pix; j++){
+                for(int j=0; j<n_pix; j++){
                         //Isotropic: Assumes all slds are real (no magnetic)
                         //Also assumes there is no polarization: No dependency on spin
@@ -183 +183 @@
                                 //full calculation
                                 //pragma omp parallel for
-                                for(size_t k=0; k<n_pix; k++){
+                                for(int k=0; k<n_pix; k++){
                                         sld_j =  sldn_val[j] * sldn_val[k] * vol_pix[j] * vol_pix[k];
                                         qr = (x_val[j]-x_val[k])*(x_val[j]-x_val[k])+

src/sas/sascalc/calculator/c_extensions/sld2i.hh

@@ -39 +39 @@
                         double s_theta);
         // compute function
-        virtual void genicomXY(int npoints, double* qx, double* qy, double *I_out);
-        virtual void genicom(int npoints, double* q, double *I_out);
+        void genicomXY(int npoints, double* qx, double* qy, double *I_out);
+        void genicom(int npoints, double* q, double *I_out);
 };
 

src/sas/sascalc/calculator/c_extensions/sld2i_module.cpp

@@ -153 +153 @@
 initsld2i(void)
 {
-    PyObject* m;
-
-    m = Py_InitModule3("sld2i", module_methods,
-                       "Sld2i module");
+    Py_InitModule3("sld2i", module_methods, "Sld2i module");
 }