source: sasmodels/sasmodels/models/flexible_cylinder.py @ c1e44e5

Last change on this file since c1e44e5 was c1e44e5, checked in by Paul Kienzle <pkienzle@…>, 5 years ago

Add local link to source files. Refs #1263.

  • Property mode set to 100644
File size: 6.6 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
66References
67----------
68
69.. [#] J S Pedersen and P Schurtenberger. *Scattering functions of semiflexible polymers with and without excluded volume effects.* Macromolecules, 29 (1996) 7602-7612
70
71Correction of the formula can be found in
72
73.. [#] 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
74
75Authorship and Verification
76----------------------------
77
78* **Author:**
79* **Last Modified by:**
80* **Last Reviewed by:** Steve King **Date:** March 26, 2019
81"""
82
83import numpy as np
84from numpy import inf
85
86name = "flexible_cylinder"
87title = "Flexible cylinder where the form factor is normalized by the volume " \
88        "of the cylinder."
89description = """Note : scale and contrast = (sld - sld_solvent) are both
90                multiplicative factors in the model and are perfectly
91                correlated. One or both of these parameters must be held fixed
92                during model fitting.
93              """
94
95category = "shape:cylinder"
96single = False  # double precision only!
97
98# pylint: disable=bad-whitespace, line-too-long
99#             ["name", "units", default, [lower, upper], "type", "description"],
100parameters = [
101    ["length",      "Ang",       1000.0, [0, inf],    "volume", "Length of the flexible cylinder"],
102    ["kuhn_length", "Ang",        100.0, [0, inf],    "volume", "Kuhn length of the flexible cylinder"],
103    ["radius",      "Ang",         20.0, [0, inf],    "volume", "Radius of the flexible cylinder"],
104    ["sld",         "1e-6/Ang^2",   1.0, [-inf, inf], "sld",    "Cylinder scattering length density"],
105    ["sld_solvent", "1e-6/Ang^2",   6.3, [-inf, inf], "sld",    "Solvent scattering length density"],
106    ]
107# pylint: enable=bad-whitespace, line-too-long
108source = ["lib/polevl.c", "lib/sas_J1.c", "lib/wrc_cyl.c", "flexible_cylinder.c"]
109
110def random():
111    """Return a random parameter set for the model."""
112    length = 10**np.random.uniform(2, 6)
113    radius = 10**np.random.uniform(1, 3)
114    kuhn_length = 10**np.random.uniform(-2, 0)*length
115    pars = dict(
116        length=length,
117        radius=radius,
118        kuhn_length=kuhn_length,
119    )
120    return pars
121
122tests = [
123    # Accuracy tests based on content in test/utest_other_models.py
124    [{'length':     1000.0,  # test T1
125      'kuhn_length': 100.0,
126      'radius':       20.0,
127      'sld':           1.0,
128      'sld_solvent':   6.3,
129      'background':    0.0001,
130     }, 0.001, 3509.2187],
131
132    # Additional tests with larger range of parameters
133    [{'length':    1000.0,  # test T2
134      'kuhn_length': 100.0,
135      'radius':       20.0,
136      'sld':           1.0,
137      'sld_solvent':   6.3,
138      'background':    0.0001,
139     }, 1.0, 0.000595345],
140    [{'length':        10.0,  # test T3
141      'kuhn_length': 800.0,
142      'radius':        2.0,
143      'sld':           6.0,
144      'sld_solvent':  12.3,
145      'background':    0.001,
146     }, 0.1, 1.55228],
147    [{'length':        100.0,  # test T4
148      'kuhn_length': 800.0,
149      'radius':       50.0,
150      'sld':           0.1,
151      'sld_solvent':   5.1,
152      'background':    0.0,
153     }, 1.0, 0.000938456]
154    ]
155
156# There are a few branches in the code that ought to have test values:
157#
158# For length > 4 * kuhn_length
159#        if length > 10 * kuhn_length then C is scaled by 3.06 (L/b)^(-0.44)
160#        q*kuhn_length <= 3.1  => Sexv_new
161#           dS/dQ < 0 has different behaviour from dS/dQ >= 0
162#  T2    q*kuhn_length > 3.1   => a_long
163#
164# For length <= 4 * kuhn_length
165#        q*kuhn_length <= max(1.9/Rg_short, 3.0)  => Sdebye((q*Rg)^2)
166#           q*Rg < 0.5 uses Pade approx, q*Rg > 1.0 uses math lib
167#  T3,T4 q*kuhn_length > max(1.9/Rg_short, 3.0)   => a_short
168#
169# Note that the transitions between branches may be abrupt.  You can see a
170# several percent change around length=10*kuhn_length and length=4*kuhn_length
171# using the following:
172#
173#    sascomp flexible_cylinder -calc=double -sets=10 length=10*kuhn_length,10.000001*kuhn_length
174#    sascomp flexible_cylinder -calc=double -sets=10 length=4*kuhn_length,4.000001*kuhn_length
175#
176# The transition between low q and high q around q*kuhn_length = 3 seems
177# to be good to 4 digits or better.  This was tested by computing the value
178# on each branches near the transition point and reporting the relative error
179# for kuhn lengths of 10, 100 and 1000 and a variety of length:kuhn_length
180# ratios.
Note: See TracBrowser for help on using the repository browser.