Changeset 463eb76e in sasview for sansmodels/src

Timestamp:

Jul 25, 2011 3:43:34 PM (13 years ago)

Author:

Jae Cho <jhjcho@…>

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:

bc7bc839

Parents:

38155f7

Message:

implementation of pearlneckclace model: need utest and doc

Location:

sansmodels/src/sans/models

Files:

• PearlNecklaceModel.py (modified) (1 diff)
• c_extensions/libmultifunc/libfunc.c (modified) (3 diffs)
• c_extensions/libmultifunc/libfunc.h (modified) (1 diff)
• c_extensions/pearlnecklace.c (modified) (4 diffs)
• c_extensions/pearlnecklace.h (modified) (1 diff)
• c_models/CPearlNecklaceModel.cpp (modified) (1 diff)
• c_models/pearlnecklace.cpp (modified) (1 diff)

sansmodels/src/sans/models/PearlNecklaceModel.py

@@ -47 +47 @@
          num_pearls      = 3.0 
          sld_pearl       = 1e-006 [1/A^(2)]
-         sld_string      = 5e-006 [1/A^(2)]
-         sld_solv        = 5e-007 [1/A^(2)]
+         sld_string      = 1e-006 [1/A^(2)]
+         sld_solv        = 6.3e-006 [1/A^(2)]
          background      = 0.0 
 

sansmodels/src/sans/models/c_extensions/libmultifunc/libfunc.c

@@ +1 @@
+// by jcho
 #include <math.h>
 #include "libmultifunc/libfunc.h"
@@ -19 +20 @@
 
 // Used in pearl nec model
