source: sasmodels/sasmodels/models/polymer_excl_volume.py @ 0507e09

core_shell_microgelsmagnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 0507e09 was 0507e09, checked in by smk78, 7 months ago

Added link to source code to each model. Closes #883

  • Property mode set to 100644
File size: 6.0 KB
Line 
1r"""
2This model describes the scattering from polymer chains subject to excluded
3volume effects and has been used as a template for describing mass fractals.
4
5Definition
6----------
7
8The form factor was originally presented in the following integral form
9(Benoit, 1957)
10
11.. math::
12
13    P(Q)=2\int_0^{1}dx(1-x)exp\left[-\frac{Q^2a^2}{6}n^{2v}x^{2v}\right]
14
15where $\nu$ is the excluded volume parameter
16(which is related to the Porod exponent $m$ as $\nu=1/m$ ),
17$a$ is the statistical segment length of the polymer chain,
18and $n$ is the degree of polymerization.
19
20This integral was put into an almost analytical form as follows
21(Hammouda, 1993)
22
23.. math::
24
25    P(Q)=\frac{1}{\nu U^{1/2\nu}}
26    \left\{
27        \gamma\left(\frac{1}{2\nu},U\right) -
28        \frac{1}{U^{1/2\nu}}\gamma\left(\frac{1}{\nu},U\right)
29    \right\}
30
31and later recast as (for example, Hore, 2013; Hammouda & Kim, 2017)
32
33.. math::
34
35    P(Q)=\frac{1}{\nu U^{1/2\nu}}\gamma\left(\frac{1}{2\nu},U\right) -
36    \frac{1}{\nu U^{1/\nu}}\gamma\left(\frac{1}{\nu},U\right)
37
38where $\gamma(x,U)$ is the incomplete gamma function
39
40.. math::
41
42    \gamma(x,U)=\int_0^{U}dt\ \exp(-t)t^{x-1}
43
44and the variable $U$ is given in terms of the scattering vector $Q$ as
45
46.. math::
47
48    U=\frac{Q^2a^2n^{2\nu}}{6} = \frac{Q^2R_{g}^2(2\nu+1)(2\nu+2)}{6}
49
50The two analytic forms are equivalent. In the 1993 paper
51
52.. math::
53
54    \frac{1}{\nu U^{1/2\nu}}
55
56has been factored out.
57
58**SasView implements the 1993 expression**.
59
60The square of the radius-of-gyration is defined as
61
62.. math::
63
64    R_{g}^2 = \frac{a^2n^{2\nu}}{(2\nu+1)(2\nu+2)}
65
66.. note::
67    This model applies only in the mass fractal range (ie, $5/3<=m<=3$ )
68    and **does not apply** to surface fractals ( $3<m<=4$ ).
69    It also does not reproduce the rigid rod limit (m=1) because it assumes chain
70    flexibility from the outset. It may cover a portion of the semi-flexible chain
71    range ( $1<m<5/3$ ).
72
73A low-Q expansion yields the Guinier form and a high-Q expansion yields the
74Porod form which is given by
75
76.. math::
77
78    P(Q\rightarrow \infty) = \frac{1}{\nu U^{1/2\nu}}\Gamma\left(
79    \frac{1}{2\nu}\right) - \frac{1}{\nu U^{1/\nu}}\Gamma\left(
80    \frac{1}{\nu}\right)
81
82Here $\Gamma(x) = \gamma(x,\infty)$ is the gamma function.
83
84The asymptotic limit is dominated by the first term
85
86.. math::
87
88    P(Q\rightarrow \infty) \sim \frac{1}{\nu U^{1/2\nu}}\Gamma\left(\frac{1}{2\nu}\right) =
89    \frac{m}{\left(QR_{g}\right)^m}\left[\frac{6}{(2\nu +1)(2\nu +2)} \right]^{m/2}
90    \Gamma (m/2)
91
92The special case when $\nu=0.5$ (or $m=1/\nu=2$ ) corresponds to Gaussian chains for
93which the form factor is given by the familiar Debye function.
94
95.. math::
96
97    P(Q) = \frac{2}{Q^4R_{g}^4} \left[\exp(-Q^2R_{g}^2) - 1 + Q^2R_{g}^2 \right]
98
99For 2D data: The 2D scattering intensity is calculated in the same way as 1D,
100where the $q$ vector is defined as
101
102.. math::
103
104    q = \sqrt{q_x^2 + q_y^2}
105
106
107References
108----------
109
110.. [#] H Benoit, *Comptes Rendus*, 245 (1957) 2244-2247
111.. [#] B Hammouda, *SANS from Homogeneous Polymer Mixtures - A Unified Overview, Advances in Polym. Sci.* 106 (1993) 87-133
112.. [#] M Hore et al, *Co-Nonsolvency of Poly(n-isopropylacrylamide) in Deuterated Water/Ethanol Mixtures* 46 (2013) 7894-7901
113.. [#] B Hammouda & M-H Kim, *The empirical core-chain model* 247 (2017) 434-440
114
115Source
116------
117
118`polymer_excl_volume.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/polymer_excl_volume.py>`_
119
120`polymer_excl_volume.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/polymer_excl_volume.c>`_
121
122Authorship and Verification
123----------------------------
124
125* **Author:**
126* **Last Modified by:**
127* **Last Reviewed by:**
128* **Source added by :** Steve King **Date:** March 25, 2019
129"""
130
131import numpy as np
132from numpy import inf, power, errstate
133from scipy.special import gammainc, gamma
134
135name = "polymer_excl_volume"
136title = "Polymer Excluded Volume model"
137description = """Compute the scattering intensity from polymers with excluded
138                volume effects.
139                rg:         radius of gyration
140                porod_exp:  Porod exponent
141              """
142category = "shape-independent"
143
144# pylint: disable=bad-whitespace, line-too-long
145#   ["name", "units", default, [lower, upper], "type", "description"],
146parameters = [
147    ["rg",        "Ang", 60.0, [0, inf], "", "Radius of Gyration"],
148    ["porod_exp", "",     3.0, [0, inf], "", "Porod exponent"],
149]
150# pylint: enable=bad-whitespace, line-too-long
151
152
153def Iq(q, rg=60.0, porod_exp=3.0):
154    """
155    :param q:         Input q-value (float or [float, float])
156    :param rg:        Radius of gyration
157    :param porod_exp: Porod exponent
158    :return:          Calculated intensity
159    """
160    usub = (q*rg)**2 * (2.0/porod_exp + 1.0) * (2.0/porod_exp + 2.0)/6.0
161    with errstate(divide='ignore', invalid='ignore'):
162        upow = power(usub, -0.5*porod_exp)
163        # Note: scipy gammainc is "regularized", being gamma(s,x)/Gamma(s),
164        # so need to scale by Gamma(s) to recover gamma(s, x).
165        result = (porod_exp*upow *
166                  (gamma(0.5*porod_exp)*gammainc(0.5*porod_exp, usub) -
167                   upow*gamma(porod_exp)*gammainc(porod_exp, usub)))
168    result[q <= 0] = 1.0
169
170    return result
171
172Iq.vectorized = True  # Iq accepts an array of q values
173
174def random():
175    """Return a random parameter set for the model."""
176    rg = 10**np.random.uniform(0, 4)
177    porod_exp = np.random.uniform(1e-3, 6)
178    scale = 10**np.random.uniform(1, 5)
179    pars = dict(
180        #background=0,
181        scale=scale,
182        rg=rg,
183        porod_exp=porod_exp,
184    )
185    return pars
186
187tests = [
188    # Accuracy tests based on content in test/polyexclvol_default_igor.txt
189    [{'rg': 60, 'porod_exp': 3.0}, 0.001, 0.999801],
190    [{'rg': 60, 'porod_exp': 3.0}, 0.105363, 0.0172751],
191    [{'rg': 60, 'porod_exp': 3.0, 'background': 0.0}, 0.665075, 6.56261e-05],
192
193    # Additional tests with larger range of parameters
194    [{'rg': 10, 'porod_exp': 4.0}, 0.1, 0.724436675809],
195    [{'rg': 2.2, 'porod_exp': 22.0, 'background': 100.0}, 5.0, 100.0],
196    [{'rg': 1.1, 'porod_exp': 1, 'background': 10.0, 'scale': 1.25},
197     20000., 10.0000712097]
198    ]
Note: See TracBrowser for help on using the repository browser.