Changeset 217590b in sasmodels for sasmodels/models/fractal.c

Timestamp:

Oct 20, 2016 4:48:48 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:

5c94f41

Parents:

2b9e63f

git-author:

Paul Kienzle <pkienzle@…> (10/20/16 16:47:25)

git-committer:

Paul Kienzle <pkienzle@…> (10/20/16 16:48:48)

Message:

fractal, fractal_core_shell: support fractal_dim as low as 0.

File:

• sasmodels/models/fractal.c (modified) (3 diffs)

sasmodels/models/fractal.c

@@ -1 @@
-#define INVALID(p) (p.fractal_dim <= 0.0)
+#define INVALID(p) (p.fractal_dim < 0.0)
 
 static double
@@ -10 +10 @@
    double sld_solvent)
 {
+    const double sq = fractal_sq(q, radius, fractal_dim, cor_length);
+
     //calculate P(q) for the spherical subunits
-    const double pq = M_4PI_3*cube(radius) * square((sld_block-sld_solvent)*sph_j1c(q*radius));
-    
-    //calculate S(q),  using Teixeira, Eq(15)
-    double sq;
-    if (q > 0. && fractal_dim > 1.) {
-        // q>0, D>0
-        const double D = fractal_dim;
-        const double Dm1 = fractal_dim - 1.0;
-        // Note: for large Dm1, sin(Dm1*atan(q*cor_length) can go negative
-        const double t1 = D*sas_gamma(Dm1)*sin(Dm1*atan(q*cor_length));
-        const double t2 = pow(q*radius, -D);
-        const double t3 = pow(1.0 + 1.0/square(q*cor_length), -0.5*Dm1);
-        sq = 1.0 + t1 * t2 * t3;
-    } else if (q > 0.) {
-        // q>0, D=1
-        sq = 1.0 + atan(q*cor_length) / (q*radius);
-    } else if (fractal_dim > 1.) {
-        // q=0, D>1
-        const double D = fractal_dim;
-        sq = 1.0 + pow(cor_length/radius, D)*sas_gamma(D+1.0);
-    } else {
-        // q=0, D=1
-        sq = 1.0 + cor_length/radius;
-    }
+    const double V = M_4PI_3*cube(radius);
+    const double pq = V * square((sld_block-sld_solvent)*sph_j1c(q*radius));
 
     // scale to units cm-1 sr-1 (assuming data on absolute scale)
@@ -41 +21 @@
     //    combined: 1e-4 * I(q)
 
-    return 1.e-4 * volfraction * pq * sq;
+    return 1.e-4 * volfraction * sq * pq;
 }