Changeset 0ce5710 in sasmodels for sasmodels/models


Ignore:
Timestamp:
Apr 14, 2016 8:32:31 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:
b151003
Parents:
a5b8477 (diff), 62cf915 (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 remote-tracking branch 'origin/master' into polydisp

Location:
sasmodels/models
Files:
23 edited

Legend:

Unmodified
Added
Removed
  • sasmodels/models/adsorbed_layer.py

    rec45c4f r62cf915  
    33#converted by Steve King, Mar 2016 
    44 
     5r""" 
     6This model describes the scattering from a layer of surfactant or polymer 
     7adsorbed on large, smooth, notionally spherical particles under the conditions 
     8that (i) the particles (cores) are contrast-matched to the dispersion medium, 
     9(ii) $S(Q) \sim 1$ (ie, the particle volume fraction is dilute), (iii) the 
     10particle radius is >> layer thickness (ie, the interface is locally flat), 
     11and (iv) scattering from excess unadsorbed adsorbate in the bulk medium is 
     12absent or has been corrected for. 
    513 
    6  
    7 r""" 
    8 This model describes the scattering from a layer of surfactant or polymer adsorbed on large, smooth, notionally spherical particles under the conditions that (i) the particles (cores) are contrast-matched to the dispersion medium, (ii) *S(Q)* ~ 1 (ie, the particle volume fraction is dilute), (iii) the particle radius is >> layer thickness (ie, the interface is locally flat), and (iv) scattering from excess unadsorbed adsorbate in the bulk medium is absent or has been corrected for. 
    9  
    10 Unlike many other core-shell models, this model does not assume any form for the density distribution of the adsorbed species normal to the interface (cf, a core-shell model normally assumes the density distribution to be a homogeneous step-function). For comparison, if the thickness of a (traditional core-shell like) step function distribution is *t*, the second moment about the mean of the density distribution (ie, the distance of the centre-of-mass of the distribution from the interface), |sigma| =  sqrt((*t* :sup:`2` )/12). 
     14Unlike many other core-shell models, this model does not assume any form 
     15for the density distribution of the adsorbed species normal to the interface 
     16(cf, a core-shell model normally assumes the density distribution to be a 
     17homogeneous step-function). For comparison, if the thickness of a (traditional 
     18core-shell like) step function distribution is $t$, the second moment about 
     19the mean of the density distribution (ie, the distance of the centre-of-mass 
     20of the distribution from the interface), $\sigma = \sqrt{t^2/12}$. 
    1121 
    1222Definition 
     
    1525.. math:: 
    1626 
    17      I(q) = \text{scale} \cdot(\rho_\text{poly}-\rho_\text{solvent})^2    \left[\frac{6\pi\phi_\text{core}}{Q^2}\frac{\Gamma^2}{\delta_\text{poly}^2R_\text{core}} \exp(-Q^2\sigma^2)\right] + \text{background} 
     27     I(q) = \text{scale} \cdot (\rho_\text{poly}-\rho_\text{solvent})^2 
     28         \left[ 
     29             \frac{6\pi\phi_\text{core}}{Q^2} 
     30             \frac{\Gamma^2}{\delta_\text{poly}^2R_\text{core}} 
     31             \exp(-Q^2\sigma^2) 
     32         \right] + \text{background} 
    1833 
    19 where *scale* is a scale factor, |rho|\ :sub:`poly` is the sld of the polymer (or surfactant) layer, |rho|\ :sub:`solv` is the sld of the solvent/medium and cores, |phi|\ :sub:`core` is the volume fraction of the core particles, |delta|\ :sub:`poly` is the bulk density of the polymer, |biggamma| is the adsorbed amount, and |sigma| is the second moment of the thickness distribution. 
     34where *scale* is a scale factor, $\rho_\text{poly}$ is the sld of the 
     35polymer (or surfactant) layer, $\rho_\text{solv}$ is the sld of the 
     36solvent/medium and cores, $\phi_\text{core}$ is the volume fraction of 
     37the core particles, $\delta_\text{poly}$ is the bulk density of the 
     38polymer, $\Gamma$ is the adsorbed amount, and $\sigma$ is the second 
     39moment of the thickness distribution. 
    2040 
    21 Note that all parameters except the |sigma| are correlated so fitting more than one of these parameters will generally fail. Also note that unlike other shape models, no volume normalization is applied to this model (the calculation is exact). 
     41Note that all parameters except $\sigma$ are correlated so fitting more 
     42than one of these parameters will generally fail. Also note that unlike 
     43other shape models, no volume normalization is applied to this model (the 
     44calculation is exact). 
    2245 
    2346References 
    2447---------- 
    2548 
    26 S King, P Griffiths, J Hone, and T Cosgrove, *SANS from Adsorbed Polymer Layers*, 
    27 *Macromol. Symp.*, 190 (2002) 33-42. 
     49S King, P Griffiths, J Hone, and T Cosgrove, 
     50*SANS from Adsorbed Polymer Layers*, *Macromol. Symp.*, 190 (2002) 33-42. 
    2851""" 
    2952 
    3053from numpy import inf, sqrt, pi, exp 
    3154 
    32 name =  "adsorbed_layer" 
    33 title =  "Scattering from an adsorbed layer on particles" 
     55name = "adsorbed_layer" 
     56title = "Scattering from an adsorbed layer on particles" 
    3457 
    35 description =  """ 
     58description = """ 
    3659    Evaluates the scattering from large particles 
    3760    with an adsorbed layer of surfactant or 
     
    3962    density distribution. 
    4063    """ 
    41 category =  "shape:sphere" 
     64category = "shape:sphere" 
    4265 
     66# pylint: disable=bad-whitespace, line-too-long 
    4367#             ["name", "units", default, [lower, upper], "type", "description"], 
    44 parameters =  [["second_moment", "Ang", 23.0, [0.0, inf], "", "Second moment of polymer distribution"], 
    45               ["adsorbed_amount", "mg/m2", 1.9, [0.0, inf], "", "Adsorbed amount of polymer"], 
    46               ["density_shell", "g/cm3", 0.7, [0.0, inf], "", "Bulk density of polymer in the shell"], 
    47               ["radius", "Ang", 500.0, [0.0, inf], "", "Core particle radius"], 
    48               ["volfraction", "None", 0.14, [0.0, inf], "", "Core particle volume fraction"], 
    49               ["sld_shell", "1e-6/Ang^2", 1.5, [-inf, inf], "sld", "Polymer shell SLD"], 
    50               ["sld_solvent", "1e-6/Ang^2", 6.3, [-inf, inf], "sld", "Solvent SLD"]] 
     68parameters = [ 
     69    ["second_moment", "Ang", 23.0, [0.0, inf], "", "Second moment of polymer distribution"], 
     70    ["adsorbed_amount", "mg/m2", 1.9, [0.0, inf], "", "Adsorbed amount of polymer"], 
     71    ["density_shell", "g/cm3", 0.7, [0.0, inf], "", "Bulk density of polymer in the shell"], 
     72    ["radius", "Ang", 500.0, [0.0, inf], "", "Core particle radius"], 
     73    ["volfraction", "None", 0.14, [0.0, inf], "", "Core particle volume fraction"], 
     74    ["sld_shell", "1e-6/Ang^2", 1.5, [-inf, inf], "sld", "Polymer shell SLD"], 
     75    ["sld_solvent", "1e-6/Ang^2", 6.3, [-inf, inf], "sld", "Solvent SLD"], 
     76] 
     77# pylint: disable=bad-whitespace, line-too-long 
    5178 
    5279# NB: Scale and Background are implicit parameters on every model 
    53 def Iq(q, second_moment, adsorbed_amount, density_shell, radius,  
    54         volfraction, sld_shell, sld_solvent): 
     80def Iq(q, second_moment, adsorbed_amount, density_shell, radius, 
     81       volfraction, sld_shell, sld_solvent): 
    5582    # pylint: disable = missing-docstring 
    56 #    deltarhosqrd =  (sld_shell - sld_solvent) * (sld_shell - sld_solvent) 
    57 #    numerator =  6.0 * pi * volfraction * (adsorbed_amount * adsorbed_amount) 
    58 #    denominator =  (q * q) * (density_shell * density_shell) * radius 
    59 #    eterm =  exp(-1.0 * (q * q) * (second_moment * second_moment)) 
    60 #    #scale by 10^-2 for units conversion to cm^-1 
    61 #    inten =  1.0e-02 * deltarhosqrd * ((numerator / denominator) * eterm) 
    62     aa =  (sld_shell - sld_solvent) * adsorbed_amount / q / density_shell  
     83    #deltarhosqrd =  (sld_shell - sld_solvent) * (sld_shell - sld_solvent) 
     84    #numerator =  6.0 * pi * volfraction * (adsorbed_amount * adsorbed_amount) 
     85    #denominator =  (q * q) * (density_shell * density_shell) * radius 
     86    #eterm =  exp(-1.0 * (q * q) * (second_moment * second_moment)) 
     87    ##scale by 10^-2 for units conversion to cm^-1 
     88    #inten =  1.0e-02 * deltarhosqrd * ((numerator / denominator) * eterm) 
     89    aa = (sld_shell - sld_solvent) * adsorbed_amount / q / density_shell 
    6390    bb = q * second_moment 
    6491    #scale by 10^-2 for units conversion to cm^-1 
    65     inten =  6.0e-02 * pi * volfraction * aa * aa * exp(-bb * bb) / radius 
     92    inten = 6.0e-02 * pi * volfraction * aa * aa * exp(-bb * bb) / radius 
    6693    return inten 
    6794Iq.vectorized =  True  # Iq accepts an array of q values 
     
    7097    # pylint: disable = missing-docstring 
    7198    return Iq(sqrt(qx ** 2 + qy ** 2), *args) 
    72 Iqxy.vectorized =  True # Iqxy accepts an array of qx, qy values 
     99Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values 
    73100 
    74 demo =  dict(scale = 1.0, 
    75             second_moment = 23.0, 
    76             adsorbed_amount = 1.9, 
    77             density_shell = 0.7, 
    78             radius = 500.0, 
    79             volfraction = 0.14, 
    80             sld_shell = 1.5, 
    81             sld_solvent = 6.3, 
    82             background = 0.0) 
     101# unit test values taken from SasView 3.1.2 
     102tests =  [ 
     103    [{'scale': 1.0, 'second_moment': 23.0, 'adsorbed_amount': 1.9, 
     104      'density_shell': 0.7, 'radius': 500.0, 'volfraction': 0.14, 
     105      'sld_shell': 1.5, 'sld_solvent': 6.3, 'background': 0.0}, 
     106     [0.0106939, 0.1], [73.741, 4.51684e-3]], 
     107] 
    83108 
    84 # these unit test values taken from SasView 3.1.2 
    85 tests =  [ 
    86     [{'scale': 1.0, 'second_moment': 23.0, 'adsorbed_amount': 1.9,  
    87      'density_shell': 0.7, 'radius': 500.0, 'volfraction': 0.14,  
    88      'sld_shell': 1.5, 'sld_solvent': 6.3, 'background': 0.0}, 
    89      [0.0106939, 0.1], [73.741, 4.51684e-3]], 
    90     ] 
    91 # ADDED by: SMK  ON: 16Mar2016  convert from sasview, check vs SANDRA, 18Mar2016 RKH some edits & renaming 
     109# 2016-03-16 SMK converted from sasview, checked vs SANDRA 
     110# 2016-03-18 RKH some edits & renaming 
     111# 2016-04-14 PAK reformatting 
  • 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.c

    rfdb1487 rce896fd  
    44    double thickness, double A) 
    55{ 
    6   const double vol = 4.0/3.0 * M_PI * r * r * r; 
     6  const double vol = M_4PI_3 * cube(r); 
    77  const double qr = q * r; 
    88  const double alpha = A * r/thickness; 
     
    1919    double sinqr, cosqr; 
    2020    SINCOS(qr, sinqr, cosqr); 
    21     fun = -3.0( 
     21    fun = -3.0*( 
    2222            ((alphasq - qrsq)*sinqr/qr - 2.0*alpha*cosqr) / sumsq 
    2323                - (alpha*sinqr/qr - cosqr) 
     
    3232f_linear(double q, double r, double sld, double slope) 
    3333{ 
    34   const double vol = 4.0/3.0 * M_PI * r * r * r; 
     34  const double vol = M_4PI_3 * cube(r); 
    3535  const double qr = q * r; 
    3636  const double bes = sph_j1c(qr); 
     
    5252{ 
    5353  const double bes = sph_j1c(q * r); 
    54   const double vol = 4.0/3.0 * M_PI * r * r * r; 
     54  const double vol = M_4PI_3 * cube(r); 
    5555  return sld * vol * bes; 
    5656} 
     
    6464    r += thickness[i]; 
    6565  } 
    66   return 4.0/3.0 * M_PI * r * r * r; 
     66  return M_4PI_3*cube(r); 
    6767} 
    6868 
    6969static double 
    70 Iq(double q, double core_sld, double core_radius, double solvent_sld, 
    71     double n, double in_sld[], double out_sld[], double thickness[], 
     70Iq(double q, double sld_core, double core_radius, double sld_solvent, 
     71    double n, double sld_in[], double sld_out[], double thickness[], 
    7272    double A[]) 
    7373{ 
    7474  int i; 
    75   r = core_radius; 
    76   f = f_constant(q, r, core_sld); 
     75  double r = core_radius; 
     76  double f = f_constant(q, r, sld_core); 
    7777  for (i=0; i<n; i++){ 
    7878    const double r0 = r; 
     
    9292    } 
    9393  } 
    94   f -= f_constant(q, r, solvent_sld); 
    95   f2 = f * f * 1.0e-4; 
     94  f -= f_constant(q, r, sld_solvent); 
     95  const double f2 = f * f * 1.0e-4; 
    9696 
    9797  return f2; 
  • 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.c

    r13ed84c rd2bb604  
    11double Iq(double q, double case_num, 
    2     double Na, double Phia, double va, double a_sld, double ba, 
    3     double Nb, double Phib, double vb, double b_sld, double bb, 
    4     double Nc, double Phic, double vc, double c_sld, double bc, 
    5     double Nd, double Phid, double vd, double d_sld, double bd, 
     2    double N[], double Phi[], double v[], double L[], double b[], 
    63    double Kab, double Kac, double Kad, 
    74    double Kbc, double Kbd, double Kcd 
    85    ); 
    96 
    10 double Iqxy(double qx, double qy, double case_num, 
    11     double Na, double Phia, double va, double a_sld, double ba, 
    12     double Nb, double Phib, double vb, double b_sld, double bb, 
    13     double Nc, double Phic, double vc, double c_sld, double bc, 
    14     double Nd, double Phid, double vd, double d_sld, double bd, 
    15     double Kab, double Kac, double Kad, 
    16     double Kbc, double Kbd, double Kcd 
    17     ); 
    18  
    19 double form_volume(void); 
    20  
    21 double form_volume(void) 
    22 { 
    23     return 1.0; 
    24 } 
    25  
    267double Iq(double q, double case_num, 
    27     double Na, double Phia, double va, double La, double ba, 
    28     double Nb, double Phib, double vb, double Lb, double bb, 
    29     double Nc, double Phic, double vc, double Lc, double bc, 
    30     double Nd, double Phid, double vd, double Ld, double bd, 
     8    double N[], double Phi[], double v[], double L[], double b[], 
    319    double Kab, double Kac, double Kad, 
    3210    double Kbc, double Kbd, double Kcd 
     
    3614#if 0  // Sasview defaults 
    3715  if (icase <= 1) { 
    38     Na=Nb=1000.0; 
    39     Phia=Phib=0.0000001; 
     16    N[0]=N[1]=1000.0; 
     17    Phi[0]=Phi[1]=0.0000001; 
    4018    Kab=Kac=Kad=Kbc=Kbd=-0.0004; 
    41     La=Lb=1e-12; 
    42     va=vb=100.0; 
    43     ba=bb=5.0; 
     19    L[0]=L[1]=1e-12; 
     20    v[0]=v[1]=100.0; 
     21    b[0]=b[1]=5.0; 
    4422  } else if (icase <= 4) { 
    45     Phia=0.0000001; 
     23    Phi[0]=0.0000001; 
    4624    Kab=Kac=Kad=-0.0004; 
    47     La=1e-12; 
    48     va=100.0; 
    49     ba=5.0; 
     25    L[0]=1e-12; 
     26    v[0]=100.0; 
     27    b[0]=5.0; 
    5028  } 
    5129#else 
    5230  if (icase <= 1) { 
    53     Na=Nb=0.0; 
    54     Phia=Phib=0.0; 
     31    N[0]=N[1]=0.0; 
     32    Phi[0]=Phi[1]=0.0; 
    5533    Kab=Kac=Kad=Kbc=Kbd=0.0; 
    56     La=Lb=Ld; 
    57     va=vb=vd; 
    58     ba=bb=0.0; 
     34    L[0]=L[1]=L[3]; 
     35    v[0]=v[1]=v[3]; 
     36    b[0]=b[1]=0.0; 
    5937  } else if (icase <= 4) { 
    60     Na = 0.0; 
    61     Phia=0.0; 
     38    N[0] = 0.0; 
     39    Phi[0]=0.0; 
    6240    Kab=Kac=Kad=0.0; 
    63     La=Ld; 
    64     va=vd; 
    65     ba=0.0; 
     41    L[0]=L[3]; 
     42    v[0]=v[3]; 
     43    b[0]=0.0; 
    6644  } 
    6745#endif 
    6846 
    69   const double Xa = q*q*ba*ba*Na/6.0; 
    70   const double Xb = q*q*bb*bb*Nb/6.0; 
    71   const double Xc = q*q*bc*bc*Nc/6.0; 
    72   const double Xd = q*q*bd*bd*Nd/6.0; 
     47  const double Xa = q*q*b[0]*b[0]*N[0]/6.0; 
     48  const double Xb = q*q*b[1]*b[1]*N[1]/6.0; 
     49  const double Xc = q*q*b[2]*b[2]*N[2]/6.0; 
     50  const double Xd = q*q*b[3]*b[3]*N[3]/6.0; 
    7351 
    7452  // limit as Xa goes to 0 is 1 
     
    9876#if 0 
    9977  const double S0aa = icase<5 
    100                       ? 1.0 : Na*Phia*va*Paa; 
     78                      ? 1.0 : N[0]*Phi[0]*v[0]*Paa; 
    10179  const double S0bb = icase<2 
    102                       ? 1.0 : Nb*Phib*vb*Pbb; 
    103   const double S0cc = Nc*Phic*vc*Pcc; 
    104   const double S0dd = Nd*Phid*vd*Pdd; 
     80                      ? 1.0 : N[1]*Phi[1]*v[1]*Pbb; 
     81  const double S0cc = N[2]*Phi[2]*v[2]*Pcc; 
     82  const double S0dd = N[3]*Phi[3]*v[3]*Pdd; 
    10583  const double S0ab = icase<8 
    106                       ? 0.0 : sqrt(Na*va*Phia*Nb*vb*Phib)*Pa*Pb; 
     84                      ? 0.0 : sqrt(N[0]*v[0]*Phi[0]*N[1]*v[1]*Phi[1])*Pa*Pb; 
    10785  const double S0ac = icase<9 
    108                       ? 0.0 : sqrt(Na*va*Phia*Nc*vc*Phic)*Pa*Pc*exp(-Xb); 
     86                      ? 0.0 : sqrt(N[0]*v[0]*Phi[0]*N[2]*v[2]*Phi[2])*Pa*Pc*exp(-Xb); 
    10987  const double S0ad = icase<9 
    110                       ? 0.0 : sqrt(Na*va*Phia*Nd*vd*Phid)*Pa*Pd*exp(-Xb-Xc); 
     88                      ? 0.0 : sqrt(N[0]*v[0]*Phi[0]*N[3]*v[3]*Phi[3])*Pa*Pd*exp(-Xb-Xc); 
    11189  const double S0bc = (icase!=4 && icase!=7 && icase!= 9) 
    112                       ? 0.0 : sqrt(Nb*vb*Phib*Nc*vc*Phic)*Pb*Pc; 
     90                      ? 0.0 : sqrt(N[1]*v[1]*Phi[1]*N[2]*v[2]*Phi[2])*Pb*Pc; 
    11391  const double S0bd = (icase!=4 && icase!=7 && icase!= 9) 
    114                       ? 0.0 : sqrt(Nb*vb*Phib*Nd*vd*Phid)*Pb*Pd*exp(-Xc); 
     92                      ? 0.0 : sqrt(N[1]*v[1]*Phi[1]*N[3]*v[3]*Phi[3])*Pb*Pd*exp(-Xc); 
    11593  const double S0cd = (icase==0 || icase==2 || icase==5) 
    116                       ? 0.0 : sqrt(Nc*vc*Phic*Nd*vd*Phid)*Pc*Pd; 
     94                      ? 0.0 : sqrt(N[2]*v[2]*Phi[2]*N[3]*v[3]*Phi[3])*Pc*Pd; 
    11795#else  // sasview equivalent 
    118 //printf("Xc=%g, S0cc=%g*%g*%g*%g\n",Xc,Nc,Phic,vc,Pcc); 
    119   double S0aa = Na*Phia*va*Paa; 
    120   double S0bb = Nb*Phib*vb*Pbb; 
    121   double S0cc = Nc*Phic*vc*Pcc; 
    122   double S0dd = Nd*Phid*vd*Pdd; 
    123   double S0ab = sqrt(Na*va*Phia*Nb*vb*Phib)*Pa*Pb; 
    124   double S0ac = sqrt(Na*va*Phia*Nc*vc*Phic)*Pa*Pc*exp(-Xb); 
    125   double S0ad = sqrt(Na*va*Phia*Nd*vd*Phid)*Pa*Pd*exp(-Xb-Xc); 
    126   double S0bc = sqrt(Nb*vb*Phib*Nc*vc*Phic)*Pb*Pc; 
    127   double S0bd = sqrt(Nb*vb*Phib*Nd*vd*Phid)*Pb*Pd*exp(-Xc); 
    128   double S0cd = sqrt(Nc*vc*Phic*Nd*vd*Phid)*Pc*Pd; 
     96//printf("Xc=%g, S0cc=%g*%g*%g*%g\n",Xc,N[2],Phi[2],v[2],Pcc); 
     97  double S0aa = N[0]*Phi[0]*v[0]*Paa; 
     98  double S0bb = N[1]*Phi[1]*v[1]*Pbb; 
     99  double S0cc = N[2]*Phi[2]*v[2]*Pcc; 
     100  double S0dd = N[3]*Phi[3]*v[3]*Pdd; 
     101  double S0ab = sqrt(N[0]*v[0]*Phi[0]*N[1]*v[1]*Phi[1])*Pa*Pb; 
     102  double S0ac = sqrt(N[0]*v[0]*Phi[0]*N[2]*v[2]*Phi[2])*Pa*Pc*exp(-Xb); 
     103  double S0ad = sqrt(N[0]*v[0]*Phi[0]*N[3]*v[3]*Phi[3])*Pa*Pd*exp(-Xb-Xc); 
     104  double S0bc = sqrt(N[1]*v[1]*Phi[1]*N[2]*v[2]*Phi[2])*Pb*Pc; 
     105  double S0bd = sqrt(N[1]*v[1]*Phi[1]*N[3]*v[3]*Phi[3])*Pb*Pd*exp(-Xc); 
     106  double S0cd = sqrt(N[2]*v[2]*Phi[2]*N[3]*v[3]*Phi[3])*Pc*Pd; 
    129107switch(icase){ 
    130108  case 0: 
     
    311289  // Note: 1e-13 to convert from fm to cm for scattering length 
    312290  const double sqrt_Nav=sqrt(6.022045e+23) * 1.0e-13; 
    313   const double Lad = icase<5 ? 0.0 : (La/va - Ld/vd)*sqrt_Nav; 
    314   const double Lbd = icase<2 ? 0.0 : (Lb/vb - Ld/vd)*sqrt_Nav; 
    315   const double Lcd = (Lc/vc - Ld/vd)*sqrt_Nav; 
     291  const double Lad = icase<5 ? 0.0 : (L[0]/v[0] - L[3]/v[3])*sqrt_Nav; 
     292  const double Lbd = icase<2 ? 0.0 : (L[1]/v[1] - L[3]/v[3])*sqrt_Nav; 
     293  const double Lcd = (L[2]/v[2] - L[3]/v[3])*sqrt_Nav; 
    316294 
    317295  const double result=Lad*Lad*S11 + Lbd*Lbd*S22 + Lcd*Lcd*S33 
     
    321299 
    322300} 
    323  
    324 double Iqxy(double qx, double qy, 
    325     double case_num, 
    326     double Na, double Phia, double va, double a_sld, double ba, 
    327     double Nb, double Phib, double vb, double b_sld, double bb, 
    328     double Nc, double Phic, double vc, double c_sld, double bc, 
    329     double Nd, double Phid, double vd, double d_sld, double bd, 
    330     double Kab, double Kac, double Kad, 
    331     double Kbc, double Kbd, double Kcd 
    332     ) 
    333 { 
    334     double q = sqrt(qx*qx + qy*qy); 
    335     return Iq(q, 
    336         case_num, 
    337         Na, Phia, va, a_sld, ba, 
    338         Nb, Phib, vb, b_sld, bb, 
    339         Nc, Phic, vc, c_sld, bc, 
    340         Nd, Phid, vd, d_sld, bd, 
    341         Kab, Kac, Kad, 
    342         Kbc, Kbd, Kcd); 
    343 } 
  • 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 rd2bb604  
    170170# pylint: disable=bad-whitespace, line-too-long 
    171171#            ["name", "units", default, [lower, upper], "type", "description"], 
    172 parameters = [["n_shells",         "",               1,      [0, 9],         "", "number of shells"], 
     172parameters = [["n",                "",               1,      [0, 9],         "", "number of shells"], 
    173173              ["radius_core",      "Ang",            50.0,   [0, inf],       "", "intern layer thickness"], 
    174174              ["sld_core",         "1e-6/Ang^2",     2.07,   [-inf, inf],    "", "sld function flat"], 
     
    192192 
    193193demo = dict( 
    194     n_shells=4, 
     194    n=4, 
    195195    scale=1.0, 
    196196    solvent_sld=1.0, 
  • 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.