source: sasmodels/sasmodels/models/stacked_disks.c @ 19f996b

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

stacked disks: remove unnecessary int→float implicit conversion

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