source: sasmodels/sasmodels/models/mass_fractal.c @ d507c3a

core_shell_microgelscostrafo411magnetic_modelrelease_v0.94release_v0.95ticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since d507c3a was 61fd21d, checked in by piotr, 8 years ago

Actively check the right value of mass_dim to avoid NaNs?.

  • Property mode set to 100644
File size: 1.5 KB
Line 
1double form_volume(double radius);
2
3double Iq(double q,
4          double radius,
5          double mass_dim,
6          double cutoff_length);
7
8double Iqxy(double qx, double qy,
9          double radius,
10          double mass_dim,
11          double cutoff_length);
12
13
14static double _mass_fractal_kernel(double q,
15          double radius,
16          double mass_dim,
17          double cutoff_length)
18{
19    // Actively check the argument.
20    if (mass_dim <= 1.0){
21       return 0.0;
22    }
23
24    //calculate P(q)
25    double pq = sph_j1c(q*radius);
26    pq = pq*pq;
27
28    //calculate S(q)
29    double mmo = mass_dim-1.0;
30    double sq = exp(lanczos_gamma(mmo))*sin((mmo)*atan(q*cutoff_length));
31    sq *= pow(cutoff_length, mmo);
32    sq /= pow((1.0 + (q*cutoff_length)*(q*cutoff_length)),(mmo/2.0));
33    sq /= q;
34
35    //combine and return
36    double result = pq * sq;
37
38    return result;
39}
40double form_volume(double radius){
41
42    return 1.333333333333333*M_PI*radius*radius*radius;
43}
44
45double Iq(double q,
46          double radius,
47          double mass_dim,
48          double cutoff_length)
49{
50    return _mass_fractal_kernel(q,
51           radius,
52           mass_dim,
53           cutoff_length);
54}
55
56// Iqxy is never called since no orientation or magnetic parameters.
57double Iqxy(double qx, double qy,
58          double radius,
59          double mass_dim,
60          double cutoff_length)
61{
62    double q = sqrt(qx*qx + qy*qy);
63    return _mass_fractal_kernel(q,
64           radius,
65           mass_dim,
66           cutoff_length);
67}
68
Note: See TracBrowser for help on using the repository browser.