source: sasmodels/sasmodels/models/sphere.py @ db3947c

ticket_1156ticket_822_more_unit_tests
Last change on this file since db3947c was db3947c, checked in by richardh, 5 years ago

more comprehensive unit tests in sphere, but still issues with beta approx

  • Property mode set to 100644
File size: 7.1 KB
Line 
1r"""
2For information about polarised and magnetic scattering, see
3the :ref:`magnetism` documentation.
4
5Definition
6----------
7
8The 1D scattering intensity is calculated in the following way (Guinier, 1955)
9
10.. math::
11
12    I(q) = \frac{\text{scale}}{V} \cdot \left[
13        3V(\Delta\rho) \cdot \frac{\sin(qr) - qr\cos(qr))}{(qr)^3}
14        \right]^2 + \text{background}
15
16where *scale* is a volume fraction, $V$ is the volume of the scatterer,
17$r$ is the radius of the sphere and *background* is the background level.
18*sld* and *sld_solvent* are the scattering length densities (SLDs) of the
19scatterer and the solvent respectively, whose difference is $\Delta\rho$.
20
21Note that if your data is in absolute scale, the *scale* should represent
22the volume fraction (which is unitless) if you have a good fit. If not,
23it should represent the volume fraction times a factor (by which your data
24might need to be rescaled).
25
26The 2D scattering intensity is the same as above, regardless of the
27orientation of $\vec q$.
28
29Validation
30----------
31
32Validation of our code was done by comparing the output of the 1D model
33to the output of the software provided by the NIST (Kline, 2006).
34
35
36References
37----------
38
39.. [#] A Guinier and G. Fournet, *Small-Angle Scattering of X-Rays*, John Wiley and Sons, New York, (1955)
40
41Source
42------
43
44`sphere.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/sphere.py>`_
45
46`sphere.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/sphere.c>`_
47
48Authorship and Verification
49----------------------------
50
51* **Author:**
52* **Last Modified by:**
53* **Last Reviewed by:** S King and P Parker **Date:** 2013/09/09 and 2014/01/06
54* **Source added by :** Steve King **Date:** March 25, 2019
55"""
56
57import numpy as np
58from numpy import inf
59
60name = "sphere"
61title = "Spheres with uniform scattering length density"
62description = """\
63P(q)=(scale/V)*[3V(sld-sld_solvent)*(sin(qr)-qr cos(qr))
64                /(qr)^3]^2 + background
65    r: radius of sphere
66    V: The volume of the scatter
67    sld: the SLD of the sphere
68    sld_solvent: the SLD of the solvent
69"""
70category = "shape:sphere"
71
72#             ["name", "units", default, [lower, upper], "type","description"],
73parameters = [["sld", "1e-6/Ang^2", 1, [-inf, inf], "sld",
74               "Layer scattering length density"],
75              ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "sld",
76               "Solvent scattering length density"],
77              ["radius", "Ang", 50, [0, inf], "volume",
78               "Sphere radius"],
79             ]
80
81source = ["lib/sas_3j1x_x.c", "sphere.c"]
82have_Fq = True
83radius_effective_modes = ["radius"]
84
85def random():
86    """Return a random parameter set for the model."""
87    radius = 10**np.random.uniform(1.3, 4)
88    pars = dict(
89        radius=radius,
90    )
91    return pars
92
93tests = [
94     [{}, 0.2, 0.726362], # each test starts with default parameter values inside { }, unless modified. Then Q and expected value of I(Q)
95     [{"scale": 1., "background": 0., "sld": 6., "sld_solvent": 1.,
96       "radius": 120.}, [0.01,0.1,0.2], [1.34836265e+04, 6.20114062e+00, 1.04733914e-01]], # each test starts with default parameter values inside { }, unless modified. Then Q and expected value of I(Q)
97     [{"scale": 1., "background": 0., "sld": 6., "sld_solvent": 1.,  # careful tests here R=120 Pd=.2,
98      #                                        then with S(Q) at default Reff=50 (but this gets changeded to 120) phi=0,2
99       "radius": 120., "radius_pd": 0.2, "radius_pd_n":45},
100      [0.01,0.1,0.2], [1.74395295e+04, 3.68016987e+00, 2.28843099e-01]],  # a list of Q values and list of expected results is also possible
101    [{"scale": 1., "background": 0., "sld": 6., "sld_solvent": 1.,"radius": 120., "radius_pd": 0.2, "radius_pd_n":45},
102      0.01, 335839.88055473, 1.41045057e+11, 120.0, 8087664.122641933, 1.0], # the longer list here checks  F1, F2, R_eff, volume, volume_ratio = call_Fq(kernel, pars)
103    [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45},
104      0.1, 482.93824329, 29763977.79867414, 120.0, 8087664.122641933, 1.0], # the longer list here checks  F1, F2, R_eff, volume, volume_ratio = call_Fq(kernel, pars)
105    [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45},
106      0.2, 1.23330406, 1850806.1197361, 120.0, 8087664.122641933, 1.0], # the longer list here checks  F1, F2, R_eff, volume, volume_ratio = call_Fq(kernel, pars)
107   #  But note P(Q) = F2/volume
108   #  F and F^2 are "unscaled", with for  n <F F*>S(q) or for beta approx I(q) = n [<F F*> + <F><F*> (S(q) - 1)]
109   #  for n the number density and <.> the orientation average, and F = integral rho(r) exp(i q . r) dr.
110   #  The number density is volume fraction divided by particle volume.
111   #  Effectively, this leaves F = V drho form, where form is the usual 3 j1(qr)/(qr) or whatever depending on the shape.
112   # NOTE the @S multiplication by S(Q) also multiplies the answer by volfraction, thus you may like to put in scale at 1/volfraction
113    [{"@S": "hardsphere",
114     "radius": 120., "radius_pd": 0.2, "radius_pd_n":45,
115     "volfraction":0.2,
116     "radius_effective":45.0,     # uses this as gives different answer to either 50 or 120 (check)
117     "structure_factor_mode": 1,  # 0 = normal decoupling approximation, 1 = beta(Q) approx
118     "radius_effective_mode": 0   # equivalent sphere, there is only one valid mode for sphere. says -this used r_eff =0 or default 50?
119     }, 0.01, 1316.2990966463444 ],
120    [{"@S": "hardsphere",
121     "radius": 120., "radius_pd": 0.2, "radius_pd_n":45,
122     "volfraction":0.2,
123     #"radius_effective":50.0,        # hard sphere structure factor
124     "structure_factor_mode": 1,  # 0 = normal decoupling approximation, 1 = beta(Q) approx
125     "radius_effective_mode": 0   # this used default 50?
126     }, [0.01,0.1,0.2], [1.32473756e+03, 7.36633631e-01, 4.67686201e-02]  ],
127    [{"@S": "hardsphere",
128     "radius": 120., "radius_pd": 0.2, "radius_pd_n":45,
129     "volfraction":0.2,
130     "radius_effective":120.0,        # hard sphere structure factor
131     "structure_factor_mode": 1,  # 0 = normal decoupling approximation, 1 = beta(Q) approx
132     "radius_effective_mode": 0   # 1 uses 120,
133     }, [0.01,0.1,0.2], [1.57928589e+03, 7.37067923e-01, 4.67686197e-02  ]],
134    [{"@S": "hardsphere",
135     "radius": 120., "radius_pd": 0.2, "radius_pd_n":45,
136     "volfraction":0.2,
137     #"radius_effective":120.0,        # hard sphere structure factor
138     "structure_factor_mode": 0,  # 0 = normal decoupling approximation, 1 = beta(Q) approx
139     "radius_effective_mode": 1   # this used 120 from the form factor
140     }, [0.01,0.1,0.2], [1.10112335e+03, 7.41366536e-01, 4.66630207e-02]],
141    [{"@S": "hardsphere",
142     "radius": 120., "radius_pd": 0.2, "radius_pd_n":45,
143     "volfraction":0.2,
144     #"radius_effective":50.0,        # hard sphere structure factor
145     "structure_factor_mode": 0,  # 0 = normal decoupling approximation, 1 = beta(Q) approx
146     "radius_effective_mode": 0   # this used 50 the default for hardsphere
147     }, [0.01,0.1,0.2], [7.82803598e+02, 6.85943611e-01, 4.71586457e-02 ]]
148]# putting None for expected result will pass the test if there are no errors from the routine, but without any check on the value of the result
Note: See TracBrowser for help on using the repository browser.