Changeset 71b751d in sasmodels for sasmodels/models/elliptical_cylinder.c

Timestamp:

Aug 14, 2018 10:09:34 AM (6 years ago)

Author:

Paul Kienzle <pkienzle@…>

Branches:

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

Children:

86aa992

Parents:

2f8cbb9

Message:

update remaining form factors to use Fq interface

File:

• sasmodels/models/elliptical_cylinder.c (modified) (4 diffs)

sasmodels/models/elliptical_cylinder.c

@@ -5 +5 @@
 }
 
-static double
-Iq(double q, double radius_minor, double r_ratio, double length,
+static void
+Fq(double q, double *F1, double *F2, double radius_minor, double r_ratio, double length,
    double sld, double solvent_sld)
 {
@@ -21 +21 @@
 
     //initialize integral
-    double outer_sum = 0.0;
+    double outer_sum_F1 = 0.0;
+    double outer_sum_F2 = 0.0;
     for(int i=0;i<GAUSS_N;i++) {
         //setup inner integral over the ellipsoidal cross-section
@@ -27 +28 @@
         const double sin_val = sqrt(1.0 - cos_val*cos_val);
         //const double arg = radius_minor*sin_val;
-        double inner_sum=0;
+        double inner_sum_F1 = 0.0;
+        double inner_sum_F2 = 0.0;
         for(int j=0;j<GAUSS_N;j++) {
             const double theta = ( GAUSS_Z[j]*(vbj-vaj) + vaj + vbj )/2.0;
             const double r = sin_val*sqrt(rA - rB*cos(theta));
             const double be = sas_2J1x_x(q*r);
-            inner_sum += GAUSS_W[j] * be * be;
+            inner_sum_F1 += GAUSS_W[j] * be;
+            inner_sum_F2 += GAUSS_W[j] * be * be;
         }
         //now calculate the value of the inner integral
-        inner_sum *= 0.5*(vbj-vaj);
+        inner_sum_F1 *= 0.5*(vbj-vaj);
+        inner_sum_F2 *= 0.5*(vbj-vaj);
 
         //now calculate outer integral
         const double si = sas_sinx_x(q*0.5*length*cos_val);
-        outer_sum += GAUSS_W[i] * inner_sum * si * si;
+        outer_sum_F1 += GAUSS_W[i] * inner_sum_F1 * si;
+        outer_sum_F2 += GAUSS_W[i] * inner_sum_F2 * si * si;
     }
-    outer_sum *= 0.5*(vb-va);
-
-    //divide integral by Pi
-    const double form = outer_sum/M_PI;
+    // correct limits and divide integral by pi
+    outer_sum_F1 *= 0.5*(vb-va)/M_PI;
+    outer_sum_F2 *= 0.5*(vb-va)/M_PI;
 
     // scale by contrast and volume, and convert to to 1/cm units
-    const double vol = form_volume(radius_minor, r_ratio, length);
-    const double delrho = sld - solvent_sld;
-    return 1.0e-4*square(delrho*vol)*form;
+    const double volume = form_volume(radius_minor, r_ratio, length);
+    const double contrast = sld - solvent_sld;
+    const double s = contrast*volume;
+    *F1 = 1.0e-2*s*outer_sum_F1;
+    *F2 = 1.0e-4*s*s*outer_sum_F2;
 }
 
@@ -63 +69 @@
     const double be = sas_2J1x_x(qr);
     const double si = sas_sinx_x(qc*0.5*length);
-    const double Aq = be * si;
-    const double delrho = sld - solvent_sld;
-    const double vol = form_volume(radius_minor, r_ratio, length);
-    return 1.0e-4 * square(delrho * vol * Aq);
+    const double fq = be * si;
+    const double contrast = sld - solvent_sld;
+    const double volume = form_volume(radius_minor, r_ratio, length);
+    return 1.0e-4 * square(contrast * volume * fq);
 }