Changeset 639c4e3 in sasmodels for sasmodels/models


Ignore:
Timestamp:
Apr 26, 2016 9:43:53 AM (8 years ago)
Author:
Paul Kienzle <pkienzle@…>
Branches:
master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
ed246ab
Parents:
cebbb5a (diff), fa227d3 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
Message:

Merge branch 'master' into polydisp

Conflicts:

sasmodels/kerneldll.py
sasmodels/models/rpa.c

Location:
sasmodels/models
Files:
3 added
25 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/core_multi_shell.c

    rf7930be rabdd01c  
    182182    if (r == r0) { 
    183183      // no thickness, so nothing to add 
    184     } else if (fabs(A[i]) < 1e-16 || sld_out[i] == sld_in[i]) { 
     184    } else if (fabs(A[i]) < 1.0e-16 || sld_out[i] == sld_in[i]) { 
    185185      f -= f_constant(q, r0, sld_in[i]); 
    186186      f += f_constant(q, r, sld_in[i]); 
    187     } else if (fabs(A[i]) < 1e-4) { 
     187    } else if (fabs(A[i]) < 1.0e-4) { 
    188188      const double slope = (sld_out[i] - sld_in[i])/thickness[i]; 
    189189      f -= f_linear(q, r0, sld_in[i], slope); 

  • sasmodels/models/core_shell_parallelepiped.c

    r44bd2be rabdd01c  
    4848    double a_scaled = a_side / b_side; 
    4949    double c_scaled = c_side / b_side; 
    50     double arim_scaled = arim_thickness / b_side; 
    51     double brim_scaled = brim_thickness / b_side; 
    52         
     50 
    5351    // DelRho values (note that drC is not used later)        
    5452        double dr0 = core_sld-solvent_sld; 

  • sasmodels/models/elliptical_cylinder.c

    r43b7eea rabdd01c  
    8787 
    8888    const double vol = form_volume(r_minor, r_ratio, length); 
    89     return answer*vol*vol*1e-4; 
     89    return answer*vol*vol*1.0e-4; 
    9090} 
    9191 

  • sasmodels/models/onion.c

    rce896fd r639c4e3  
    8080    if (r == r0) { 
    8181      // no thickness, so nothing to add 
    82     } else if (fabs(A[i]) < 1e-16 || sld_out[i] == sld_in[i]) { 
     82    } else if (fabs(A[i]) < 1.0e-16 || sld_out[i] == sld_in[i]) { 
    8383      f -= f_constant(q, r0, sld_in[i]); 
    8484      f += f_constant(q, r, sld_in[i]); 
    85     } else if (fabs(A[i]) < 1e-4) { 
     85    } else if (fabs(A[i]) < 1.0e-4) { 
    8686      const double slope = (sld_out[i] - sld_in[i])/thickness[i]; 
    8787      f -= f_linear(q, r0, sld_in[i], slope); 

  • sasmodels/models/rpa.c

    rd2bb604 r639c4e3  
    1717    Phi[0]=Phi[1]=0.0000001; 
    1818    Kab=Kac=Kad=Kbc=Kbd=-0.0004; 
    19     L[0]=L[1]=1e-12; 
    20     v[0]=v[1]=100.0; 
    21     b[0]=b[1]=5.0; 
     19    La=Lb=1.0e-12; 
     20    va=vb=100.0; 
     21    ba=bb=5.0; 
    2222  } else if (icase <= 4) { 
    2323    Phi[0]=0.0000001; 
    2424    Kab=Kac=Kad=-0.0004; 
    25     L[0]=1e-12; 
    26     v[0]=100.0; 
    27     b[0]=5.0; 
     25    La=1.0e-12; 
     26    va=100.0; 
     27    ba=5.0; 
    2828  } 
    2929#else 

  • sasmodels/models/cylinder.c

    r26141cb re9b1663d  
    33double Iqxy(double qx, double qy, double sld, double solvent_sld, 
    44    double radius, double length, double theta, double phi); 
     5 
     6#define INVALID(v) (v.radius<0 || v.length<0) 
    57 
    68double form_volume(double radius, double length) 
    … …  
    1517    double length) 
    1618{ 
    17     // TODO: return NaN if radius<0 or length<0? 
    1819    // precompute qr and qh to save time in the loop 
    1920    const double qr = q*radius; 
    … …  
    4748    double phi) 
    4849{ 
    49     // TODO: return NaN if radius<0 or length<0? 
    5050    double sn, cn; // slots to hold sincos function output 
    5151 

  • sasmodels/models/cylinder.py

    rf247314 r7ae2b7f  
    8282""" 
    8383 
    84 import numpy as np 
    85 from numpy import pi, inf 
     84import numpy as np  # type: ignore 
     85from numpy import pi, inf  # type: ignore 
    8686 
    8787name = "cylinder" 

  • sasmodels/models/flexible_cylinder.c

    re6408d0 r4937980  
    1 double form_volume(double length, double kuhn_length, double radius); 
    2 double Iq(double q, double length, double kuhn_length, double radius, 
    3           double sld, double solvent_sld); 
    4 double Iqxy(double qx, double qy, double length, double kuhn_length, 
    5             double radius, double sld, double solvent_sld); 
    6 double flexible_cylinder_kernel(double q, double length, double kuhn_length, 
    7                                 double radius, double sld, double solvent_sld); 
    8  
    9  
    10 double form_volume(double length, double kuhn_length, double radius) 
     1static double 
     2form_volume(length, kuhn_length, radius) 
    113{ 
    124    return 1.0; 
    135} 
    146 
    15 double flexible_cylinder_kernel(double q, 
    16           double length, 
    17           double kuhn_length, 
    18           double radius, 
    19           double sld, 
    20           double solvent_sld) 
     7static double 
     8Iq(double q, 
     9   double length, 
     10   double kuhn_length, 
     11   double radius, 
     12   double sld, 
     13   double solvent_sld) 
    2114{ 
    22  
    23     const double cont = sld-solvent_sld; 
    24     const double qr = q*radius; 
    25     //const double crossSect = (2.0*J1(qr)/qr)*(2.0*J1(qr)/qr); 
    26     const double crossSect = sas_J1c(qr); 
    27     double flex = Sk_WR(q,length,kuhn_length); 
    28     flex *= crossSect*crossSect; 
    29     flex *= M_PI*radius*radius*length; 
    30     flex *= cont*cont; 
    31     flex *= 1.0e-4; 
    32     return flex; 
     15    const double contrast = sld - solvent_sld; 
     16    const double cross_section = sas_J1c(q*radius); 
     17    const double volume = M_PI*radius*radius*length; 
     18    const double flex = Sk_WR(q, length, kuhn_length); 
     19    return 1.0e-4 * volume * square(contrast*cross_section) * flex; 
    3320} 
    34  
    35 double Iq(double q, 
    36           double length, 
    37           double kuhn_length, 
    38           double radius, 
    39           double sld, 
    40           double solvent_sld) 
    41 { 
    42  
    43     double result = flexible_cylinder_kernel(q, length, kuhn_length, radius, sld, solvent_sld); 
    44     return result; 
    45 } 
    46  
    47 double Iqxy(double qx, double qy, 
    48             double length, 
    49             double kuhn_length, 
    50             double radius, 
    51             double sld, 
    52             double solvent_sld) 
    53 { 
    54     double q; 
    55     q = sqrt(qx*qx+qy*qy); 
    56     double result = flexible_cylinder_kernel(q, length, kuhn_length, radius, sld, solvent_sld); 
    57  
    58     return result; 
    59 } 

  • sasmodels/models/gel_fit.c

    r30b4ddf r03cac08  
    1 double form_volume(void); 
    2  
    3 double Iq(double q, 
    4           double guinier_scale, 
    5           double lorentzian_scale, 
    6           double gyration_radius, 
    7           double fractal_exp, 
    8           double cor_length); 
    9  
    10 double Iqxy(double qx, double qy, 
    11           double guinier_scale, 
    12           double lorentzian_scale, 
    13           double gyration_radius, 
    14           double fractal_exp, 
    15           double cor_length); 
    16  
    17 static double _gel_fit_kernel(double q, 
     1static double Iq(double q, 
    182          double guinier_scale, 
    193          double lorentzian_scale, 
    … …  
    248    // Lorentzian Term 
    259    ////////////////////////double a(x[i]*x[i]*zeta*zeta); 
    26     double lorentzian_term = q*q*cor_length*cor_length; 
     10    double lorentzian_term = square(q*cor_length); 
    2711    lorentzian_term = 1.0 + ((fractal_exp + 1.0)/3.0)*lorentzian_term; 
    2812    lorentzian_term = pow(lorentzian_term, (fractal_exp/2.0) ); 
    … …  
    3014    // Exponential Term 
    3115    ////////////////////////double d(x[i]*x[i]*rg*rg); 
    32     double exp_term = q*q*gyration_radius*gyration_radius; 
     16    double exp_term = square(q*gyration_radius); 
    3317    exp_term = exp(-1.0*(exp_term/3.0)); 
    3418 
    … …  
    3721    return result; 
    3822} 
    39 double form_volume(void){ 
    40     // Unused, so free to return garbage. 
    41     return NAN; 
    42 } 
    43  
    44 double Iq(double q, 
    45           double guinier_scale, 
    46           double lorentzian_scale, 
    47           double gyration_radius, 
    48           double fractal_exp, 
    49           double cor_length) 
    50 { 
    51     return _gel_fit_kernel(q, 
    52                           guinier_scale, 
    53                           lorentzian_scale, 
    54                           gyration_radius, 
    55                           fractal_exp, 
    56                           cor_length); 
    57 } 
    58  
    59 // Iqxy is never called since no orientation or magnetic parameters. 
    60 double Iqxy(double qx, double qy, 
    61           double guinier_scale, 
    62           double lorentzian_scale, 
    63           double gyration_radius, 
    64           double fractal_exp, 
    65           double cor_length) 
    66 { 
    67     double q = sqrt(qx*qx + qy*qy); 
    68     return _gel_fit_kernel(q, 
    69                           guinier_scale, 
    70                           lorentzian_scale, 
    71                           gyration_radius, 
    72                           fractal_exp, 
    73                           cor_length); 
    74 } 
    75  

  • sasmodels/models/hardsphere.py

    rec45c4f rd2bb604  
    149149   """ 
    150150 
    151 Iqxy = """ 
    152     // never called since no orientation or magnetic parameters. 
    153     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    154     """ 
    155  
    156151# ER defaults to 0.0 
    157152# VR defaults to 1.0 

  • sasmodels/models/hayter_msa.py

    rec45c4f rd2bb604  
    8787    return 1.0; 
    8888    """ 
    89 Iqxy = """ 
    90     // never called since no orientation or magnetic parameters. 
    91     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    92     """ 
    9389# ER defaults to 0.0 
    9490# VR defaults to 1.0 

  • sasmodels/models/lamellar.py

    rec45c4f rd2bb604  
    8282    """ 
    8383 
    84 Iqxy = """ 
    85     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    86     """ 
    87  
    8884# ER defaults to 0.0 
    8985# VR defaults to 1.0 

  • sasmodels/models/lamellar_hg.py

    rec45c4f rd2bb604  
    101101    """ 
    102102 
    103 Iqxy = """ 
    104     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    105     """ 
    106  
    107103# ER defaults to 0.0 
    108104# VR defaults to 1.0 

  • sasmodels/models/lamellar_hg_stack_caille.py

    rec45c4f rd2bb604  
    120120    """ 
    121121 
    122 Iqxy = """ 
    123     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    124     """ 
    125  
    126122# ER defaults to 0.0 
    127123# VR defaults to 1.0 

  • sasmodels/models/lamellar_stack_caille.py

    rec45c4f rd2bb604  
    104104    """ 
    105105 
    106 Iqxy = """ 
    107     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    108     """ 
    109  
    110106# ER defaults to 0.0 
    111107# VR defaults to 1.0 

  • sasmodels/models/lamellar_stack_paracrystal.py

    rec45c4f rd2bb604  
    132132    """ 
    133133 
    134 Iqxy = """ 
    135     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    136     """ 
    137  
    138134# ER defaults to 0.0 
    139135# VR defaults to 1.0 

  • sasmodels/models/lib/sas_JN.c

    re6408d0 r4937980  
    4848*/ 
    4949 
    50 static double 
    51 sas_JN( int n, double x ) { 
     50double sas_JN( int n, double x ); 
     51 
     52double sas_JN( int n, double x ) { 
    5253 
    5354    const double MACHEP = 1.11022302462515654042E-16; 

  • sasmodels/models/lib/sph_j1c.c

    re6f1410 rba32cdd  
    77* using double precision that are the source. 
    88*/ 
     9double sph_j1c(double q); 
    910 
    1011// The choice of the number of terms in the series and the cutoff value for 
    … …  
    4344#endif 
    4445 
    45 double sph_j1c(double q); 
    4646double sph_j1c(double q) 
    4747{ 

  • sasmodels/models/lib/sphere_form.c

    rad90df9 rba32cdd  
    1 inline double 
    2 sphere_volume(double radius) 
     1double sphere_volume(double radius); 
     2double sphere_form(double q, double radius, double sld, double solvent_sld); 
     3 
     4double sphere_volume(double radius) 
    35{ 
    46    return M_4PI_3*cube(radius); 
    57} 
    68 
    7 inline double 
    8 sphere_form(double q, double radius, double sld, double solvent_sld) 
     9double sphere_form(double q, double radius, double sld, double solvent_sld) 
    910{ 
    1011    const double fq = sphere_volume(radius) * sph_j1c(q*radius); 

  • sasmodels/models/lib/wrc_cyl.c

    re7678b2 rba32cdd  
    22    Functions for WRC implementation of flexible cylinders 
    33*/ 
     4double Sk_WR(double q, double L, double b); 
     5 
     6 
    47static double 
    58AlphaSquare(double x) 
    … …  
    363366} 
    364367 
    365 double Sk_WR(double q, double L, double b); 
    366368double Sk_WR(double q, double L, double b) 
    367369{ 

  • sasmodels/models/onion.py

    r416609b rfa5fd8d  
    293293 
    294294#             ["name", "units", default, [lower, upper], "type","description"], 
    295 parameters = [["core_sld", "1e-6/Ang^2", 1.0, [-inf, inf], "", 
     295parameters = [["sld_core", "1e-6/Ang^2", 1.0, [-inf, inf], "", 
    296296               "Core scattering length density"], 
    297297              ["core_radius", "Ang", 200., [0, inf], "volume", 
    298298               "Radius of the core"], 
    299               ["solvent_sld", "1e-6/Ang^2", 6.4, [-inf, inf], "", 
     299              ["sld_solvent", "1e-6/Ang^2", 6.4, [-inf, inf], "", 
    300300               "Solvent scattering length density"], 
    301               ["n", "", 1, [0, 10], "volume", 
     301              ["n_shells", "", 1, [0, 10], "volume", 
    302302               "number of shells"], 
    303               ["in_sld[n]", "1e-6/Ang^2", 1.7, [-inf, inf], "", 
     303              ["sld_in[n_shells]", "1e-6/Ang^2", 1.7, [-inf, inf], "", 
    304304               "scattering length density at the inner radius of shell k"], 
    305               ["out_sld[n]", "1e-6/Ang^2", 2.0, [-inf, inf], "", 
     305              ["sld_out[n_shells]", "1e-6/Ang^2", 2.0, [-inf, inf], "", 
    306306               "scattering length density at the outer radius of shell k"], 
    307               ["thickness[n]", "Ang", 40., [0, inf], "volume", 
     307              ["thickness[n_shells]", "Ang", 40., [0, inf], "volume", 
    308308               "Thickness of shell k"], 
    309               ["A[n]", "", 1.0, [-inf, inf], "", 
     309              ["A[n_shells]", "", 1.0, [-inf, inf], "", 
    310310               "Decay rate of shell k"], 
    311311              ] 
    312312 
    313 #source = ["lib/sph_j1c.c", "onion.c"] 
    314  
    315 def Iq(q, *args, **kw): 
    316     return q 
    317  
    318 def Iqxy(qx, *args, **kw): 
    319     return qx 
    320  
    321  
    322 def shape(core_sld, core_radius, solvent_sld, n, in_sld, out_sld, thickness, A): 
     313source = ["lib/sph_j1c.c", "onion.c"] 
     314 
     315#def Iq(q, *args, **kw): 
     316#    return q 
     317 
     318profile_axes = ['Radius (A)', 'SLD (1e-6/A^2)'] 
     319def profile(core_sld, core_radius, solvent_sld, n_shells, 
     320            in_sld, out_sld, thickness, A): 
    323321    """ 
    324     SLD profile 
     322    Returns shape profile with x=radius, y=SLD. 
    325323    """ 
    326324 
    327     total_radius = 1.25*(sum(thickness[:n]) + core_radius + 1) 
     325    total_radius = 1.25*(sum(thickness[:n_shells]) + core_radius + 1) 
    328326    dr = total_radius/400  # 400 points for a smooth plot 
    329327 
    … …  
    338336 
    339337    # add in the shells 
    340     for k in range(n): 
     338    for k in range(n_shells): 
    341339        # Left side of each shells 
    342340        r0 = r[-1] 
    … …  
    365363    beta.append(solvent_sld) 
    366364 
    367     return np.asarray(r), np.asarray(beta) 
     365    return np.asarray(r), np.asarray(beta)*1e-6 
    368366 
    369367def ER(core_radius, n, thickness): 
    … …  
    374372 
    375373demo = { 
    376     "solvent_sld": 2.2, 
    377     "core_sld": 1.0, 
     374    "sld_solvent": 2.2, 
     375    "sld_core": 1.0, 
    378376    "core_radius": 100, 
    379377    "n": 4, 
    380     "in_sld": [0.5, 1.5, 0.9, 2.0], 
    381     "out_sld": [nan, 0.9, 1.2, 1.6], 
     378    "sld_in": [0.5, 1.5, 0.9, 2.0], 
     379    "sld_out": [nan, 0.9, 1.2, 1.6], 
    382380    "thickness": [50, 75, 150, 75], 
    383381    "A": [0, -1, 1e-4, 1], 

  • sasmodels/models/rpa.py

    rec45c4f ra5b8477  
    8686#   ["name", "units", default, [lower, upper], "type","description"], 
    8787parameters = [ 
    88     ["case_num", CASES, 0, [0, 10], "", "Component organization"], 
     88    ["case_num", "", 1, [CASES], "", "Component organization"], 
    8989 
    90     ["Na", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
    91     ["Phia", "", 0.25, [0, 1], "", "volume fraction"], 
    92     ["va", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
    93     ["La", "fm", 10.0, [-inf, inf], "", "scattering length"], 
    94     ["ba", "Ang", 5.0, [0, inf], "", "segment length"], 
    95  
    96     ["Nb", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
    97     ["Phib", "", 0.25, [0, 1], "", "volume fraction"], 
    98     ["vb", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
    99     ["Lb", "fm", 10.0, [-inf, inf], "", "scattering length"], 
    100     ["bb", "Ang", 5.0, [0, inf], "", "segment length"], 
    101  
    102     ["Nc", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
    103     ["Phic", "", 0.25, [0, 1], "", "volume fraction"], 
    104     ["vc", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
    105     ["Lc", "fm", 10.0, [-inf, inf], "", "scattering length"], 
    106     ["bc", "Ang", 5.0, [0, inf], "", "segment length"], 
    107  
    108     ["Nd", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
    109     ["Phid", "", 0.25, [0, 1], "", "volume fraction"], 
    110     ["vd", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
    111     ["Ld", "fm", 10.0, [-inf, inf], "", "scattering length"], 
    112     ["bd", "Ang", 5.0, [0, inf], "", "segment length"], 
     90    ["N[4]", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
     91    ["Phi[4]", "", 0.25, [0, 1], "", "volume fraction"], 
     92    ["v[4]", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
     93    ["L[4]", "fm", 10.0, [-inf, inf], "", "scattering length"], 
     94    ["b[4]", "Ang", 5.0, [0, inf], "", "segment length"], 
    11395 
    11496    ["Kab", "", -0.0004, [-inf, inf], "", "Interaction parameter"], 

  • sasmodels/models/spherical_sld.py

    rec45c4f rdd71228  
    11r""" 
    2 This model calculates an empirical functional form for SAS data using SpericalSLD profile 
    3  
    4 Similarly to the OnionExpShellModel, this model provides the form factor, P(q), for a multi-shell sphere, 
    5 where the interface between the each neighboring shells can be described by one of a number of functions 
    6 including error, power-law, and exponential functions. This model is to calculate the scattering intensity 
    7 by building a continuous custom SLD profile against the radius of the particle. The SLD profile is composed 
    8 of a flat core, a flat solvent, a number (up to 9 ) flat shells, and the interfacial layers between 
    9 the adjacent flat shells (or core, and solvent) (see below). 
     2This model calculates an empirical functional form for SAS data using 
     3SpericalSLD profile 
     4 
     5Similarly to the OnionExpShellModel, this model provides the form factor, 
     6P(q), for a multi-shell sphere, where the interface between the each neighboring 
     7shells can be described by one of a number of functions including error, 
     8power-law, and exponential functions. 
     9This model is to calculate the scattering intensity by building a continuous 
     10custom SLD profile against the radius of the particle. 
     11The SLD profile is composed of a flat core, a flat solvent, a number (up to 9 ) 
     12flat shells, and the interfacial layers between the adjacent flat shells 
     13(or core, and solvent) (see below). 
    1014 
    1115.. figure:: img/spherical_sld_profile.gif 
    … …  
    1317    Exemplary SLD profile 
    1418 
    15 Unlike the <onion> model (using an analytical integration), 
    16 the interfacial layers here are sub-divided and numerically integrated assuming each of the sub-layers are described 
    17 by a line function. The number of the sub-layer can be given by users by setting the integer values of npts_inter. 
    18 The form factor is normalized by the total volume of the sphere. 
     19Unlike the <onion> model (using an analytical integration), the interfacial 
     20layers here are sub-divided and numerically integrated assuming each of the 
     21sub-layers are described by a line function. 
     22The number of the sub-layer can be given by users by setting the integer values 
     23of npts_inter. The form factor is normalized by the total volume of the sphere. 
    1924 
    2025Definition 
    … …  
    2934    \sum_{\text{flat}_i=0}^N f_{\text{flat}_i} +f_\text{solvent} 
    3035 
    31 For a spherically symmetric particle with a particle density $\rho_x(r)$ the sld function can be defined as: 
     36For a spherically symmetric particle with a particle density $\rho_x(r)$ 
     37the sld function can be defined as: 
    3238 
    3339.. math:: 
    … …  
    3945 
    4046.. math:: 
    41     f_\text{core} = 4 \pi \int_{0}^{r_\text{core}} \rho_\text{core} \frac{\sin(qr)} {qr} r^2 dr = 
     47    f_\text{core} = 4 \pi \int_{0}^{r_\text{core}} \rho_\text{core} 
     48    \frac{\sin(qr)} {qr} r^2 dr = 
    4249    3 \rho_\text{core} V(r_\text{core}) 
    43     \Big[ \frac{\sin(qr_\text{core}) - qr_\text{core} \cos(qr_\text{core})} {qr_\text{core}^3} \Big] 
    44  
    45     f_{\text{inter}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } \rho_{ \text{inter}_i } \frac{\sin(qr)} {qr} r^2 dr 
    46  
    47     f_{\text{shell}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } \rho_{ \text{flat}_i } \frac{\sin(qr)} {qr} r^2 dr = 
    48     3 \rho_{ \text{flat}_i } V ( r_{ \text{inter}_i } + \Delta t_{ \text{inter}_i } ) 
    49     \Big[ \frac{\sin(qr_{\text{inter}_i} + \Delta t_{ \text{inter}_i } ) - q (r_{\text{inter}_i} + 
    50     \Delta t_{ \text{inter}_i }) \cos(q( r_{\text{inter}_i} + \Delta t_{ \text{inter}_i } ) ) } 
     50    \Big[ \frac{\sin(qr_\text{core}) - qr_\text{core} \cos(qr_\text{core})} 
     51    {qr_\text{core}^3} \Big] 
     52 
     53    f_{\text{inter}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } 
     54    \rho_{ \text{inter}_i } \frac{\sin(qr)} {qr} r^2 dr 
     55 
     56    f_{\text{shell}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } 
     57    \rho_{ \text{flat}_i } \frac{\sin(qr)} {qr} r^2 dr = 
     58    3 \rho_{ \text{flat}_i } V ( r_{ \text{inter}_i } + 
     59    \Delta t_{ \text{inter}_i } ) 
     60    \Big[ \frac{\sin(qr_{\text{inter}_i} + \Delta t_{ \text{inter}_i } ) 
     61    - q (r_{\text{inter}_i} + \Delta t_{ \text{inter}_i }) 
     62    \cos(q( r_{\text{inter}_i} + \Delta t_{ \text{inter}_i } ) ) } 
    5163    {q ( r_{\text{inter}_i} + \Delta t_{ \text{inter}_i } )^3 }  \Big] 
    5264    -3 \rho_{ \text{flat}_i } V(r_{ \text{inter}_i }) 
    53     \Big[ \frac{\sin(qr_{\text{inter}_i}) - qr_{\text{flat}_i} \cos(qr_{\text{inter}_i}) } {qr_{\text{inter}_i}^3} \Big] 
    54  
    55     f_\text{solvent} = 4 \pi \int_{r_N}^{\infty} \rho_\text{solvent} \frac{\sin(qr)} {qr} r^2 dr = 
     65    \Big[ \frac{\sin(qr_{\text{inter}_i}) - qr_{\text{flat}_i} 
     66    \cos(qr_{\text{inter}_i}) } {qr_{\text{inter}_i}^3} \Big] 
     67 
     68    f_\text{solvent} = 4 \pi \int_{r_N}^{\infty} \rho_\text{solvent} 
     69    \frac{\sin(qr)} {qr} r^2 dr = 
    5670    3 \rho_\text{solvent} V(r_N) 
    5771    \Big[ \frac{\sin(qr_N) - qr_N \cos(qr_N)} {qr_N^3} \Big] 
    … …  
    6680.. math:: 
    6781    \rho_{{inter}_i} (r) = \begin{cases} 
    68     B \exp\Big( \frac {\pm A(r - r_{\text{flat}_i})} {\Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A \neq 0 \\ 
    69     B \Big( \frac {(r - r_{\text{flat}_i})} {\Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A = 0 \\ 
     82    B \exp\Big( \frac {\pm A(r - r_{\text{flat}_i})} 
     83    {\Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A \neq 0 \\ 
     84    B \Big( \frac {(r - r_{\text{flat}_i})} 
     85    {\Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A = 0 \\ 
    7086    \end{cases} 
    7187 
    … …  
    7490.. math:: 
    7591    \rho_{{inter}_i} (r) = \begin{cases} 
    76     \pm B \Big( \frac {(r - r_{\text{flat}_i} )} {\Delta t_{ \text{inter}_i }} \Big) ^A  +C  & \text{for} A \neq 0 \\ 
     92    \pm B \Big( \frac {(r - r_{\text{flat}_i} )} {\Delta t_{ \text{inter}_i }} 
     93    \Big) ^A  +C  & \text{for} A \neq 0 \\ 
    7794    \rho_{\text{flat}_{i+1}}  & \text{for} A = 0 \\ 
    7895    \end{cases} 
    … …  
    8299.. math:: 
    83100    \rho_{{inter}_i} (r) = \begin{cases} 
    84     B \text{erf} \Big( \frac { A(r - r_{\text{flat}_i})} {\sqrt{2} \Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A \neq 0 \\ 
    85     B \Big( \frac {(r - r_{\text{flat}_i} )} {\Delta t_{ \text{inter}_i }} \Big)  +C  & \text{for} A = 0 \\ 
     101    B \text{erf} \Big( \frac { A(r - r_{\text{flat}_i})} 
     102    {\sqrt{2} \Delta t_{ \text{inter}_i }} \Big) +C  & \text{for} A \neq 0 \\ 
     103    B \Big( \frac {(r - r_{\text{flat}_i} )} {\Delta t_{ \text{inter}_i }} 
     104    \Big)  +C  & \text{for} A = 0 \\ 
    86105    \end{cases} 
    87106 
    88 The functions are normalized so that they vary between 0 and 1, and they are constrained such that the SLD 
    89 is continuous at the boundaries of the interface as well as each sub-layers. Thus B and C are determined. 
    90  
    91 Once $\rho_{\text{inter}_i}$ is found at the boundary of the sub-layer of the interface, we can find its contribution 
    92 to the form factor $P(q)$ 
    93  
    94 .. math:: 
    95     f_{\text{inter}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } \rho_{ \text{inter}_i } \frac{\sin(qr)} {qr} r^2 dr = 
     107The functions are normalized so that they vary between 0 and 1, and they are 
     108constrained such that the SLD is continuous at the boundaries of the interface 
     109as well as each sub-layers. Thus B and C are determined. 
     110 
     111Once $\rho_{\text{inter}_i}$ is found at the boundary of the sub-layer of the 
     112interface, we can find its contribution to the form factor $P(q)$ 
     113 
     114.. math:: 
     115    f_{\text{inter}_i} = 4 \pi \int_{\Delta t_{ \text{inter}_i } } 
     116    \rho_{ \text{inter}_i } \frac{\sin(qr)} {qr} r^2 dr = 
    96117    4 \pi \sum_{j=0}^{npts_{\text{inter}_i} -1 } 
    97     \int_{r_j}^{r_{j+1}} \rho_{ \text{inter}_i } (r_j) \frac{\sin(qr)} {qr} r^2 dr \approx 
     118    \int_{r_j}^{r_{j+1}} \rho_{ \text{inter}_i } (r_j) 
     119    \frac{\sin(qr)} {qr} r^2 dr \approx 
    98120 
    99121    4 \pi \sum_{j=0}^{npts_{\text{inter}_i} -1 } \Big[ 
    100     3 ( \rho_{ \text{inter}_i } ( r_{j+1} ) - \rho_{ \text{inter}_i } ( r_{j} ) V ( r_{ \text{sublayer}_j } ) 
    101     \Big[ \frac {r_j^2 \beta_\text{out}^2 \sin(\beta_\text{out}) - (\beta_\text{out}^2-2) \cos(\beta_\text{out}) } 
     122    3 ( \rho_{ \text{inter}_i } ( r_{j+1} ) - \rho_{ \text{inter}_i } 
     123    ( r_{j} ) V ( r_{ \text{sublayer}_j } ) 
     124    \Big[ \frac {r_j^2 \beta_\text{out}^2 \sin(\beta_\text{out}) 
     125    - (\beta_\text{out}^2-2) \cos(\beta_\text{out}) } 
    102126    {\beta_\text{out}^4 } \Big] 
    103127 
    104     - 3 ( \rho_{ \text{inter}_i } ( r_{j+1} ) - \rho_{ \text{inter}_i } ( r_{j} ) V ( r_{ \text{sublayer}_j-1 } ) 
    105     \Big[ \frac {r_{j-1}^2 \sin(\beta_\text{in}) - (\beta_\text{in}^2-2) \cos(\beta_\text{in}) } 
     128    - 3 ( \rho_{ \text{inter}_i } ( r_{j+1} ) - \rho_{ \text{inter}_i } 
     129    ( r_{j} ) V ( r_{ \text{sublayer}_j-1 } ) 
     130    \Big[ \frac {r_{j-1}^2 \sin(\beta_\text{in}) 
     131    - (\beta_\text{in}^2-2) \cos(\beta_\text{in}) } 
    106132    {\beta_\text{in}^4 } \Big] 
    107133 
    … …  
    120146    V(a) = \frac {4\pi}{3}a^3 
    121147 
    122     a_\text{in} ~ \frac{r_j}{r_{j+1} -r_j} \text{, } a_\text{out} ~ \frac{r_{j+1}}{r_{j+1} -r_j} 
     148    a_\text{in} ~ \frac{r_j}{r_{j+1} -r_j} \text{, } a_\text{out} 
     149    ~ \frac{r_{j+1}}{r_{j+1} -r_j} 
    123150 
    124151    \beta_\text{in} = qr_j \text{, } \beta_\text{out} = qr_{j+1} 
    125152 
    126153 
    127 We assume the $\rho_{\text{inter}_i} (r)$ can be approximately linear within a sub-layer $j$ 
     154We assume the $\rho_{\text{inter}_i} (r)$ can be approximately linear 
     155within a sub-layer $j$ 
    128156 
    129157Finally form factor can be calculated by 
    … …  
    131159.. math:: 
    132160 
    133     P(q) = \frac{[f]^2} {V_\text{particle}} \text{where} V_\text{particle} = V(r_{\text{shell}_N}) 
     161    P(q) = \frac{[f]^2} {V_\text{particle}} \text{where} V_\text{particle} 
     162    = V(r_{\text{shell}_N}) 
    134163 
    135164For 2D data the scattering intensity is calculated in the same way as 1D, 
    … …  
    150179 
    151180.. note:: 
    152     The outer most radius is used as the effective radius for S(Q) when $P(Q) * S(Q)$ is applied. 
     181    The outer most radius is used as the effective radius for S(Q) 
     182    when $P(Q) * S(Q)$ is applied. 
    153183 
    154184References 
    155185---------- 
    156 L A Feigin and D I Svergun, Structure Analysis by Small-Angle X-Ray and Neutron Scattering, Plenum Press, New York, (1987) 
     186L A Feigin and D I Svergun, Structure Analysis by Small-Angle X-Ray 
     187and Neutron Scattering, Plenum Press, New York, (1987) 
    157188 
    158189""" 
    … …  
    170201# pylint: disable=bad-whitespace, line-too-long 
    171202#            ["name", "units", default, [lower, upper], "type", "description"], 
    172 parameters = [["n_shells",         "",               1,      [0, 9],         "", "number of shells"], 
    173               ["radius_core",      "Ang",            50.0,   [0, inf],       "", "intern layer thickness"], 
    174               ["sld_core",         "1e-6/Ang^2",     2.07,   [-inf, inf],    "", "sld function flat"], 
    175               ["sld_flat[n]",      "1e-6/Ang^2",     4.06,   [-inf, inf],    "", "sld function flat"], 
    176               ["thick_flat[n]",    "Ang",            100.0,  [0, inf],       "", "flat layer_thickness"], 
    177               ["func_inter[n]",    "",               0,      [0, 4],         "", "Erf:0, RPower:1, LPower:2, RExp:3, LExp:4"], 
    178               ["thick_inter[n]",   "Ang",            50.0,   [0, inf],       "", "intern layer thickness"], 
    179               ["inter_nu[n]",      "",               2.5,    [-inf, inf],    "", "steepness parameter"], 
    180               ["npts_inter",       "",               35,     [0, 35],        "", "number of points in each sublayer Must be odd number"], 
    181               ["sld_solvent",      "1e-6/Ang^2",     1.0,    [-inf, inf],    "", "sld function solvent"], 
     203parameters = [["n_shells",          "",               1,      [0, 9],         "volume", "number of shells"], 
     204              ["npts_inter",        "",               35,     [0, inf],        "", "number of points in each sublayer Must be odd number"], 
     205              ["radius_core",       "Ang",            50.0,   [0, inf],       "volume", "intern layer thickness"], 
     206              ["sld_core",          "1e-6/Ang^2",     2.07,   [-inf, inf],    "", "sld function flat"], 
     207              ["sld_solvent",       "1e-6/Ang^2",     1.0,    [-inf, inf],    "", "sld function solvent"], 
     208              ["func_inter0",       "",               0,      [0, 4],         "", "Erf:0, RPower:1, LPower:2, RExp:3, LExp:4"], 
     209              ["thick_inter0",      "Ang",            50.0,   [0, inf],       "volume", "intern layer thickness for core layer"], 
     210              ["nu_inter0",         "",               2.5,    [-inf, inf],    "", "steepness parameter for core layer"], 
     211              ["sld_flat[n_shells]",      "1e-6/Ang^2",     4.06,   [-inf, inf],    "", "sld function flat"], 
     212              ["thick_flat[n_shells]",    "Ang",            100.0,  [0, inf],       "volume", "flat layer_thickness"], 
     213              ["func_inter[n_shells]",    "",               0,      [0, 4],         "", "Erf:0, RPower:1, LPower:2, RExp:3, LExp:4"], 
     214              ["thick_inter[n_shells]",   "Ang",            50.0,   [0, inf],       "volume", "intern layer thickness"], 
     215              ["nu_inter[n_shells]",      "",               2.5,    [-inf, inf],    "", "steepness parameter"], 
    182216              ] 
    183217# pylint: enable=bad-whitespace, line-too-long 
    184 #source = ["lib/librefl.c",  "lib/sph_j1c.c", "spherical_sld.c"] 
    185  
    186 def Iq(q, *args, **kw): 
    187     return q 
    188  
    189 def Iqxy(qx, *args, **kw): 
    190     return qx 
    191  
    192  
    193 demo = dict( 
    194     n_shells=4, 
    195     scale=1.0, 
    196     solvent_sld=1.0, 
    197     background=0.0, 
    198     npts_inter=35.0, 
    199     ) 
     218source = ["lib/librefl.c",  "lib/sph_j1c.c", "spherical_sld.c"] 
     219 
     220profile_axes = ['Radius (A)', 'SLD (1e-6/A^2)'] 
     221def profile(n_shells, radius_core,  sld_core,  sld_solvent, sld_flat, 
     222            thick_flat, func_inter, thick_inter, nu_inter, npts_inter): 
     223    """ 
     224    Returns shape profile with x=radius, y=SLD. 
     225    """ 
     226 
     227    z = [] 
     228    beta = [] 
     229    z0 = 0 
     230    # two sld points for core 
     231    z.append(0) 
     232    beta.append(sld_core) 
     233    z.append(radius_core) 
     234    beta.append(sld_core) 
     235    z0 += radius_core 
     236 
     237    for i in range(1, n_shells+2): 
     238        dz = thick_inter[i-1]/npts_inter 
     239        # j=0 for interface, j=1 for flat layer 
     240        for j in range(0, 2): 
     241            # interation for sub-layers 
     242            for n_s in range(0, npts_inter+1): 
     243                if j == 1: 
     244                    if i == n_shells+1: 
     245                        break 
     246                    # shift half sub thickness for the first point 
     247                    z0 -= dz#/2.0 
     248                    z.append(z0) 
     249                    #z0 -= dz/2.0 
     250                    z0 += thick_flat[i] 
     251                    sld_i = sld_flat[i] 
     252                    beta.append(sld_flat[i]) 
     253                    dz = 0 
     254                else: 
     255                    nu = nu_inter[i-1] 
     256                    # decide which sld is which, sld_r or sld_l 
     257                    if i == 1: 
     258                        sld_l = sld_core 
     259                    else: 
     260                        sld_l = sld_flat[i-1] 
     261                    if i == n_shells+1: 
     262                        sld_r = sld_solvent 
     263                    else: 
     264                        sld_r = sld_flat[i] 
     265                    # get function type 
     266                    func_idx = func_inter[i-1] 
     267                    # calculate the sld 
     268                    sld_i = intersldfunc(func_idx, npts_inter, n_s, nu, 
     269                                            sld_l, sld_r) 
     270                # append to the list 
     271                z.append(z0) 
     272                beta.append(sld_i) 
     273                z0 += dz 
     274                if j == 1: 
     275                    break 
     276    z.append(z0) 
     277    beta.append(sld_solvent) 
     278    z_ext = z0/5.0 
     279    z.append(z0+z_ext) 
     280    beta.append(sld_solvent) 
     281    # return sld profile (r, beta) 
     282    return np.asarray(z), np.asarray(beta)*1e-6 
     283 
     284def ER(n_shells, radius_core, thick_inter0, thick_inter, thick_flat): 
     285    total_thickness = thick_inter0 
     286    total_thickness += np.sum(thick_inter[:n_shells], axis=0) 
     287    total_thickness += np.sum(thick_flat[:n_shells], axis=0) 
     288    return total_thickness + radius_core 
     289 
     290 
     291demo = { 
     292    "n_shells": 4, 
     293    "npts_inter": 35.0, 
     294    "radius_core": 50.0, 
     295    "sld_core": 2.07, 
     296    "sld_solvent": 1.0, 
     297    "thick_inter0": 50.0, 
     298    "func_inter0": 0, 
     299    "nu_inter0": 2.5, 
     300    "sld_flat":[4.0,3.5,4.0,3.5], 
     301    "thick_flat":[100.0,100.0,100.0,100.0], 
     302    "func_inter":[0,0,0,0], 
     303    "thick_inter":[50.0,50.0,50.0,50.0], 
     304    "nu_inter":[2.5,2.5,2.5,2.5], 
     305    } 
    200306 
    201307#TODO: Not working yet 
    202308tests = [ 
    203309    # Accuracy tests based on content in test/utest_extra_models.py 
    204     [{'npts_iter':35, 
    205         'sld_solv':1, 
    206         'radius_core':50.0, 
    207         'sld_core':2.07, 
    208         'func_inter2':0.0, 
    209         'thick_inter2':50, 
    210         'nu_inter2':2.5, 
    211         'sld_flat2':4, 
    212         'thick_flat2':100, 
    213         'func_inter1':0.0, 
    214         'thick_inter1':50, 
    215         'nu_inter1':2.5, 
    216         'background': 0.0, 
     310    [{"n_shells":4, 
     311        'npts_inter':35, 
     312        "radius_core":50.0, 
     313        "sld_core":2.07, 
     314        "sld_solvent": 1.0, 
     315        "sld_flat":[4.0,3.5,4.0,3.5], 
     316        "thick_flat":[100.0,100.0,100.0,100.0], 
     317        "func_inter":[0,0,0,0], 
     318        "thick_inter":[50.0,50.0,50.0,50.0], 
     319        "nu_inter":[2.5,2.5,2.5,2.5] 
    217320    }, 0.001, 0.001], 
    218321] 

  • sasmodels/models/squarewell.py

    rec45c4f rd2bb604  
    123123""" 
    124124 
    125 Iqxy = """ 
    126     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    127     """ 
    128  
    129125# ER defaults to 0.0 
    130126# VR defaults to 1.0 

  • sasmodels/models/stickyhardsphere.py

    rec45c4f rd2bb604  
    171171""" 
    172172 
    173 Iqxy = """ 
    174     return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS); 
    175     """ 
    176  
    177173# ER defaults to 0.0 
    178174# VR defaults to 1.0 
Note: See TracChangeset for help on using the changeset viewer.