Changeset 7e0b281 in sasmodels


Ignore:
Timestamp:
Apr 17, 2017 6:34:52 PM (8 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:
64ca163
Parents:
cb038a2
Message:

adjust paracrystal equations to better match terminology in Matsuoka paper

Files:
5 edited

Legend:

Unmodified
Added
Removed

  • explore/bccpy.py

    rcb038a2 r7e0b281  
    4747SLD_SOLVENT = 6.3 
    4848 
     49# Note: using Matsuoka 1990; this is different from what 
     50# is in the sasmodels/models code  (see bcc vs bcc_old). 
     51# The difference is that the sign of phi and theta seem to be 
     52# negative in the old vs. the new, yielding a pattern that is 
     53# swapped left to right and top to bottom. 
     54def sc(qa, qb, qc): 
     55    return qa, qb, qc 
     56 
     57def bcc(qa, qb, qc): 
     58    a1 = (+qa + qb + qc)/2 
     59    a2 = (-qa - qb + qc)/2 
     60    a3 = (-qa + qb - qc)/2 
     61    return a1, a2, a3 
     62 
     63def bcc_old(qa, qb, qc): 
     64    a1 = (+qa + qb - qc)/2.0 
     65    a2 = (+qa - qb + qc)/2.0 
     66    a3 = (-qa + qb + qc)/2.0 
     67    return a1, a2, a3 
     68 
     69def fcc(qa, qb, qc): 
     70    a1 = ( 0. + qb + qc)/2 
     71    a2 = (-qa + 0. + qc)/2 
     72    a3 = (-qa + qb + 0.)/2 
     73    return a1, a2, a3 
     74 
     75def fcc_old(qa, qb, qc): 
     76    a1 = ( qa + qb + 0.)/2 
     77    a2 = ( qa + 0. + qc)/2 
     78    a3 = ( 0. + qb + qc)/2 
     79    return a1, a2, a3 
     80 
     81KERNEL = bcc 
     82 
    4983def kernel(q, dnn, d_factor, theta, phi): 
    5084    """ 
    … …  
    5690    qc = q*cos(theta) 
    5791 
    58     if 0: # sc 
    59         a1, a2, a3 = qa, qb, qc 
    60         dcos = dnn 
    61     if 1: # bcc 
    62         a1 = +qa - qc + qb 
    63         a2 = +qa + qc - qb 
    64         a3 = -qa + qc + qb 
    65         dcos = dnn/2 
    66     if 0: # fcc 
    67         a1 = qb + qa 
    68         a2 = qa + qc 
    69         a3 = qb + qc 
    70         dcos = dnn/2 
    71  
    72     arg = 0.5*square(dnn*d_factor)*(a1**2 + a2**2 + a3**2) 
    73     exp_arg = exp(-arg) 
    74     den = [((exp_arg - 2*cos(dcos*a))*exp_arg + 1.0) for a in (a1, a2, a3)] 
    75     Sq = -expm1(-2*arg)**3/np.prod(den, axis=0) 
    76     return Sq 
     92    a1, a2, a3 = KERNEL(qa, qb, qc) 
     93 
     94    # Note: paper says that different directions can have different distortion 
     95    # factors.  Easy enough to add to the code.  This would definitely break 
     96    # 8-fold symmetry. 
     97    arg = -0.5*square(dnn*d_factor)*(a1**2 + a2**2 + a3**2) 
     98    exp_arg = exp(arg) 
     99    den = [((exp_arg - 2*cos(dnn*a))*exp_arg + 1.0) for a in (a1, a2, a3)] 
     100    Zq = -expm1(2*arg)**3/np.prod(den, axis=0) 
     101    return Zq 
    77102 
    78103 
    … …  
    85110    def integrand(theta, phi): 
    86111        evals[0] += 1 
    87         Sq = kernel(q=q, dnn=dnn, d_factor=d_factor, theta=theta, phi=phi) 
    88         return Sq*sin(theta) 
     112        Zq = kernel(q=q, dnn=dnn, d_factor=d_factor, theta=theta, phi=phi) 
     113        return Zq*sin(theta) 
    89114    ans = dblquad(integrand, 0, pi/2, lambda x: 0, lambda x: pi/2)[0]*8/(4*pi) 
    90115    print("dblquad evals =", evals[0]) 
    … …  
    134159    Aw = w[None, :] * w[:, None] 
    135160    sin_theta = np.fmax(abs(sin(Atheta)), 1e-6) 
    136     Sq = kernel(q=q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
     161    Zq = kernel(q=q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
    137162    print("gauss %d evals ="%n, n**2) 
    138     return np.sum(Sq*Aw*sin_theta)*8/(4*pi) 
     163    return np.sum(Zq*Aw*sin_theta)*8/(4*pi) 
    139164 
    140165 
    … …  
    148173    phi = np.linspace(0, pi/2, n) 
    149174    Atheta, Aphi = np.meshgrid(theta, phi) 
    150     Sq = kernel(q=Q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
    151     Sq *= abs(sin(Atheta)) 
     175    Zq = kernel(q=Q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
     176    Zq *= abs(sin(Atheta)) 
    152177    dx, dy = theta[1]-theta[0], phi[1]-phi[0] 
    153     print("rect", n, np.sum(Sq)*dx*dy*8/(4*pi)) 
    154     print("trapz", n, np.trapz(np.trapz(Sq, dx=dx), dx=dy)*8/(4*pi)) 
    155     print("simpson", n, simps(simps(Sq, dx=dx), dx=dy)*8/(4*pi)) 
    156     print("romb", n, romb(romb(Sq, dx=dx), dx=dy)*8/(4*pi)) 
     178    print("rect", n, np.sum(Zq)*dx*dy*8/(4*pi)) 
     179    print("trapz", n, np.trapz(np.trapz(Zq, dx=dx), dx=dy)*8/(4*pi)) 
     180    print("simpson", n, simps(simps(Zq, dx=dx), dx=dy)*8/(4*pi)) 
     181    print("romb", n, romb(romb(Zq, dx=dx), dx=dy)*8/(4*pi)) 
    157182    print("gridded %d evals ="%n, n**2) 
    158183 
    … …  
    168193    #phi = np.linspace(0, pi/2, n) 
    169194    Atheta, Aphi = np.meshgrid(theta, phi) 
    170     Sq = kernel(q=Q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
    171     Sq *= abs(sin(Atheta)) 
    172     pylab.pcolor(degrees(theta), degrees(phi), log10(np.fmax(Sq, 1.e-6))) 
     195    Zq = kernel(q=Q, dnn=dnn, d_factor=d_factor, theta=Atheta, phi=Aphi) 
     196    Zq *= abs(sin(Atheta)) 
     197    pylab.pcolor(degrees(theta), degrees(phi), log10(np.fmax(Zq, 1.e-6))) 
    173198    pylab.axis('tight') 
    174     pylab.title("BCC S(q) for q=%g, dnn=%g d_factor=%g" % (q, dnn, d_factor)) 
     199    pylab.title("%s Z(q) for q=%g, dnn=%g d_factor=%g" 
     200                % (KERNEL.__name__, q, dnn, d_factor)) 
    175201    pylab.xlabel("theta (degrees)") 
    176202    pylab.ylabel("phi (degrees)") 

  • explore/sc.c

    rfdd56a1 r7e0b281  
    11static double 
    2 _sq_sc(double qa, double qb, double qc, double dnn, double d_factor) 
     2sc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    44    // Rewriting equations for efficiency, accuracy and readability, and so 
    … …  
    88    const double a3 = qc; 
    99 
    10     const double arg = 0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     10    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    1111 
    1212    // Numerator: (1 - exp(a)^2)^3 
    … …  
    1616    //         => exp(a)^2 - 2 cos(xk) exp(a) + 1 
    1717    //         => (exp(a) - 2 cos(xk)) * exp(a) + 1 
    18     const double exp_arg = exp(-arg); 
    19     const double Sq = -cube(expm1(-2.0*arg)) 
     18    const double exp_arg = exp(arg); 
     19    const double Zq = -cube(expm1(2.0*arg)) 
    2020        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
    2121          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
    2222          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2323 
    24     return Sq; 
     24    return Zq; 
    2525} 
    2626 
    2727// occupied volume fraction calculated from lattice symmetry and sphere radius 
    2828static double 
    29 _sc_volume_fraction(double radius, double dnn) 
     29sc_volume_fraction(double radius, double dnn) 
    3030{ 
    3131    return sphere_volume(radius/dnn); 
    … …  
    9898            const double qa = qab*cos_phi; 
    9999            const double qb = qab*sin_phi; 
    100             const double fq = _sq_sc(qa, qb, qc, dnn, d_factor); 
    101             inner_sum += fq; 
     100            const double form = sc_Zq(qa, qb, qc, dnn, d_factor); 
     101            inner_sum += form; 
    102102        } 
    103103        inner_sum *= phi_m;  // sum(f(x)dx) = sum(f(x)) dx 
    … …  
    106106    outer_sum *= theta_m/(n*n); 
    107107#endif 
    108 double Sq; 
     108double Zq; 
    109109if (sym > 0.) { 
    110     Sq = outer_sum/M_PI_2; 
     110    Zq = outer_sum/M_PI_2; 
    111111} else { 
    112     Sq = outer_sum/(4.0*M_PI); 
     112    Zq = outer_sum/(4.0*M_PI); 
    113113} 
    114114 
     115    return Zq; 
    115116    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    116  
    117     return _sc_volume_fraction(radius, dnn) * Pq * Sq; 
     117    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    118118} 
    119119 
    … …  
    133133    q = sqrt(qa*qa + qb*qb + qc*qc); 
    134134    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    135     const double Sq = _sq_sc(qa, qb, qc, dnn, d_factor); 
    136     return _sc_volume_fraction(radius, dnn) * Pq * Sq; 
     135    const double Zq = sc_Zq(qa, qb, qc, dnn, d_factor); 
     136    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    137137} 

  • sasmodels/models/bcc_paracrystal.c

    r2a0b2b1 r7e0b281  
    11static double 
    2 _sq_bcc(double qa, double qb, double qc, double dnn, double d_factor) 
     2bcc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    4     // Rewriting equations for efficiency, accuracy and readability, and so 
    5     // code is reusable between 1D and 2D models. 
    6     const double a1 = +qa - qc + qb; 
    7     const double a2 = +qa + qc - qb; 
    8     const double a3 = -qa + qc + qb; 
    9  
    10     const double half_dnn = 0.5*dnn; 
    11     const double arg = 0.5*square(half_dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     4#if 0  // Equations as written in Matsuoka 
     5    const double a1 = (+qa + qb + qc)/2.0; 
     6    const double a2 = (-qa - qb + qc)/2.0; 
     7    const double a3 = (-qa + qb - qc)/2.0; 
     8#else 
     9    const double a1 = (+qa + qb - qc)/2.0; 
     10    const double a2 = (+qa - qb + qc)/2.0; 
     11    const double a3 = (-qa + qb + qc)/2.0; 
     12#endif 
    1213 
    1314#if 1 
    … …  
    1516    //         => (-(exp(2a) - 1))^3 
    1617    //         => -expm1(2a)^3 
    17     // Denominator: prod(1 - 2 cos(xk) exp(a) + exp(a)^2) 
    18     //         => exp(a)^2 - 2 cos(xk) exp(a) + 1 
    19     //         => (exp(a) - 2 cos(xk)) * exp(a) + 1 
    20     const double exp_arg = exp(-arg); 
    21     const double Sq = -cube(expm1(-2.0*arg)) 
    22         / ( ((exp_arg - 2.0*cos(half_dnn*a1))*exp_arg + 1.0) 
    23           * ((exp_arg - 2.0*cos(half_dnn*a2))*exp_arg + 1.0) 
    24           * ((exp_arg - 2.0*cos(half_dnn*a3))*exp_arg + 1.0)); 
     18    // Denominator: prod(1 - 2 cos(d ak) exp(a) + exp(2a)) 
     19    //         => prod(exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
     20    //         => prod((exp(a) - 2 cos(d ak)) * exp(a) + 1) 
     21    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     22    const double exp_arg = exp(arg); 
     23    const double Zq = -cube(expm1(2.0*arg)) 
     24        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
     25          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
     26          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2527#else 
    2628    // Alternate form, which perhaps is more approachable 
     29    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    2730    const double sinh_qd = sinh(arg); 
    2831    const double cosh_qd = cosh(arg); 
    29     const double Sq = sinh_qd/(cosh_qd - cos(half_dnn*a1)) 
    30                     * sinh_qd/(cosh_qd - cos(half_dnn*a2)) 
    31                     * sinh_qd/(cosh_qd - cos(half_dnn*a3)); 
     32    const double Zq = sinh_qd/(cosh_qd - cos(dnn*a1)) 
     33                    * sinh_qd/(cosh_qd - cos(dnn*a2)) 
     34                    * sinh_qd/(cosh_qd - cos(dnn*a3)); 
    3235#endif 
    3336 
    34     return Sq; 
     37    return Zq; 
    3538} 
    3639 
    … …  
    3841// occupied volume fraction calculated from lattice symmetry and sphere radius 
    3942static double 
    40 _bcc_volume_fraction(double radius, double dnn) 
     43bcc_volume_fraction(double radius, double dnn) 
    4144{ 
    4245    return 2.0*sphere_volume(sqrt(0.75)*radius/dnn); 
    … …  
    7578            const double qa = qab*cos_phi; 
    7679            const double qb = qab*sin_phi; 
    77             const double fq = _sq_bcc(qa, qb, qc, dnn, d_factor); 
    78             inner_sum += Gauss150Wt[j] * fq; 
     80            const double form = bcc_Zq(qa, qb, qc, dnn, d_factor); 
     81            inner_sum += Gauss150Wt[j] * form; 
    7982        } 
    8083        inner_sum *= phi_m;  // sum(f(x)dx) = sum(f(x)) dx 
    … …  
    8285    } 
    8386    outer_sum *= theta_m; 
    84     const double Sq = outer_sum/(4.0*M_PI); 
     87    const double Zq = outer_sum/(4.0*M_PI); 
    8588    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    86  
    87     return _bcc_volume_fraction(radius, dnn) * Pq * Sq; 
     89    return bcc_volume_fraction(radius, dnn) * Pq * Zq; 
    8890} 
    8991 
    … …  
    101103 
    102104    q = sqrt(qa*qa + qb*qb + qc*qc); 
     105    const double Zq = bcc_Zq(qa, qb, qc, dnn, d_factor); 
    103106    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    104     const double Sq = _sq_bcc(qa, qb, qc, dnn, d_factor); 
    105     return _bcc_volume_fraction(radius, dnn) * Pq * Sq; 
     107    return bcc_volume_fraction(radius, dnn) * Pq * Zq; 
    106108} 

  • sasmodels/models/fcc_paracrystal.c

    r2a0b2b1 r7e0b281  
    11static double 
    2 _sq_fcc(double qa, double qb, double qc, double dnn, double d_factor) 
     2fcc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    4     // Rewriting equations for efficiency, accuracy and readability, and so 
    5     // code is reusable between 1D and 2D models. 
    6     const double a1 = qb + qa; 
    7     const double a2 = qa + qc; 
    8     const double a3 = qb + qc; 
    9  
    10     const double half_dnn = 0.5*dnn; 
    11     const double arg = 0.5*square(half_dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     4#if 0  // Equations as written in Matsuoka 
     5    const double a1 = ( qa + qb)/2.0; 
     6    const double a2 = (-qa + qc)/2.0; 
     7    const double a3 = (-qa + qb)/2.0; 
     8#else 
     9    const double a1 = ( qa + qb)/2.0; 
     10    const double a2 = ( qa + qc)/2.0; 
     11    const double a3 = ( qb + qc)/2.0; 
     12#endif 
    1213 
    1314    // Numerator: (1 - exp(a)^2)^3 
    1415    //         => (-(exp(2a) - 1))^3 
    1516    //         => -expm1(2a)^3 
    16     // Denominator: prod(1 - 2 cos(xk) exp(a) + exp(a)^2) 
    17     //         => exp(a)^2 - 2 cos(xk) exp(a) + 1 
    18     //         => (exp(a) - 2 cos(xk)) * exp(a) + 1 
    19     const double exp_arg = exp(-arg); 
    20     const double Sq = -cube(expm1(-2.0*arg)) 
    21         / ( ((exp_arg - 2.0*cos(half_dnn*a1))*exp_arg + 1.0) 
    22           * ((exp_arg - 2.0*cos(half_dnn*a2))*exp_arg + 1.0) 
    23           * ((exp_arg - 2.0*cos(half_dnn*a3))*exp_arg + 1.0)); 
     17    // Denominator: prod(1 - 2 cos(d ak) exp(a) + exp(2a)) 
     18    //         => prod(exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
     19    //         => prod((exp(a) - 2 cos(d ak)) * exp(a) + 1) 
     20    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     21    const double exp_arg = exp(arg); 
     22    const double Zq = -cube(expm1(2.0*arg)) 
     23        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
     24          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
     25          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2426 
    25     return Sq; 
     27    return Zq; 
    2628} 
    2729 
    … …  
    2931// occupied volume fraction calculated from lattice symmetry and sphere radius 
    3032static double 
    31 _fcc_volume_fraction(double radius, double dnn) 
     33fcc_volume_fraction(double radius, double dnn) 
    3234{ 
    3335    return 4.0*sphere_volume(M_SQRT1_2*radius/dnn); 
    … …  
    6668            const double qa = qab*cos_phi; 
    6769            const double qb = qab*sin_phi; 
    68             const double fq = _sq_fcc(qa, qb, qc, dnn, d_factor); 
    69             inner_sum += Gauss150Wt[j] * fq; 
     70            const double form = fcc_Zq(qa, qb, qc, dnn, d_factor); 
     71            inner_sum += Gauss150Wt[j] * form; 
    7072        } 
    7173        inner_sum *= phi_m;  // sum(f(x)dx) = sum(f(x)) dx 
    … …  
    7375    } 
    7476    outer_sum *= theta_m; 
    75     const double Sq = outer_sum/(4.0*M_PI); 
     77    const double Zq = outer_sum/(4.0*M_PI); 
    7678    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    7779 
    78     return _fcc_volume_fraction(radius, dnn) * Pq * Sq; 
     80    return fcc_volume_fraction(radius, dnn) * Pq * Zq; 
    7981} 
    8082 
    … …  
    9395    q = sqrt(qa*qa + qb*qb + qc*qc); 
    9496    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    95     const double Sq = _sq_fcc(qa, qb, qc, dnn, d_factor); 
    96     return _fcc_volume_fraction(radius, dnn) * Pq * Sq; 
     97    const double Zq = fcc_Zq(qa, qb, qc, dnn, d_factor); 
     98    return fcc_volume_fraction(radius, dnn) * Pq * Zq; 
    9799} 

  • sasmodels/models/sc_paracrystal.c

    r2a0b2b1 r7e0b281  
    11static double 
    2 _sq_sc(double qa, double qb, double qc, double dnn, double d_factor) 
     2sc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    4     // Rewriting equations for efficiency, accuracy and readability, and so 
    5     // code is reusable between 1D and 2D models. 
    64    const double a1 = qa; 
    75    const double a2 = qb; 
    86    const double a3 = qc; 
    97 
    10     const double arg = 0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    11  
    128    // Numerator: (1 - exp(a)^2)^3 
    139    //         => (-(exp(2a) - 1))^3 
    1410    //         => -expm1(2a)^3 
    15     // Denominator: prod(1 - 2 cos(xk) exp(a) + exp(a)^2) 
    16     //         => exp(a)^2 - 2 cos(xk) exp(a) + 1 
    17     //         => (exp(a) - 2 cos(xk)) * exp(a) + 1 
    18     const double exp_arg = exp(-arg); 
    19     const double Sq = -cube(expm1(-2.0*arg)) 
     11    // Denominator: prod(1 - 2 cos(d ak) exp(a) + exp(2a)) 
     12    //         => prod(exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
     13    //         => prod((exp(a) - 2 cos(d ak)) * exp(a) + 1) 
     14    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
     15    const double exp_arg = exp(arg); 
     16    const double Zq = -cube(expm1(2.0*arg)) 
    2017        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
    2118          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
    2219          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2320 
    24     return Sq; 
     21    return Zq; 
    2522} 
    2623 
    2724// occupied volume fraction calculated from lattice symmetry and sphere radius 
    2825static double 
    29 _sc_volume_fraction(double radius, double dnn) 
     26sc_volume_fraction(double radius, double dnn) 
    3027{ 
    3128    return sphere_volume(radius/dnn); 
    … …  
    6562            const double qa = qab*cos_phi; 
    6663            const double qb = qab*sin_phi; 
    67             const double fq = _sq_sc(qa, qb, qc, dnn, d_factor); 
    68             inner_sum += Gauss150Wt[j] * fq; 
     64            const double form = sc_Zq(qa, qb, qc, dnn, d_factor); 
     65            inner_sum += Gauss150Wt[j] * form; 
    6966        } 
    7067        inner_sum *= phi_m;  // sum(f(x)dx) = sum(f(x)) dx 
    … …  
    7269    } 
    7370    outer_sum *= theta_m; 
    74     const double Sq = outer_sum/M_PI_2; 
     71    const double Zq = outer_sum/M_PI_2; 
    7572    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    7673 
    77     return _sc_volume_fraction(radius, dnn) * Pq * Sq; 
     74    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    7875} 
    7976 
    … …  
    9289    q = sqrt(qa*qa + qb*qb + qc*qc); 
    9390    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    94     const double Sq = _sq_sc(qa, qb, qc, dnn, d_factor); 
    95     return _sc_volume_fraction(radius, dnn) * Pq * Sq; 
     91    const double Zq = sc_Zq(qa, qb, qc, dnn, d_factor); 
     92    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    9693} 
Note: See TracChangeset for help on using the changeset viewer.