source: sasmodels/sasmodels/models/stacked_disks.py @ ef07e95

core_shell_microgelsmagnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since ef07e95 was ef07e95, checked in by Paul Kienzle <pkienzle@…>, 5 years ago

use consistent formatting for 'Last Reviewed'

  • Property mode set to 100644
File size: 10.7 KB
Line 
1r"""
2Definition
3----------
4
5This model provides the form factor, $P(q)$, for stacked discs (tactoids)
6with a core/layer structure which is constructed itself as $P(q) S(Q)$
7multiplying a $P(q)$ for individual core/layer disks by a structure factor
8$S(q)$ proposed by Kratky and Porod in 1949\ [#CIT1949]_ assuming the next
9neighbor distance (d-spacing) in the stack of parallel discs obeys a Gaussian
10distribution. As such the normalization of this "composite" form factor is
11relative to the individual disk volume, not the volume of the stack of disks.
12This model is appropriate for example for non non exfoliated clay particles
13such as Laponite.
14
15.. figure:: img/stacked_disks_geometry.png
16
17   Geometry of a single core/layer disk
18
19The scattered intensity $I(q)$ is calculated as
20
21.. math::
22
23    I(q) = N\int_{0}^{\pi /2}\left[ \Delta \rho_t
24    \left( V_t f_t(q,\alpha) - V_c f_c(q,\alpha)\right) + \Delta
25    \rho_c V_c f_c(q,\alpha)\right]^2 S(q,\alpha)\sin{\alpha}\ d\alpha
26    + \text{background}
27
28where the contrast
29
30.. math::
31
32    \Delta \rho_i = \rho_i - \rho_\text{solvent}
33
34and $N$ is the number of individual (single) discs per unit volume, $\alpha$
35is the angle between the axis of the disc and $q$, and $V_t$ and $V_c$ are the
36total volume and the core volume of a single disc, respectively, and
37
38.. math::
39
40    f_t(q,\alpha) =
41    \left(\frac{\sin(q(d+h)\cos{\alpha})}{q(d+h)\cos{\alpha}}\right)
42    \left(\frac{2J_1(qR\sin{\alpha})}{qR\sin{\alpha}} \right)
43
44    f_c(q,\alpha) =
45    \left(\frac{\sin(qh)\cos{\alpha})}{qh\cos{\alpha}}\right)
46    \left(\frac{2J_1(qR\sin{\alpha})}{qR\sin{\alpha}}\right)
47
48where $d$ = thickness of the layer (*thick_layer*),
49$2h$ = core thickness (*thick_core*), and $R$ = radius of the disc (*radius*).
50
51.. math::
52
53    S(q,\alpha) = 1 + \frac{1}{2}\sum_{k=1}^n(n-k)\cos{(kDq\cos{\alpha})}
54    \exp\left[ -k(q)^2(D\cos{\alpha}~\sigma_d)^2/2\right]
55
56where $n$ is the total number of the disc stacked (*n_stacking*),
57$D = 2(d+h)$ is the next neighbor center-to-center distance (d-spacing),
58and $\sigma_d$ = the Gaussian standard deviation of the d-spacing (*sigma_d*).
59Note that $D\cos(\alpha)$ is the component of $D$ parallel to $q$ and the last
60term in the equation above is effectively a Debye-Waller factor term.
61
62.. note::
63
64    1. Each assembly in the stack is layer/core/layer, so the spacing of the
65    cores is core plus two layers. The 2nd virial coefficient of the cylinder
66    is calculated based on the *radius* and *length*
67    = *n_stacking* * (*thick_core* + 2 * *thick_layer*)
68    values, and used as the effective radius for $S(Q)$ when $P(Q) * S(Q)$
69    is applied.
70
71    2. the resolution smearing calculation uses 76 Gaussian quadrature points
72    to properly smear the model since the function is HIGHLY oscillatory,
73    especially around the q-values that correspond to the repeat distance of
74    the layers.
75
762d scattering from oriented stacks is calculated in the same way as for
77cylinders, for further details of the calculation and angular dispersions
78see :ref:`orientation`.
79
80.. figure:: img/cylinder_angle_definition.png
81
82    Angles $\theta$ and $\phi$ orient the stack of discs relative
83    to the beam line coordinates, where the beam is along the $z$ axis.
84    Rotation $\theta$, initially in the $xz$ plane, is carried out first,
85    then rotation $\phi$ about the $z$ axis. Orientation distributions are
86    described as rotations about two perpendicular axes $\delta_1$ and
87    $\delta_2$ in the frame of the cylinder itself, which when
88    $\theta = \phi = 0$ are parallel to the $Y$ and $X$ axes.
89
90
91Our model is derived from the form factor calculations implemented in a
92c-library provided by the NIST Center for Neutron Research\ [#CIT_Kline]_
93
94References
95----------
96
97.. [#CIT1949] O Kratky and G Porod, *J. Colloid Science*, 4, (1949) 35
98.. [#CIT_Kline] S R Kline, *J Appl. Cryst.*, 39 (2006) 895
99.. [#] J S Higgins and H C Benoit, *Polymers and Neutron Scattering*,
100   Clarendon, Oxford, 1994
101.. [#] A Guinier and G Fournet, *Small-Angle Scattering of X-Rays*,
102   John Wiley and Sons, New York, 1955
103
104Authorship and Verification
105----------------------------
106
107* **Author:** NIST IGOR/DANSE **Date:** pre 2010
108* **Last Modified by:** Paul Butler and Paul Kienzle **Date:** November 26, 2016
109* **Last Reviewed by:** Paul Butler and Paul Kienzle **Date:** November 26, 2016
110"""
111
112import numpy as np
113from numpy import inf, sin, cos, pi
114
115name = "stacked_disks"
116title = "Form factor for a stacked set of non exfoliated core/shell disks"
117description = """\
118    One layer of disk consists of a core, a top layer, and a bottom layer.
119    radius =  the radius of the disk
120    thick_core = thickness of the core
121    thick_layer = thickness of a layer
122    sld_core = the SLD of the core
123    sld_layer = the SLD of the layers
124    n_stacking = the number of the disks
125    sigma_d =  Gaussian STD of d-spacing
126    sld_solvent = the SLD of the solvent
127    """
128category = "shape:cylinder"
129
130# pylint: disable=bad-whitespace, line-too-long
131#   ["name", "units", default, [lower, upper], "type","description"],
132parameters = [
133    ["thick_core",  "Ang",        10.0, [0, inf],    "volume",      "Thickness of the core disk"],
134    ["thick_layer", "Ang",        10.0, [0, inf],    "volume",      "Thickness of layer each side of core"],
135    ["radius",      "Ang",        15.0, [0, inf],    "volume",      "Radius of the stacked disk"],
136    ["n_stacking",  "",            1.0, [1, inf],    "volume",      "Number of stacked layer/core/layer disks"],
137    ["sigma_d",     "Ang",         0,   [0, inf],    "",            "Sigma of nearest neighbor spacing"],
138    ["sld_core",    "1e-6/Ang^2",  4,   [-inf, inf], "sld",         "Core scattering length density"],
139    ["sld_layer",   "1e-6/Ang^2",  0.0, [-inf, inf], "sld",         "Layer scattering length density"],
140    ["sld_solvent", "1e-6/Ang^2",  5.0, [-inf, inf], "sld",         "Solvent scattering length density"],
141    ["theta",       "degrees",     0,   [-360, 360], "orientation", "Orientation of the stacked disk axis w/respect incoming beam"],
142    ["phi",         "degrees",     0,   [-360, 360], "orientation", "Rotation about beam"],
143    ]
144# pylint: enable=bad-whitespace, line-too-long
145
146source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "stacked_disks.c"]
147
148def random():
149    radius = 10**np.random.uniform(1, 4.7)
150    total_stack = 10**np.random.uniform(1, 4.7)
151    n_stacking = int(10**np.random.uniform(0, np.log10(total_stack)-1) + 0.5)
152    d = total_stack/n_stacking
153    thick_core = np.random.uniform(0, d-2)  # at least 1 A for each layer
154    thick_layer = (d - thick_core)/2
155    # Let polydispersity peak around 15%; 95% < 0.4; max=100%
156    sigma_d = d * np.random.beta(1.5, 7)
157    pars = dict(
158        thick_core=thick_core,
159        thick_layer=thick_layer,
160        radius=radius,
161        n_stacking=n_stacking,
162        sigma_d=sigma_d,
163    )
164    return pars
165
166demo = dict(background=0.001,
167            scale=0.01,
168            thick_core=10.0,
169            thick_layer=10.0,
170            radius=15.0,
171            n_stacking=1,
172            sigma_d=0,
173            sld_core=4,
174            sld_layer=0.0,
175            sld_solvent=5.0,
176            theta=90,
177            phi=0)
178# After redefinition of spherical coordinates -
179# tests had in old coords theta=0, phi=0; new coords theta=90, phi=0
180q = 0.1
181# april 6 2017, rkh added a 2d unit test, assume correct!
182qx = q*cos(pi/6.0)
183qy = q*sin(pi/6.0)
184# Accuracy tests based on content in test/utest_extra_models.py.
185# Added 2 tests with n_stacked = 5 using SasView 3.1.2 - PDB;
186# for which alas q=0.001 values seem closer to n_stacked=1 not 5,
187# changed assuming my 4.1 code OK, RKH
188tests = [
189    [{'thick_core': 10.0,
190      'thick_layer': 15.0,
191      'radius': 3000.0,
192      'n_stacking': 1.0,
193      'sigma_d': 0.0,
194      'sld_core': 4.0,
195      'sld_layer': -0.4,
196      'sld_solvent': 5.0,
197      'theta': 90.0,
198      'phi': 0.0,
199      'scale': 0.01,
200      'background': 0.001,
201     }, 0.001, 5075.12],
202    [{'thick_core': 10.0,
203      'thick_layer': 15.0,
204      'radius': 3000.0,
205      'n_stacking': 5,
206      'sigma_d': 0.0,
207      'sld_core': 4.0,
208      'sld_layer': -0.4,
209      'sld_solvent': 5.0,
210      'theta': 90.0,
211      'phi': 0.0,
212      'scale': 0.01,
213      'background': 0.001,
214      # n_stacking=1 not 5 ? slight change in value here 11jan2017,
215      # check other cpu types
216      #}, 0.001, 5065.12793824],
217      #}, 0.001, 5075.11570493],
218     }, 0.001, 25325.635693],
219    [{'thick_core': 10.0,
220      'thick_layer': 15.0,
221      'radius': 100.0,
222      'n_stacking': 5,
223      'sigma_d': 0.0,
224      'sld_core': 4.0,
225      'sld_layer': -0.4,
226      'sld_solvent': 5.0,
227      'theta': 90.0,
228      'phi': 20.0,
229      'scale': 0.01,
230      'background': 0.001,
231     }, (qx, qy), 0.0491167089952],
232    [{'thick_core': 10.0,
233      'thick_layer': 15.0,
234      'radius': 3000.0,
235      'n_stacking': 5,
236      'sigma_d': 0.0,
237      'sld_core': 4.0,
238      'sld_layer': -0.4,
239      'sld_solvent': 5.0,
240      'theta': 90.0,
241      'phi': 0.0,
242      'scale': 0.01,
243      'background': 0.001,
244      # n_stacking=1 not 5 ?  slight change in value here 11jan2017,
245      # check other cpu types
246      #}, 0.164, 0.0127673597265],
247      #}, 0.164, 0.01480812968],
248     }, 0.164, 0.0598367986509],
249
250    [{'thick_core': 10.0,
251      'thick_layer': 15.0,
252      'radius': 3000.0,
253      'n_stacking': 1.0,
254      'sigma_d': 0.0,
255      'sld_core': 4.0,
256      'sld_layer': -0.4,
257      'sld_solvent': 5.0,
258      'theta': 90.0,
259      'phi': 0.0,
260      'scale': 0.01,
261      'background': 0.001,
262      # second test here was at q=90, changed it to q=5,
263      # note I(q) is then just value of flat background
264     }, [0.001, 5.0], [5075.12, 0.001]],
265
266    [{'thick_core': 10.0,
267      'thick_layer': 15.0,
268      'radius': 3000.0,
269      'n_stacking': 1.0,
270      'sigma_d': 0.0,
271      'sld_core': 4.0,
272      'sld_layer': -0.4,
273      'sld_solvent': 5.0,
274      'theta': 90.0,
275      'phi': 0.0,
276      'scale': 0.01,
277      'background': 0.001,
278     }, ([0.4, 0.5]), [0.00105074, 0.00121761]],
279    #[{'thick_core': 10.0,
280    #  'thick_layer': 15.0,
281    #  'radius': 3000.0,
282    #  'n_stacking': 1.0,
283    #  'sigma_d': 0.0,
284    #  'sld_core': 4.0,
285    #  'sld_layer': -0.4,
286    #  'sld_solvent': 5.0,
287    #  'theta': 90.0,
288    #  'phi': 20.0,
289    #  'scale': 0.01,
290    #  'background': 0.001,
291    # 2017-05-18 PAK temporarily suppress output of qx,qy test; j1 is
292    #     not accurate for large qr
293    # }, (qx, qy), 0.0341738733124],
294    # }, (qx, qy), None],
295
296    [{'thick_core': 10.0,
297      'thick_layer': 15.0,
298      'radius': 3000.0,
299      'n_stacking': 1.0,
300      'sigma_d': 0.0,
301      'sld_core': 4.0,
302      'sld_layer': -0.4,
303      'sld_solvent': 5.0,
304      'theta': 90.0,
305      'phi': 0.0,
306      'scale': 0.01,
307      'background': 0.001,
308     }, ([1.3, 1.57]), [0.0010039, 0.0010038]],
309    ]
310# 11Jan2017   RKH checking unit test again, note they are all 1D, no 2D
Note: See TracBrowser for help on using the repository browser.