source: sasmodels/sasmodels/models/stacked_disks.c @ becded3

core_shell_microgelsmagnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since becded3 was becded3, checked in by Paul Kienzle <pkienzle@…>, 6 years ago

update oriented models to new interface (which will be in the next commit)

  • Property mode set to 100644
File size: 4.7 KB
Line 
1static double
2stacked_disks_kernel(
3    double qab,
4    double qc,
5    double halfheight,
6    double thick_layer,
7    double radius,
8    int n_stacking,
9    double sigma_dnn,
10    double core_sld,
11    double layer_sld,
12    double solvent_sld,
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 = radius*qab;
23    //const double besarg2 = radius*qab;
24
25    const double sinarg1 = halfheight*qc;
26    const double sinarg2 = (halfheight+thick_layer)*qc;
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 = d*qc;
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
64stacked_disks_1d(
65    double q,
66    double thick_core,
67    double thick_layer,
68    double radius,
69    int n_stacking,
70    double sigma_dnn,
71    double core_sld,
72    double layer_sld,
73    double solvent_sld)
74{
75/*    StackedDiscsX  :  calculates the form factor of a stacked "tactoid" of core shell disks
76like clay platelets that are not exfoliated
77*/
78    double summ = 0.0;    //initialize integral
79
80    double d = 2.0*thick_layer+thick_core;
81    double halfheight = 0.5*thick_core;
82
83    for(int i=0; i<N_POINTS_76; i++) {
84        double zi = (Gauss76Z[i] + 1.0)*M_PI_4;
85        double sin_alpha, cos_alpha; // slots to hold sincos function output
86        SINCOS(zi, sin_alpha, cos_alpha);
87        double yyy = stacked_disks_kernel(q*sin_alpha, q*cos_alpha,
88                           halfheight,
89                           thick_layer,
90                           radius,
91                           n_stacking,
92                           sigma_dnn,
93                           core_sld,
94                           layer_sld,
95                           solvent_sld,
96                           d);
97        summ += Gauss76Wt[i] * yyy * sin_alpha;
98    }
99
100    double answer = M_PI_4*summ;
101
102    //Convert to [cm-1]
103    return 1.0e-4*answer;
104}
105
106static double
107form_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
119Iq(
120    double q,
121    double thick_core,
122    double thick_layer,
123    double radius,
124    double fp_n_stacking,
125    double sigma_dnn,
126    double core_sld,
127    double layer_sld,
128    double solvent_sld)
129{
130    int n_stacking = (int)(fp_n_stacking + 0.5);
131    return stacked_disks_1d(q,
132                    thick_core,
133                    thick_layer,
134                    radius,
135                    n_stacking,
136                    sigma_dnn,
137                    core_sld,
138                    layer_sld,
139                    solvent_sld);
140}
141
142
143static double
144Iqxy(double qab, double qc,
145    double thick_core,
146    double thick_layer,
147    double radius,
148    double fp_n_stacking,
149    double sigma_dnn,
150    double core_sld,
151    double layer_sld,
152    double solvent_sld)
153{
154    int n_stacking = (int)(fp_n_stacking + 0.5);
155    double d = 2.0 * thick_layer + thick_core;
156    double halfheight = 0.5*thick_core;
157    double answer = stacked_disks_kernel(qab, qc,
158                     halfheight,
159                     thick_layer,
160                     radius,
161                     n_stacking,
162                     sigma_dnn,
163                     core_sld,
164                     layer_sld,
165                     solvent_sld,
166                     d);
167
168    //convert to [cm-1]
169    answer *= 1.0e-4;
170
171    return answer;
172}
173
Note: See TracBrowser for help on using the repository browser.