1 | |
---|
2 | static double |
---|
3 | f_exp(double q, double r, double sld_in, double sld_out, |
---|
4 | double thickness, double A, double side) |
---|
5 | { |
---|
6 | const double vol = M_4PI_3 * cube(r); |
---|
7 | const double qr = q * r; |
---|
8 | const double bes = sas_3j1x_x(qr); |
---|
9 | const double alpha = A * r/thickness; |
---|
10 | double result; |
---|
11 | if (qr == 0.0) { |
---|
12 | result = 1.0; |
---|
13 | } else if (fabs(A) > 0.0) { |
---|
14 | const double qrsq = qr * qr; |
---|
15 | const double alphasq = alpha * alpha; |
---|
16 | const double sumsq = alphasq + qrsq; |
---|
17 | double sinqr, cosqr; |
---|
18 | SINCOS(qr, sinqr, cosqr); |
---|
19 | const double t1 = (alphasq - qrsq)*sinqr/qr - 2.0*alpha*cosqr; |
---|
20 | const double t2 = alpha*sinqr/qr - cosqr; |
---|
21 | const double fun = -3.0*(t1/sumsq - t2)/sumsq; |
---|
22 | const double slope = (sld_out - sld_in)/expm1(A); |
---|
23 | const double contrast = slope*exp(A*side); |
---|
24 | result = contrast*fun + (sld_in-slope)*bes; |
---|
25 | } else { |
---|
26 | result = sld_in*bes; |
---|
27 | } |
---|
28 | return vol * result; |
---|
29 | } |
---|
30 | |
---|
31 | static double |
---|
32 | form_volume(double radius_core, double n_shells, double thickness[]) |
---|
33 | { |
---|
34 | int n = (int)(n_shells+0.5); |
---|
35 | double r = radius_core; |
---|
36 | for (int i=0; i < n; i++) { |
---|
37 | r += thickness[i]; |
---|
38 | } |
---|
39 | return M_4PI_3*cube(r); |
---|
40 | } |
---|
41 | |
---|
42 | static void |
---|
43 | Fq(double q, double *F1, double *F2, double sld_core, double radius_core, double sld_solvent, |
---|
44 | double n_shells, double sld_in[], double sld_out[], double thickness[], |
---|
45 | double A[]) |
---|
46 | { |
---|
47 | int n = (int)(n_shells+0.5); |
---|
48 | double r_out = radius_core; |
---|
49 | double f = f_exp(q, r_out, sld_core, 0.0, 0.0, 0.0, 0.0); |
---|
50 | for (int i=0; i < n; i++){ |
---|
51 | const double r_in = r_out; |
---|
52 | r_out += thickness[i]; |
---|
53 | f -= f_exp(q, r_in, sld_in[i], sld_out[i], thickness[i], A[i], 0.0); |
---|
54 | f += f_exp(q, r_out, sld_in[i], sld_out[i], thickness[i], A[i], 1.0); |
---|
55 | } |
---|
56 | f -= f_exp(q, r_out, sld_solvent, 0.0, 0.0, 0.0, 0.0); |
---|
57 | |
---|
58 | *F1 = 1e-2 * f; |
---|
59 | *F2 = 1e-4 * f * f; |
---|
60 | } |
---|