Changeset 71b751d in sasmodels for sasmodels/models/spherical_sld.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/spherical_sld.c (modified) (1 diff)

sasmodels/models/spherical_sld.c

@@ -101 +101 @@
 }
 
+static void Fq(
+    double q,
+    double *F1,
+    double *F2,
+    double fp_n_shells,
+    double sld_solvent,
+    double sld[],
+    double thickness[],
+    double interface[],
+    double shape[],
+    double nu[],
+    double fp_n_steps)
+{
+    // iteration for # of shells + core + solvent
+    int n_shells = (int)(fp_n_shells + 0.5);
+    int n_steps = (int)(fp_n_steps + 0.5);
+    double f=0.0;
+    double r=0.0;
+    for (int shell=0; shell<n_shells; shell++){
+        const double sld_l = sld[shell];
+
+        // uniform shell; r=0 => r^3=0 => f=0, so works for core as well.
+        f -= M_4PI_3 * cube(r) * sld_l * sas_3j1x_x(q*r);
+        r += thickness[shell];
+        f += M_4PI_3 * cube(r) * sld_l * sas_3j1x_x(q*r);
+
+        // iterate over sub_shells in the interface
+        const double dr = interface[shell]/n_steps;
+        const double delta = (shell==n_shells-1 ? sld_solvent : sld[shell+1]) - sld_l;
+        const double nu_shell = fmax(fabs(nu[shell]), 1.e-14);
+        const int shape_shell = (int)(shape[shell]);
+
+        // if there is no interface the equations don't work
+        if (dr == 0.) continue;
+
+        double sld_in = sld_l;
+        for (int step=1; step <= n_steps; step++) {
+            // find sld_i at the outer boundary of sub-shell step
+            //const double z = (double)step/(double)n_steps;
+            const double z = (double)step/(double)n_steps;
+            const double fraction = blend(shape_shell, nu_shell, z);
+            const double sld_out = fraction*delta + sld_l;
+            // calculate slope
+            const double slope = (sld_out - sld_in)/dr;
+            const double contrast = sld_in - slope*r;
+
+            // iteration for the left and right boundary of the shells
+            f -= f_linear(q, r, contrast, slope);
+            r += dr;
+            f += f_linear(q, r, contrast, slope);
+            sld_in = sld_out;
+        }
+    }
+    // add in solvent effect
+    f -= M_4PI_3 * cube(r) * sld_solvent * sas_3j1x_x(q*r);
+
+    *F1 = 1e-2*f;
+    *F2 = 1e-4*f*f;
+}