source: sasmodels/sasmodels/models/squarewell.py @ 348557a

core_shell_microgelscostrafo411magnetic_modelrelease_v0.94release_v0.95ticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 348557a was 348557a, checked in by richardh, 8 years ago

rkh new hayter_msa, new squarewell, unit test to lamellar

  • Property mode set to 100644
File size: 4.4 KB
Line 
1# Note: model title and parameter table are inserted automatically
2r"""
3This calculates the interparticle structure factor for a square well fluid spherical particles. The mean spherical
4approximation (MSA) closure was used for this calculation, and is not the most appropriate closure for an attractive
5interparticle potential. This solution has been compared to Monte Carlo simulations for a square well fluid, showing
6this calculation to be limited in applicability to well depths |epsilon| < 1.5 kT and volume fractions |phi| < 0.08.
7
8Positive well depths correspond to an attractive potential well. Negative well depths correspond to a potential
9"shoulder", which may or may not be physically reasonable. The stickyhardsphere model may be a better choice in
10some circumstances. Computed values may behave badly at extremely small Q.
11
12The well width (*l*\ ) is defined as multiples of the particle diameter (2\*\ *R*\ )
13
14The interaction potential is:
15
16.. image:: img/image225.PNG
17
18.. math::
19
20    U(r) = \begin{cases}
21    \infty & r < 2R \\
22    -\epsilon , 2R \leq < 2R\lambda
23    0 & r \geq 2R
24    \end{cases}
25
26where $r$ is the distance from the center of the sphere of a radius $R$.
27
28
29For 2D data: The 2D scattering intensity is calculated in the same way as 1D, where the *q* vector is defined as
30
31.. math::
32
33    q = \sqrt{q_x^2 + q_y^2}
34
35.. figure:: img/squarewell_226.jpg
36
37    1D plot using the default values (in linear scale).*
38
39REFERENCE
40
41R V Sharma, K C Sharma, *Physica*, 89A (1977) 213
42
43"""
44from numpy import inf
45
46name = "squarewell"
47title = "Square well structure factor, with MSA closure"
48description = """\
49    [Square well structure factor, with MSA closure]
50        Interparticle structure factor S(Q)for a hard sphere fluid with
51        a narrow attractive well. Fits are prone to deliver non-physical
52        parameters, use with care and read the references in the full manual.
53        In sasview the effective radius will be calculated from the
54        parameters used in P(Q).
55"""
56category = "structure-factor"
57
58#single = False
59#             ["name", "units", default, [lower, upper], "type","description"],
60parameters = [
61    #   [ "name", "units", default, [lower, upper], "type",
62    #     "description" ],
63    ["effect_radius", "Ang", 50.0, [0, inf], "volume",
64     "effective radius of hard sphere"],
65    ["volfraction", "", 0.04, [0, 0.08], "",
66     "volume fraction of spheres"],
67    ["welldepth", "kT", 1.5, [0.0, 1.5], "",
68     "depth of well, epsilon"],
69    ["wellwidth", "diameters", 1.2, [0, inf], "",
70     "width of well in diameters (=2R) units"],
71    ]
72
73# No volume normalization despite having a volume parameter
74# This should perhaps be volume normalized?
75form_volume = """
76    return 1.0;
77    """
78
79Iq = """
80// single precision is very poor at extreme small Q, would need a Taylor series
81        double req,phis,edibkb,lambda,struc;
82        double sigma,eta,eta2,eta3,eta4,etam1,etam14,alpha,beta,gamm;
83        double x,sk,sk2,sk3,sk4,t1,t2,t3,t4,ck;
84        double S,C,SL,CL;
85        x= q;
86       
87        req = effect_radius;
88        phis = volfraction;
89        edibkb = welldepth;
90        lambda = wellwidth;
91       
92        sigma = req*2.;
93        eta = phis;
94        eta2 = eta*eta;
95        eta3 = eta*eta2;
96        eta4 = eta*eta3;
97        etam1 = 1. - eta;
98        etam14 = etam1*etam1*etam1*etam1;
99        // temp borrow sk for an intermediate calc
100        sk = 1.0 +2.0*eta;
101        sk *= sk;
102        alpha = (  sk + eta3*( eta-4.0 )  )/etam14;
103        beta = -(eta/3.0) * ( 18. + 20.*eta - 12.*eta2 + eta4 )/etam14;
104        gamm = 0.5*eta*( sk + eta3*(eta-4.) )/etam14;
105       
106        //  calculate the structure factor
107       
108        sk = x*sigma;
109        sk2 = sk*sk;
110        sk3 = sk*sk2;
111        sk4 = sk3*sk;
112        SINCOS(sk,S,C);
113        SINCOS(lambda*sk,SL,CL);
114        t1 = alpha * sk3 * ( S - sk * C );
115        t2 = beta * sk2 * 2.0*( sk*S - (0.5*sk2 - 1.)*C - 1.0 );
116        t3 = gamm*( ( 4.0*sk3 - 24.*sk ) * S - ( sk4 - 12.0*sk2 + 24.0 )*C + 24.0 );
117        t4 = -edibkb*sk3*(SL +sk*(C - lambda*CL) - S );
118        ck =  -24.0*eta*( t1 + t2 + t3 + t4 )/sk3/sk3;
119        struc  = 1./(1.-ck);
120       
121        return(struc);
122"""
123
124Iqxy = """
125    return Iq(sqrt(qx*qx+qy*qy), IQ_PARAMETERS);
126    """
127
128# ER defaults to 0.0
129# VR defaults to 1.0
130
131oldname = 'SquareWellStructure'
132oldpars = dict()
133demo = dict(effect_radius=50, volfraction=0.04, welldepth=1.5,
134            wellwidth=1.2, effect_radius_pd=0, effect_radius_pd_n=0)
135#
136tests = [
137        [ {'scale': 1.0, 'background' : 0.0, 'effect_radius' : 50.0, 'volfraction' : 0.04,'welldepth' : 1.5, 'wellwidth' : 1.2, 
138           'effect_radius_pd' : 0}, [0.001], [0.97665742]]
139        ]
140# ADDED by: converting from sasview RKH  ON: 16Mar2016 - in progress
141
Note: See TracBrowser for help on using the repository browser.