source: sasmodels/sasmodels/models/core_shell_bicelle.c @ 2a0b2b1

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

restructure all 2D models to work with (qa,qb,qc) = rotate(qx,qy) rather than working with angles directly in preparation for revised jitter algorithm

  • Property mode set to 100644
File size: 2.6 KB
Line 
1static double
2form_volume(double radius, double thick_rim, double thick_face, double length)
3{
4    return M_PI*square(radius+thick_rim)*(length+2.0*thick_face);
5}
6
7static double
8bicelle_kernel(double qab,
9    double qc,
10    double radius,
11    double thick_radius,
12    double thick_face,
13    double halflength,
14    double sld_core,
15    double sld_face,
16    double sld_rim,
17    double sld_solvent)
18{
19    const double dr1 = sld_core-sld_face;
20    const double dr2 = sld_rim-sld_solvent;
21    const double dr3 = sld_face-sld_rim;
22    const double vol1 = M_PI*square(radius)*2.0*(halflength);
23    const double vol2 = M_PI*square(radius+thick_radius)*2.0*(halflength+thick_face);
24    const double vol3 = M_PI*square(radius)*2.0*(halflength+thick_face);
25
26    const double be1 = sas_2J1x_x((radius)*qab);
27    const double be2 = sas_2J1x_x((radius+thick_radius)*qab);
28    const double si1 = sas_sinx_x((halflength)*qc);
29    const double si2 = sas_sinx_x((halflength+thick_face)*qc);
30
31    const double t = vol1*dr1*si1*be1 +
32                     vol2*dr2*si2*be2 +
33                     vol3*dr3*si2*be1;
34
35    return t;
36}
37
38static double
39Iq(double q,
40    double radius,
41    double thick_radius,
42    double thick_face,
43    double length,
44    double sld_core,
45    double sld_face,
46    double sld_rim,
47    double sld_solvent)
48{
49    // set up the integration end points
50    const double uplim = M_PI_4;
51    const double halflength = 0.5*length;
52
53    double total = 0.0;
54    for(int i=0;i<N_POINTS_76;i++) {
55        double theta = (Gauss76Z[i] + 1.0)*uplim;
56        double sin_theta, cos_theta; // slots to hold sincos function output
57        SINCOS(theta, sin_theta, cos_theta);
58        double fq = bicelle_kernel(q*sin_theta, q*cos_theta, radius, thick_radius, thick_face,
59                                   halflength, sld_core, sld_face, sld_rim, sld_solvent);
60        total += Gauss76Wt[i]*fq*fq*sin_theta;
61    }
62
63    // calculate value of integral to return
64    double answer = total*uplim;
65    return 1.0e-4*answer;
66}
67
68static double
69Iqxy(double qx, double qy,
70    double radius,
71    double thick_rim,
72    double thick_face,
73    double length,
74    double core_sld,
75    double face_sld,
76    double rim_sld,
77    double solvent_sld,
78    double theta,
79    double phi)
80{
81    double q, sin_alpha, cos_alpha;
82    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha);
83    const double qab = q*sin_alpha;
84    const double qc = q*cos_alpha;
85
86    double fq = bicelle_kernel(qab, qc, radius, thick_rim, thick_face,
87                           0.5*length, core_sld, face_sld, rim_sld,
88                           solvent_sld);
89    return 1.0e-4*fq*fq;
90}
Note: See TracBrowser for help on using the repository browser.