source: sasmodels/sasmodels/models/stacked_disks.c @ 6831fa0

core_shell_microgelscostrafo411magnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 6831fa0 was 6831fa0, checked in by Paul Kienzle <pkienzle@…>, 8 years ago

stacked_disks: enable 2D modeling for stacked disks
  • Property mode set to 100644
File size: 5.4 KB
Line 
1double form_volume(double thick_core,
2                   double thick_layer,
3                   double radius,
4                   double n_stacking);
5
6double Iq(double q,
7          double thick_core,
8          double thick_layer,
9          double radius,
10          double n_stacking,
11          double sigma_dnn,
12          double core_sld,
13          double layer_sld,
14          double solvent_sld);
15
16double Iqxy(double qx, double qy,
17          double thick_core,
18          double thick_layer,
19          double radius,
20          double n_stacking,
21          double sigma_dnn,
22          double core_sld,
23          double layer_sld,
24          double solvent_sld,
25          double theta,
26          double phi);
27
28static
29double _kernel(double qq,
30               double radius,
31               double core_sld,
32               double layer_sld,
33               double solvent_sld,
34               double halfheight,
35               double thick_layer,
36               double sin_alpha,
37               double cos_alpha,
38               double sigma_dnn,
39               double d,
40               double n_stacking)
41
42{
43        // qq is the q-value for the calculation (1/A)
44        // radius is the core radius of the cylinder (A)
45        // *_sld are the respective SLD's
46        // halfheight is the *Half* CORE-LENGTH of the cylinder = L (A)
47        // zi is the dummy variable for the integration (x in Feigin's notation)
48
49        const double besarg1 = qq*radius*sin_alpha;
50        //const double besarg2 = qq*radius*sin_alpha;
51
52        const double sinarg1 = qq*halfheight*cos_alpha;
53        const double sinarg2 = qq*(halfheight+thick_layer)*cos_alpha;
54
55        const double be1 = sas_J1c(besarg1);
56        //const double be2 = sas_J1c(besarg2);
57        const double be2 = be1;
58        const double si1 = sinc(sinarg1);
59        const double si2 = sinc(sinarg2);
60
61        const double dr1 = (core_sld-solvent_sld);
62        const double dr2 = (layer_sld-solvent_sld);
63        const double area = M_PI*radius*radius;
64        const double totald = 2.0*(thick_layer+halfheight);
65
66        const double t1 = area*(2.0*halfheight)*dr1*(si1)*(be1);
67        const double t2 = area*dr2*(totald*si2-2.0*halfheight*si1)*(be2);
68
69
70        double retval =((t1+t2)*(t1+t2));
71
72        // loop for the structure facture S(q)
73        double sqq=0.0;
74        for(int kk=1;kk<n_stacking;kk+=1) {
75                double qd_cos_alpha = qq*d*cos_alpha;
76                double dexpt=square(qd_cos_alpha*sigma_dnn)*kk/2.0;
77                sqq += (n_stacking-kk)*cos(qd_cos_alpha*kk)*exp(-1.*dexpt);
78        }
79        // end of loop for S(q)
80        sqq=1.0+2.0*sqq/n_stacking;
81
82        return retval * sqq;
83}
84
85
86static
87double stacked_disks_kernel(double q,
88                            double thick_core,
89                            double thick_layer,
90                            double radius,
91                            double n_stacking,
92                            double sigma_dnn,
93                            double core_sld,
94                            double layer_sld,
95                            double solvent_sld)
96{
97/*      StackedDiscsX  :  calculates the form factor of a stacked "tactoid" of core shell disks
98like clay platelets that are not exfoliated
99*/
100        double summ = 0.0;      //initialize integral
101
102        double d = 2.0*thick_layer+thick_core;
103        double halfheight = 0.5*thick_core;
104
105        for(int i=0; i<N_POINTS_76; i++) {
106                double zi = (Gauss76Z[i] + 1.0)*M_PI_4;
107                double sin_alpha, cos_alpha; // slots to hold sincos function output
108                SINCOS(zi, sin_alpha, cos_alpha);
109                double yyy = _kernel(q,
110                                   radius,
111                                   core_sld,
112                                   layer_sld,
113                                   solvent_sld,
114                                   halfheight,
115                                   thick_layer,
116                                   sin_alpha,
117                                   cos_alpha,
118                                   sigma_dnn,
119                                   d,
120                                   n_stacking);
121                summ += Gauss76Wt[i] * yyy * sin_alpha;
122        }
123
124        double answer = M_PI_4*summ;
125
126        //Convert to [cm-1]
127        answer *= 1.0e-4;
128
129        return answer;
130}
131
132double form_volume(double thick_core,
133                   double thick_layer,
134                   double radius,
135                   double n_stacking){
136    double d = 2.0 * thick_layer + thick_core;
137    return M_PI * radius * radius * d * n_stacking;
138}
139
140double Iq(double q,
141          double thick_core,
142          double thick_layer,
143          double radius,
144          double n_stacking,
145          double sigma_dnn,
146          double core_sld,
147          double layer_sld,
148          double solvent_sld)
149{
150    return stacked_disks_kernel(q,
151                    thick_core,
152                    thick_layer,
153                    radius,
154                    n_stacking,
155                    sigma_dnn,
156                    core_sld,
157                    layer_sld,
158                    solvent_sld);
159}
160
161
162double
163Iqxy(double qx, double qy,
164     double thick_core,
165     double thick_layer,
166     double radius,
167     double n_stacking,
168     double sigma_dnn,
169     double core_sld,
170     double layer_sld,
171     double solvent_sld,
172     double theta,
173     double phi)
174{
175    double q, sin_alpha, cos_alpha;
176    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha);
177
178    double d = 2.0 * thick_layer + thick_core;
179    double halfheight = 0.5*thick_core;
180    double answer = _kernel(q,
181                     radius,
182                     core_sld,
183                     layer_sld,
184                     solvent_sld,
185                     halfheight,
186                     thick_layer,
187                     sin_alpha,
188                     cos_alpha,
189                     sigma_dnn,
190                     d,
191                     n_stacking);
192
193    //convert to [cm-1]
194    answer *= 1.0e-4;
195
196    return answer;
197}
198
Note: See TracBrowser for help on using the repository browser.