-//integral of sin(x)/x up to namx term nmax>10 looks fine.
-double Si(double x, int nmax )
+// Sine integral function: approximated within 1%!!!
+// integral of sin(x)/x up to namx term nmax=6 looks the best.
+double Si(double x)
 {
         int i;
+        int nmax=6;
         double out;
         long double power;
-        if (x > 5.0){
-                double pi = 4.0*atan(1.0);
-                return pi/2.0;
+        double pi = 4.0*atan(1.0);
+        if (x >= pi*6.2/4.0){
+
+                double out_sin = 0.0;
+                double out_cos = 0.0;
+                out = pi/2.0;
+                for (i=0; i<nmax-2; i+=1){
+                        out_cos += pow(-1.0, i) * (double)factorial(2*i) / pow(x, 2*i+1);
+                        out_sin += pow(-1.0, i) * (double)factorial(2*i+1) / pow(x, 2*i+2);
+                }
+                out -= cos(x) * out_cos;
+                out -= sin(x) * out_sin;
+                return out;
         }
         out = 0.0;
@@ -43 +56 @@
 }
 
-double sin_x(double x)
+double sinc(double x)
 {
         if (x==0){

sansmodels/src/sans/models/c_extensions/libmultifunc/libfunc.h

@@ -4 +4 @@
 int factorial(int i);
 
-double Si(double x, int nmax );
+double Si(double x);
 
-double sin_x(double x);
+double sinc(double x);
 
 #endif

sansmodels/src/sans/models/c_extensions/pearlnecklace.c

@@ -53 +53 @@
         // Note take only 20 terms in Si series: 10 terms may be enough though.
         double gamma;
-        gamma = Si(q* edge_separation, 6);
+        gamma = Si(q* edge_separation);
         gamma /= (q* edge_separation);
         double beta;
-        beta = Si(q * (edge_separation + radius), 6);
-        beta -= Si(q * radius, 6);
+        beta = Si(q * (edge_separation + radius));
+        beta -= Si(q * radius);
         beta /= (q* edge_separation);
 
@@ -68 +68 @@
 
         // form factor for num_pearls
-        sss = 1.0 - pow(sin_x(q*A_s), (int)num_pearls );
-        sss /= pow((1.0-sin_x(q*A_s)), 2);
-        sss *= sin_x(q*A_s);
+        sss = 1.0 - pow(sinc(q*A_s), num_pearls );
+        sss /= pow((1.0-sinc(q*A_s)), 2);
+        sss *= -sinc(q*A_s);
         sss -= num_pearls/2.0;
-        sss += num_pearls/(1.0-sin_x(q*A_s));
+        sss += num_pearls/(1.0-sinc(q*A_s));
         sss *= 2.0 * pow((m_s*psi), 2);
 
         // form factor for num_strings (like thin rods)
         double srr_1;
-        srr_1 = -pow(sin_x(q*edge_separation/2.0), 2);
-        //printf ("sss %g, %g, %g\n",sss, srr_1, gamma);
+        srr_1 = -pow(sinc(q*edge_separation/2.0), 2);
+
         srr_1 += 2.0 * gamma;
         srr_1 *= num_strings;
         double srr_2;
-        srr_2 = 2.0/(1.0-sin_x(q*A_s));
+        srr_2 = 2.0/(1.0-sinc(q*A_s));
         srr_2 *= num_strings;
         srr_2 *= pow(beta, 2);
         double srr_3;
-        srr_3 = 1.0 - pow(sin_x(q*A_s), (int)num_strings);
-        srr_3 /= pow((1.0-sin_x(q*A_s)), 2);
+        srr_3 = 1.0 - pow(sinc(q*A_s), num_strings);
+        srr_3 /= pow((1.0-sinc(q*A_s)), 2);
         srr_3 *= pow(beta, 2);
         srr_3 *= -2.0;
@@ -96 +96 @@
 
         // form factor for correlations
-        srs = -1.0;
-        srs -= pow(sin_x(q*A_s), (int)num_strings);
-        srs /= pow((1.0-sin_x(q*A_s)), 2);
-        srs *= sin_x(q*A_s);
-        srs += (num_strings/(1.0-sin_x(q*A_s)));
+        srs = 1.0;
+        srs -= pow(sinc(q*A_s), num_strings);
+        srs /= pow((1.0-sinc(q*A_s)), 2);
+        srs *= -sinc(q*A_s);
+        srs += (num_strings/(1.0-sinc(q*A_s)));
         srs *= 4.0;
         srs *= (m_r * m_s * beta * psi);
@@ -111 +111 @@
         form_factor *= scale;
         form_factor += background;
-
     return (form_factor);
 }

sansmodels/src/sans/models/c_extensions/pearlnecklace.h

@@ -44 +44 @@
         double sld_pearl;
         ///     sld_string
-        //  [DEFAULT]=sld_string= 5.0e-06 [1/A^(2)]
+        //  [DEFAULT]=sld_string= 1.0e-06 [1/A^(2)]
         double sld_string;
         ///     sld_solv
-        //  [DEFAULT]=sld_solv= 0.5e-06 [1/A^(2)]
+        //  [DEFAULT]=sld_solv= 6.3e-06 [1/A^(2)]
         double sld_solv;
         /// Background

sansmodels/src/sans/models/c_models/CPearlNecklaceModel.cpp

@@ -93 +93 @@
         PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000000000));
         PyDict_SetItemString(self->params,"thick_string",Py_BuildValue("d",2.500000000000));
-        PyDict_SetItemString(self->params,"sld_string",Py_BuildValue("d",0.000005000000));
+        PyDict_SetItemString(self->params,"sld_string",Py_BuildValue("d",0.000001000000));
         PyDict_SetItemString(self->params,"edge_separation",Py_BuildValue("d",350.000000000000));
         PyDict_SetItemString(self->params,"sld_pearl",Py_BuildValue("d",0.000001000000));
         PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",80.000000000000));
         PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000000000));
-        PyDict_SetItemString(self->params,"sld_solv",Py_BuildValue("d",0.000000500000));
+        PyDict_SetItemString(self->params,"sld_solv",Py_BuildValue("d",0.000006300000));
         // Initialize dispersion / averaging parameter dict
         DispersionVisitor* visitor = new DispersionVisitor();

sansmodels/src/sans/models/c_models/pearlnecklace.cpp

@@ -21 +21 @@
         num_pearls.set_min(0.0);
         sld_pearl = Parameter(1.0e-06);
-        sld_string = Parameter(5.0e-06);
-        sld_solv = Parameter(0.5e-06);
+        sld_string = Parameter(1.0e-06);
+        sld_solv = Parameter(6.3e-06);
     background = Parameter(0.0);