source: sasmodels/sasmodels/models/flexible_cylinder.py @ 830cf6b

ticket-1257-vesicle-productticket_1156ticket_822_more_unit_tests
Last change on this file since 830cf6b was 830cf6b, checked in by Paul Kienzle <pkienzle@…>, 5 years ago

Merge branch 'master' into ticket_822_v5_unit_tests

  • Property mode set to 100644
File size: 7.9 KB
Line 
1r"""
2This model provides the form factor, $P(q)$, for a flexible cylinder
3where the form factor is normalized by the volume of the cylinder.
4**Inter-cylinder interactions are NOT provided for.**
5
6.. math::
7
8    P(q) = \text{scale} \left<F^2\right>/V + \text{background}
9
10where the averaging $\left<\ldots\right>$ is applied only for the 1D
11calculation
12
13The 2D scattering intensity is the same as 1D, regardless of the orientation of
14the q vector which is defined as
15
16.. math::
17
18    q = \sqrt{q_x^2 + q_y^2}
19
20Definitions
21-----------
22
23.. figure:: img/flexible_cylinder_geometry.jpg
24
25
26The chain of contour length, $L$, (the total length) can be described as a
27chain of some number of locally stiff segments of length $l_p$, the persistence
28length (the length along the cylinder over which the flexible cylinder can be
29considered a rigid rod).
30The Kuhn length $(b = 2*l_p)$ is also used to describe the stiffness of a chain.
31
32In the parameters, the sld and sld\_solvent represent the SLD of the cylinder
33and solvent respectively.
34
35Our model uses the form factor calculations in reference [1] as implemented in a
36c-library provided by the NIST Center for Neutron Research (Kline, 2006). This states:
37
38    'Method 3 With Excluded Volume' is used.
39    The model is a parametrization of simulations of a discrete representation
40    of the worm-like chain model of Kratky and Porod applied in the
41    pseudocontinuous limit.
42    See equations (13,26-27) in the original reference for the details.
43
44.. note::
45
46    There are several typos in the original reference that have been corrected
47    by WRC [2]. Details of the corrections are in the reference below. Most notably
48
49    - Equation (13): the term $(1 - w(QR))$ should swap position with $w(QR)$
50
51    - Equations (23) and (24) are incorrect; WRC has entered these into
52      Mathematica and solved analytically. The results were then converted to
53      code.
54
55    - Equation (27) should be $q0 = max(a3/(Rg^2)^{1/2},3)$ instead of
56      $max(a3*b(Rg^2)^{1/2},3)$
57
58    - The scattering function is negative for a range of parameter values and
59      q-values that are experimentally accessible. A correction function has been
60      added to give the proper behavior.
61
62
63**This is a model with complex behaviour depending on the ratio of** $L/b$ **and the
64reader is strongly encouraged to read reference [1] before use.**
65
66.. note::
67
68    There are several typos in the original reference that have been corrected
69    by WRC [2]. Details of the corrections are in the reference below. Most notably
70
71    - Equation (13): the term $(1 - w(QR))$ should swap position with $w(QR)$
72
73    - Equations (23) and (24) are incorrect; WRC has entered these into
74      Mathematica and solved analytically. The results were then converted to
75      code.
76
77    - Equation (27) should be $q0 = max(a3/(Rg^2)^{1/2},3)$ instead of
78      $max(a3*b(Rg^2)^{1/2},3)$
79
80    - The scattering function is negative for a range of parameter values and
81      q-values that are experimentally accessible. A correction function has been
82      added to give the proper behavior.
83
84
85**This is a model with complex behaviour depending on the ratio of** $L/b$ **and the
86reader is strongly encouraged to read reference [1] before use.**
87
88References
89----------
90
91.. [#] J S Pedersen and P Schurtenberger. *Scattering functions of semiflexible polymers with and without excluded volume effects.* Macromolecules, 29 (1996) 7602-7612
92
93Correction of the formula can be found in
94
95.. [#] W R Chen, P D Butler and L J Magid, *Incorporating Intermicellar Interactions in the Fitting of SANS Data from Cationic Wormlike Micelles.* Langmuir, 22(15) 2006 6539-6548
96
97Source
98------
99
100`flexible_cylinder.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/flexible_cylinder.py>`_
101
102`flexible_cylinder.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/flexible_cylinder.c>`_
103
104`wrc_cyl.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/lib/wrc_cyl.c>`_
105
106Authorship and Verification
107----------------------------
108
109* **Author:**
110* **Last Modified by:**
111* **Last Reviewed by:** Steve King **Date:** March 26, 2019
112* **Source added by :** Steve King **Date:** March 25, 2019
113"""
114
115import numpy as np
116from numpy import inf
117
118name = "flexible_cylinder"
119title = "Flexible cylinder where the form factor is normalized by the volume " \
120        "of the cylinder."
121description = """Note : scale and contrast = (sld - sld_solvent) are both
122                multiplicative factors in the model and are perfectly
123                correlated. One or both of these parameters must be held fixed
124                during model fitting.
125              """
126
127category = "shape:cylinder"
128single = False  # double precision only!
129
130# pylint: disable=bad-whitespace, line-too-long
131#             ["name", "units", default, [lower, upper], "type", "description"],
132parameters = [
133    ["length",      "Ang",       1000.0, [0, inf],    "volume", "Length of the flexible cylinder"],
134    ["kuhn_length", "Ang",        100.0, [0, inf],    "volume", "Kuhn length of the flexible cylinder"],
135    ["radius",      "Ang",         20.0, [0, inf],    "volume", "Radius of the flexible cylinder"],
136    ["sld",         "1e-6/Ang^2",   1.0, [-inf, inf], "sld",    "Cylinder scattering length density"],
137    ["sld_solvent", "1e-6/Ang^2",   6.3, [-inf, inf], "sld",    "Solvent scattering length density"],
138    ]
139# pylint: enable=bad-whitespace, line-too-long
140source = ["lib/polevl.c", "lib/sas_J1.c", "lib/wrc_cyl.c", "flexible_cylinder.c"]
141
142def random():
143    """Return a random parameter set for the model."""
144    length = 10**np.random.uniform(2, 6)
145    radius = 10**np.random.uniform(1, 3)
146    kuhn_length = 10**np.random.uniform(-2, 0)*length
147    pars = dict(
148        length=length,
149        radius=radius,
150        kuhn_length=kuhn_length,
151    )
152    return pars
153
154tests = [
155    # Accuracy tests based on content in test/utest_other_models.py
156    [{'length':     1000.0,  # test T1
157      'kuhn_length': 100.0,
158      'radius':       20.0,
159      'sld':           1.0,
160      'sld_solvent':   6.3,
161      'background':    0.0001,
162     }, 0.001, 3509.2187],
163
164    # Additional tests with larger range of parameters
165    [{'length':    1000.0,  # test T2
166      'kuhn_length': 100.0,
167      'radius':       20.0,
168      'sld':           1.0,
169      'sld_solvent':   6.3,
170      'background':    0.0001,
171     }, 1.0, 0.000595345],
172    [{'length':        10.0,  # test T3
173      'kuhn_length': 800.0,
174      'radius':        2.0,
175      'sld':           6.0,
176      'sld_solvent':  12.3,
177      'background':    0.001,
178     }, 0.1, 1.55228],
179    [{'length':        100.0,  # test T4
180      'kuhn_length': 800.0,
181      'radius':       50.0,
182      'sld':           0.1,
183      'sld_solvent':   5.1,
184      'background':    0.0,
185     }, 1.0, 0.000938456]
186    ]
187
188# There are a few branches in the code that ought to have test values:
189#
190# For length > 4 * kuhn_length
191#        if length > 10 * kuhn_length then C is scaled by 3.06 (L/b)^(-0.44)
192#        q*kuhn_length <= 3.1  => Sexv_new
193#           dS/dQ < 0 has different behaviour from dS/dQ >= 0
194#  T2    q*kuhn_length > 3.1   => a_long
195#
196# For length <= 4 * kuhn_length
197#        q*kuhn_length <= max(1.9/Rg_short, 3.0)  => Sdebye((q*Rg)^2)
198#           q*Rg < 0.5 uses Pade approx, q*Rg > 1.0 uses math lib
199#  T3,T4 q*kuhn_length > max(1.9/Rg_short, 3.0)   => a_short
200#
201# Note that the transitions between branches may be abrupt.  You can see a
202# several percent change around length=10*kuhn_length and length=4*kuhn_length
203# using the following:
204#
205#    sascomp flexible_cylinder -calc=double -sets=10 length=10*kuhn_length,10.000001*kuhn_length
206#    sascomp flexible_cylinder -calc=double -sets=10 length=4*kuhn_length,4.000001*kuhn_length
207#
208# The transition between low q and high q around q*kuhn_length = 3 seems
209# to be good to 4 digits or better.  This was tested by computing the value
210# on each branches near the transition point and reporting the relative error
211# for kuhn lengths of 10, 100 and 1000 and a variety of length:kuhn_length
212# ratios.
Note: See TracBrowser for help on using the repository browser.