source:sasmodels/sasmodels/models/core_shell_parallelepiped.py@2d81cfe

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

lint

• Property mode set to 100644
File size: 9.0 KB
Line
1r"""
2Definition
3----------
4
5Calculates the form factor for a rectangular solid with a core-shell structure.
6The thickness and the scattering length density of the shell or
7"rim" can be different on each (pair) of faces. However at this time the 1D
8calculation does **NOT** actually calculate a c face rim despite the presence
9of the parameter. Some other aspects of the 1D calculation may be wrong.
10
11.. note::
12   This model was originally ported from NIST IGOR macros. However, it is not
13   yet fully understood by the SasView developers and is currently under review.
14
15The form factor is normalized by the particle volume $V$ such that
16
17.. math::
18
19    I(q) = \text{scale}\frac{\langle f^2 \rangle}{V} + \text{background}
20
21where $\langle \ldots \rangle$ is an average over all possible orientations
22of the rectangular solid.
23
24
25The function calculated is the form factor of the rectangular solid below.
26The core of the solid is defined by the dimensions $A$, $B$, $C$ such that
27$A < B < C$.
28
29.. image:: img/core_shell_parallelepiped_geometry.jpg
30
31There are rectangular "slabs" of thickness $t_A$ that add to the $A$ dimension
32(on the $BC$ faces). There are similar slabs on the $AC$ $(=t_B)$ and $AB$
33$(=t_C)$ faces. The projection in the $AB$ plane is then
34
35.. image:: img/core_shell_parallelepiped_projection.jpg
36
37The volume of the solid is
38
39.. math::
40
41    V = ABC + 2t_ABC + 2t_BAC + 2t_CAB
42
43**meaning that there are "gaps" at the corners of the solid.**  Again note that
44$t_C = 0$ currently.
45
46The intensity calculated follows the :ref:parallelepiped model, with the
47core-shell intensity being calculated as the square of the sum of the
48amplitudes of the core and shell, in the same manner as a core-shell model.
49
50.. math::
51
52    F_{a}(Q,\alpha,\beta)=
53    \left[\frac{\sin(\tfrac{1}{2}Q(L_A+2t_A)\sin\alpha \sin\beta)
54               }{\tfrac{1}{2}Q(L_A+2t_A)\sin\alpha\sin\beta}
55    - \frac{\sin(\tfrac{1}{2}QL_A\sin\alpha \sin\beta)
56           }{\tfrac{1}{2}QL_A\sin\alpha \sin\beta} \right]
57    \left[\frac{\sin(\tfrac{1}{2}QL_B\sin\alpha \sin\beta)
58               }{\tfrac{1}{2}QL_B\sin\alpha \sin\beta} \right]
59    \left[\frac{\sin(\tfrac{1}{2}QL_C\sin\alpha \sin\beta)
60               }{\tfrac{1}{2}QL_C\sin\alpha \sin\beta} \right]
61
62.. note::
63
64    Why does t_B not appear in the above equation?
65    For the calculation of the form factor to be valid, the sides of the solid
66    MUST (perhaps not any more?) be chosen such that** $A < B < C$.
67    If this inequality is not satisfied, the model will not report an error,
68    but the calculation will not be correct and thus the result wrong.
69
70FITTING NOTES
71If the scale is set equal to the particle volume fraction, $\phi$, the returned
72value is the scattered intensity per unit volume, $I(q) = \phi P(q)$.
73However, **no interparticle interference effects are included in this
74calculation.**
75
76There are many parameters in this model. Hold as many fixed as possible with
77known values, or you will certainly end up at a solution that is unphysical.
78
79Constraints must be applied during fitting to ensure that the inequality
80$A < B < C$ is not violated. The calculation will not report an error,
81but the results will not be correct.
82
83The returned value is in units of |cm^-1|, on absolute scale.
84
85NB: The 2nd virial coefficient of the core_shell_parallelepiped is calculated
86based on the the averaged effective radius $(=\sqrt{(A+2t_A)(B+2t_B)/\pi})$
87and length $(C+2t_C)$ values, after appropriately
88sorting the three dimensions to give an oblate or prolate particle, to give an
89effective radius, for $S(Q)$ when $P(Q) * S(Q)$ is applied.
90
91For 2d data the orientation of the particle is required, described using
92angles $\theta$, $\phi$ and $\Psi$ as in the diagrams below, for further details
93of the calculation and angular dispersions see :ref:orientation .
94The angle $\Psi$ is the rotational angle around the *long_c* axis. For example,
95$\Psi = 0$ when the *short_b* axis is parallel to the *x*-axis of the detector.
96
97.. figure:: img/parallelepiped_angle_definition.png
98
99    Definition of the angles for oriented core-shell parallelepipeds.
100    Note that rotation $\theta$, initially in the $xz$ plane, is carried
101    out first, then rotation $\phi$ about the $z$ axis, finally rotation
102    $\Psi$ is now around the axis of the cylinder. The neutron or X-ray
103    beam is along the $z$ axis.
104
105.. figure:: img/parallelepiped_angle_projection.png
106
107    Examples of the angles for oriented core-shell parallelepipeds against the
108    detector plane.
109
110References
111----------
112
113.. [#] P Mittelbach and G Porod, *Acta Physica Austriaca*, 14 (1961) 185-211
114    Equations (1), (13-14). (in German)
115.. [#] D Singh (2009). *Small angle scattering studies of self assembly in
116   lipid mixtures*, John's Hopkins University Thesis (2009) 223-225. Available
117   from Proquest <http://search.proquest.com/docview/304915826?accountid
118   =26379>_
119
120Authorship and Verification
121----------------------------
122
123* **Author:** NIST IGOR/DANSE **Date:** pre 2010
124* **Converted to sasmodels by:** Miguel Gonzales **Date:** February 26, 2016
126* **Currently Under review by:** Paul Butler
127"""
128
129import numpy as np
130from numpy import pi, inf, sqrt, cos, sin
131
132name = "core_shell_parallelepiped"
133title = "Rectangular solid with a core-shell structure."
134description = """
135     P(q)=
136"""
137category = "shape:parallelepiped"
138
139#             ["name", "units", default, [lower, upper], "type","description"],
140parameters = [["sld_core", "1e-6/Ang^2", 1, [-inf, inf], "sld",
141               "Parallelepiped core scattering length density"],
142              ["sld_a", "1e-6/Ang^2", 2, [-inf, inf], "sld",
143               "Parallelepiped A rim scattering length density"],
144              ["sld_b", "1e-6/Ang^2", 4, [-inf, inf], "sld",
145               "Parallelepiped B rim scattering length density"],
146              ["sld_c", "1e-6/Ang^2", 2, [-inf, inf], "sld",
147               "Parallelepiped C rim scattering length density"],
148              ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "sld",
149               "Solvent scattering length density"],
150              ["length_a", "Ang", 35, [0, inf], "volume",
151               "Shorter side of the parallelepiped"],
152              ["length_b", "Ang", 75, [0, inf], "volume",
153               "Second side of the parallelepiped"],
154              ["length_c", "Ang", 400, [0, inf], "volume",
155               "Larger side of the parallelepiped"],
156              ["thick_rim_a", "Ang", 10, [0, inf], "volume",
157               "Thickness of A rim"],
158              ["thick_rim_b", "Ang", 10, [0, inf], "volume",
159               "Thickness of B rim"],
160              ["thick_rim_c", "Ang", 10, [0, inf], "volume",
161               "Thickness of C rim"],
162              ["theta", "degrees", 0, [-360, 360], "orientation",
163               "c axis to beam angle"],
164              ["phi", "degrees", 0, [-360, 360], "orientation",
166              ["psi", "degrees", 0, [-360, 360], "orientation",
168             ]
169
170source = ["lib/gauss76.c", "core_shell_parallelepiped.c"]
171
172
173def ER(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c):
174    """
176    """
177
178    # surface average radius (rough approximation)
179    surf_rad = sqrt((length_a + 2.0*thick_rim_a) * (length_b + 2.0*thick_rim_b) / pi)
180
181    height = length_c + 2.0*thick_rim_c
182
184    return 0.5 * (ddd) ** (1. / 3.)
185
186# VR defaults to 1.0
187
188def random():
189    outer = 10**np.random.uniform(1, 4.7, size=3)
190    thick = np.random.beta(0.5, 0.5, size=3)*(outer-2) + 1
191    length = outer - thick
192    pars = dict(
193        length_a=length[0],
194        length_b=length[1],
195        length_c=length[2],
196        thick_rim_a=thick[0],
197        thick_rim_b=thick[1],
198        thick_rim_c=thick[2],
199    )
200    return pars
201
202# parameters for demo
203demo = dict(scale=1, background=0.0,
204            sld_core=1, sld_a=2, sld_b=4, sld_c=2, sld_solvent=6,
205            length_a=35, length_b=75, length_c=400,
206            thick_rim_a=10, thick_rim_b=10, thick_rim_c=10,
207            theta=0, phi=0, psi=0,
208            length_a_pd=0.1, length_a_pd_n=1,
209            length_b_pd=0.1, length_b_pd_n=1,
210            length_c_pd=0.1, length_c_pd_n=1,
211            thick_rim_a_pd=0.1, thick_rim_a_pd_n=1,
212            thick_rim_b_pd=0.1, thick_rim_b_pd_n=1,
213            thick_rim_c_pd=0.1, thick_rim_c_pd_n=1,
214            theta_pd=10, theta_pd_n=1,
215            phi_pd=10, phi_pd_n=1,
216            psi_pd=10, psi_pd_n=1)
217
218# rkh 7/4/17 add random unit test for 2d, note make all params different, 2d values not tested against other codes or models
219if 0:  # pak: model rewrite; need to update tests
220    qx, qy = 0.2 * cos(pi/6.), 0.2 * sin(pi/6.)
221    tests = [[{}, 0.2, 0.533149288477],
222             [{}, [0.2], [0.533149288477]],
223             [{'theta':10.0, 'phi':20.0}, (qx, qy), 0.0853299803222],
224             [{'theta':10.0, 'phi':20.0}, [(qx, qy)], [0.0853299803222]],
225            ]
226    del qx, qy  # not necessary to delete, but cleaner
Note: See TracBrowser for help on using the repository browser.