# source:sasmodels/sasmodels/models/ellipsoid.py@0d6e865

core_shell_microgelscostrafo411magnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 0d6e865 was 0d6e865, checked in by dirk, 6 years ago

Rewriting selected models in spherical coordinates. This fixes ticket #492 for these models.

• Property mode set to 100644
File size: 5.7 KB
Line
1# ellipsoid model
2# Note: model title and parameter table are inserted automatically
3r"""
4The form factor is normalized by the particle volume
5
6Definition
7----------
8
9The output of the 2D scattering intensity function for oriented ellipsoids
10is given by (Feigin, 1987)
11
12.. math::
13
14    P(q,\alpha) = \frac{\text{scale}}{V} F^2(q,\alpha) + \text{background}
15
16where
17
18.. math::
19
20    F(q,\alpha) = \frac{3 \Delta \rho V (\sin[qr(R_p,R_e,\alpha)]
21                - \cos[qr(R_p,R_e,\alpha)])}
22                {[qr(R_p,R_e,\alpha)]^3}
23
24and
25
26.. math::
27
28    r(R_p,R_e,\alpha) = \left[ R_e^2 \sin^2 \alpha
29        + R_p^2 \cos^2 \alpha \right]^{1/2}
30
31
32$\alpha$ is the angle between the axis of the ellipsoid and $\vec q$,
33$V = (4/3)\pi R_pR_e^2$ is the volume of the ellipsoid , $R_p$ is the polar radius along the
34rotational axis of the ellipsoid, $R_e$ is the equatorial radius perpendicular
35to the rotational axis of the ellipsoid and $\Delta \rho$ (contrast) is the
36scattering length density difference between the scatterer and the solvent.
37
38For randomly oriented particles:
39
40.. math::
41
42   F^2(q)=\int_{0}^{\pi/2}{F^2(q,\alpha)\sin(\alpha)d\alpha}
43
44
45To provide easy access to the orientation of the ellipsoid, we define
46the rotation axis of the ellipsoid using two angles $\theta$ and $\phi$.
47These angles are defined in the
48:ref:cylinder orientation figure <cylinder-angle-definition>.
49For the ellipsoid, $\theta$ is the angle between the rotational axis
50and the $z$ -axis.
51
52NB: The 2nd virial coefficient of the solid ellipsoid is calculated based
53on the $R_p$ and $R_e$ values, and used as the effective radius for
54$S(q)$ when $P(q) \cdot S(q)$ is applied.
55
56
57The $\theta$ and $\phi$ parameters are not used for the 1D output.
58
59
60
61Validation
62----------
63
64Validation of the code was done by comparing the output of the 1D model
65to the output of the software provided by the NIST (Kline, 2006).
66
67The implementation of the intensity for fully oriented ellipsoids was
68validated by averaging the 2D output using a uniform distribution
69$p(\theta,\phi) = 1.0$ and comparing with the output of the 1D calculation.
70
71
72.. _ellipsoid-comparison-2d:
73
74.. figure:: img/ellipsoid_comparison_2d.jpg
75
76    Comparison of the intensity for uniformly distributed ellipsoids
77    calculated from our 2D model and the intensity from the NIST SANS
78    analysis software. The parameters used were: *scale* = 1.0,
80    *contrast* = 3e-6 |Ang^-2|, and *background* = 0.0 |cm^-1|.
81
82The discrepancy above $q$ = 0.3 |cm^-1| is due to the way the form factors
83are calculated in the c-library provided by NIST. A numerical integration
84has to be performed to obtain $P(q)$ for randomly oriented particles.
85The NIST software performs that integration with a 76-point Gaussian
86quadrature rule, which will become imprecise at high $q$ where the amplitude
87varies quickly as a function of $q$. The SasView result shown has been
88obtained by summing over 501 equidistant points. Our result was found
89to be stable over the range of $q$ shown for a number of points higher
90than 500.
91
92References
93----------
94
95L A Feigin and D I Svergun.
96*Structure Analysis by Small-Angle X-Ray and Neutron Scattering*,
97Plenum Press, New York, 1987.
98"""
99
100from numpy import inf
101
102name = "ellipsoid"
103title = "Ellipsoid of revolution with uniform scattering length density."
104
105description = """\
106P(q.alpha)= scale*f(q)^2 + background, where f(q)= 3*(sld
107        - sld_solvent)*V*[sin(q*r(Rp,Re,alpha))
108        -q*r*cos(qr(Rp,Re,alpha))]
109        /[qr(Rp,Re,alpha)]^3"
110
111     r(Rp,Re,alpha)= [Re^(2)*(sin(alpha))^2
112        + Rp^(2)*(cos(alpha))^2]^(1/2)
113
114        sld: SLD of the ellipsoid
115        sld_solvent: SLD of the solvent
116        V: volume of the ellipsoid
117        Rp: polar radius of the ellipsoid
118        Re: equatorial radius of the ellipsoid
119"""
120category = "shape:ellipsoid"
121
122#             ["name", "units", default, [lower, upper], "type","description"],
123parameters = [["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld",
124               "Ellipsoid scattering length density"],
125              ["sld_solvent", "1e-6/Ang^2", 1, [-inf, inf], "sld",
126               "Solvent scattering length density"],
127              ["radius_polar", "Ang", 20, [0, inf], "volume",
129              ["radius_equatorial", "Ang", 400, [0, inf], "volume",
131              ["theta", "degrees", 60, [-inf, inf], "orientation",
132               "In plane angle"],
133              ["phi", "degrees", 60, [-inf, inf], "orientation",
134               "Out of plane angle"],
135             ]
136
137source = ["lib/sph_j1c.c", "lib/gauss76.c", "ellipsoid.c"]
138
140    import numpy as np
141
148    ee[idx] = 2 * radius_polar[idx]
150    bd = 1.0 - ee[valid]
151    e1 = np.sqrt(ee[valid])
152    b1 = 1.0 + np.arcsin(e1) / (e1 * np.sqrt(bd))
153    bL = (1.0 + e1) / (1.0 - e1)
154    b2 = 1.0 + bd / 2 / e1 * np.log(bL)
155    delta = 0.75 * b1 * b2
156