Changes in / [bf94e6e:9616dfe] in sasmodels


Ignore:
Location:
sasmodels/models
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/core_shell_parallelepiped.c

    re077231 rdbf1a60  
    5959 
    6060    // outer integral (with gauss points), integration limits = 0, 1 
     61    // substitute d_cos_alpha for sin_alpha d_alpha 
    6162    double outer_sum = 0; //initialize integral 
    6263    for( int i=0; i<GAUSS_N; i++) { 
    6364        const double cos_alpha = 0.5 * ( GAUSS_Z[i] + 1.0 ); 
    6465        const double mu = half_q * sqrt(1.0-cos_alpha*cos_alpha); 
    65  
    66         // inner integral (with gauss points), integration limits = 0, pi/2 
    6766        const double siC = length_c * sas_sinx_x(length_c * cos_alpha * half_q); 
    6867        const double siCt = tC * sas_sinx_x(tC * cos_alpha * half_q); 
     68 
     69        // inner integral (with gauss points), integration limits = 0, 1 
     70        // substitute beta = PI/2 u (so 2/PI * d_(PI/2 * beta) = d_beta) 
    6971        double inner_sum = 0.0; 
    7072        for(int j=0; j<GAUSS_N; j++) { 
    71             const double beta = 0.5 * ( GAUSS_Z[j] + 1.0 ); 
     73            const double u = 0.5 * ( GAUSS_Z[j] + 1.0 ); 
    7274            double sin_beta, cos_beta; 
    73             SINCOS(M_PI_2*beta, sin_beta, cos_beta); 
     75            SINCOS(M_PI_2*u, sin_beta, cos_beta); 
    7476            const double siA = length_a * sas_sinx_x(length_a * mu * sin_beta); 
    7577            const double siB = length_b * sas_sinx_x(length_b * mu * cos_beta); 
    … …  
    9193            inner_sum += GAUSS_W[j] * f * f; 
    9294        } 
     95        // now complete change of inner integration variable (1-0)/(1-(-1))= 0.5 
    9396        inner_sum *= 0.5; 
    9497        // now sum up the outer integral 
    9598        outer_sum += GAUSS_W[i] * inner_sum; 
    9699    } 
     100    // now complete change of outer integration variable (1-0)/(1-(-1))= 0.5 
    97101    outer_sum *= 0.5; 
    98102 

  • sasmodels/models/core_shell_parallelepiped.py

    r97be877 rdbf1a60  
    1111.. math:: 
    1212 
    13     I(q) = \text{scale}\frac{\langle f^2 \rangle}{V} + \text{background} 
     13    I(q) = \frac{\text{scale}}{V} \langle P(q,\alpha,\beta) \rangle  
     14    + \text{background} 
    1415 
    1516where $\langle \ldots \rangle$ is an average over all possible orientations 
    16 of the rectangular solid. 
    17  
    18 The function calculated is the form factor of the rectangular solid below. 
     17of the rectangular solid, and the usual $\Delta \rho^2 \ V^2$ term cannot be 
     18pulled out of the form factor term due to the multiple slds in the model. 
     19 
    1920The core of the solid is defined by the dimensions $A$, $B$, $C$ such that 
    2021$A < B < C$. 
    2122 
    22 .. image:: img/core_shell_parallelepiped_geometry.jpg 
     23.. figure:: img/parallelepiped_geometry.jpg 
     24 
     25   Core of the core shell Parallelepiped with the corresponding definition 
     26   of sides. 
     27 
    2328 
    2429There are rectangular "slabs" of thickness $t_A$ that add to the $A$ dimension 
    2530(on the $BC$ faces). There are similar slabs on the $AC$ $(=t_B)$ and $AB$ 
    26 $(=t_C)$ faces. The projection in the $AB$ plane is then 
    27  
    28 .. image:: img/core_shell_parallelepiped_projection.jpg 
    29  
    30 The volume of the solid is 
     31$(=t_C)$ faces. The projection in the $AB$ plane is 
     32 
     33.. figure:: img/core_shell_parallelepiped_projection.jpg 
     34 
     35   AB cut through the core-shell parllelipiped showing the cross secion of 
     36   four of the six shell slabs. As can be seen This model leaves **"gaps"** 
     37   at the corners of the solid. 
     38 
     39 
     40The total volume of the solid is thus given as 
    3141 
    3242.. math:: 
    3343 
    3444    V = ABC + 2t_ABC + 2t_BAC + 2t_CAB 
    35  
    36 **meaning that there are "gaps" at the corners of the solid.** 
    3745 
    3846The intensity calculated follows the :ref:`parallelepiped` model, with the 
    3947core-shell intensity being calculated as the square of the sum of the 
    40 amplitudes of the core and the slabs on the edges. 
    41  
    42 the scattering amplitude is computed for a particular orientation of the 
    43 core-shell parallelepiped with respect to the scattering vector and then 
    44 averaged over all possible orientations, where $\alpha$ is the angle between 
    45 the $z$ axis and the $C$ axis of the parallelepiped, $\beta$ is 
    46 the angle between projection of the particle in the $xy$ detector plane 
    47 and the $y$ axis. 
    48  
    49 .. math:: 
    50  
    51     F(Q) 
     48amplitudes of the core and the slabs on the edges. The scattering amplitude is 
     49computed for a particular orientation of the core-shell parallelepiped with 
     50respect to the scattering vector and then averaged over all possible 
     51orientations, where $\alpha$ is the angle between the $z$ axis and the $C$ axis 
     52of the parallelepiped, and $\beta$ is the angle between the projection of the 
     53particle in the $xy$ detector plane and the $y$ axis. 
     54 
     55.. math:: 
     56 
     57    P(q)=\frac {\int_{0}^{\pi/2}\int_{0}^{\pi/2}F^2(q,\alpha,\beta) \ sin\alpha 
     58    \ d\alpha \ d\beta} {\int_{0}^{\pi/2} \ sin\alpha \ d\alpha \ d\beta} 
     59 
     60and 
     61 
     62.. math:: 
     63 
     64    F(q,\alpha,\beta) 
    5265    &= (\rho_\text{core}-\rho_\text{solvent}) 
    5366       S(Q_A, A) S(Q_B, B) S(Q_C, C) \\ 
    5467    &+ (\rho_\text{A}-\rho_\text{solvent}) 
    55         \left[S(Q_A, A+2t_A) - S(Q_A, Q)\right] S(Q_B, B) S(Q_C, C) \\ 
     68        \left[S(Q_A, A+2t_A) - S(Q_A, A)\right] S(Q_B, B) S(Q_C, C) \\ 
    5669    &+ (\rho_\text{B}-\rho_\text{solvent}) 
    5770        S(Q_A, A) \left[S(Q_B, B+2t_B) - S(Q_B, B)\right] S(Q_C, C) \\ 
    … …  
    6376.. math:: 
    6477 
    65     S(Q, L) = L \frac{\sin \tfrac{1}{2} Q L}{\tfrac{1}{2} Q L} 
     78    S(Q_X, L) = L \frac{\sin \tfrac{1}{2} Q_X L}{\tfrac{1}{2} Q_X L} 
    6679 
    6780and 
    … …  
    6982.. math:: 
    7083 
    71     Q_A &= \sin\alpha \sin\beta \\ 
    72     Q_B &= \sin\alpha \cos\beta \\ 
    73     Q_C &= \cos\alpha 
     84    Q_A &= q \sin\alpha \sin\beta \\ 
     85    Q_B &= q \sin\alpha \cos\beta \\ 
     86    Q_C &= q \cos\alpha 
    7487 
    7588 
    … …  
    7992is the scattering length of the solvent. 
    8093 
     94.. note::  
     95 
     96   the code actually implements two substitutions: $d(cos\alpha)$ is 
     97   substituted for -$sin\alpha \ d\alpha$ (note that in the 
     98   :ref:`parallelepiped` code this is explicitly implemented with 
     99   $\sigma = cos\alpha$), and $\beta$ is set to $\beta = u \pi/2$ so that 
     100   $du = \pi/2 \ d\beta$.  Thus both integrals go from 0 to 1 rather than 0 
     101   to $\pi/2$. 
     102 
    81103FITTING NOTES 
    82104~~~~~~~~~~~~~ 
    … …  
    95117and length $(C+2t_C)$ values, after appropriately sorting the three dimensions 
    96118to give an oblate or prolate particle, to give an effective radius, 
    97 for $S(Q)$ when $P(Q) * S(Q)$ is applied. 
     119for $S(q)$ when $P(q) * S(q)$ is applied. 
    98120 
    99121For 2d data the orientation of the particle is required, described using 
    … …  
    103125$\Psi = 0$ when the *short_b* axis is parallel to the *x*-axis of the detector. 
    104126 
    105 For 2d, constraints must be applied during fitting to ensure that the 
    106 inequality $A < B < C$ is not violated, and hence the correct definition 
    107 of angles is preserved. The calculation will not report an error, 
    108 but the results may be not correct. 
     127.. note:: For 2d, constraints must be applied during fitting to ensure that the 
     128   inequality $A < B < C$ is not violated, and hence the correct definition 
     129   of angles is preserved. The calculation will not report an error, 
     130   but the results may be not correct. 
    109131 
    110132.. figure:: img/parallelepiped_angle_definition.png 

  • sasmodels/models/parallelepiped.c

    r108e70e rdbf1a60  
    3838            inner_total += GAUSS_W[j] * square(si1 * si2); 
    3939        } 
     40        // now complete change of inner integration variable (1-0)/(1-(-1))= 0.5 
    4041        inner_total *= 0.5; 
    4142 
    … …  
    4344        outer_total += GAUSS_W[i] * inner_total * si * si; 
    4445    } 
     46    // now complete change of outer integration variable (1-0)/(1-(-1))= 0.5 
    4547    outer_total *= 0.5; 
    4648 

  • sasmodels/models/parallelepiped.py

    ref07e95 rdbf1a60  
    1010 
    1111 This model calculates the scattering from a rectangular parallelepiped 
    12  (\:numref:`parallelepiped-image`\). 
     12 (:numref:`parallelepiped-image`). 
    1313 If you need to apply polydispersity, see also :ref:`rectangular-prism`. 
    1414 
    … …  
    3939 
    4040    I(q) = \frac{\text{scale}}{V} (\Delta\rho \cdot V)^2 
    41            \left< P(q, \alpha) \right> + \text{background} 
     41           \left< P(q, \alpha, \beta) \right> + \text{background} 
    4242 
    4343where the volume $V = A B C$, the contrast is defined as 
    44 $\Delta\rho = \rho_\text{p} - \rho_\text{solvent}$, 
    45 $P(q, \alpha)$ is the form factor corresponding to a parallelepiped oriented 
    46 at an angle $\alpha$ (angle between the long axis C and $\vec q$), 
    47 and the averaging $\left<\ldots\right>$ is applied over all orientations. 
     44$\Delta\rho = \rho_\text{p} - \rho_\text{solvent}$, $P(q, \alpha, \beta)$ 
     45is the form factor corresponding to a parallelepiped oriented 
     46at an angle $\alpha$ (angle between the long axis C and $\vec q$), and $\beta$ 
     47( the angle between the projection of the particle in the $xy$ detector plane 
     48and the $y$ axis) and the averaging $\left<\ldots\right>$ is applied over all 
     49orientations. 
    4850 
    4951Assuming $a = A/B < 1$, $b = B /B = 1$, and $c = C/B > 1$, the 
    50 form factor is given by (Mittelbach and Porod, 1961) 
     52form factor is given by (Mittelbach and Porod, 1961 [#Mittelbach]_) 
    5153 
    5254.. math:: 
    … …  
    6668    \mu &= qB 
    6769 
    68 The scattering intensity per unit volume is returned in units of |cm^-1|. 
     70where substitution of $\sigma = cos\alpha$ and $\beta = \pi/2 \ u$ have been 
     71applied. 
    6972 
    7073NB: The 2nd virial coefficient of the parallelepiped is calculated based on 
    … …  
    120123.. math:: 
    121124 
    122     P(q_x, q_y) = \left[\frac{\sin(\tfrac{1}{2}qA\cos\alpha)}{(\tfrac{1}{2}qA\cos\alpha)}\right]^2 
    123                   \left[\frac{\sin(\tfrac{1}{2}qB\cos\beta)}{(\tfrac{1}{2}qB\cos\beta)}\right]^2 
    124                   \left[\frac{\sin(\tfrac{1}{2}qC\cos\gamma)}{(\tfrac{1}{2}qC\cos\gamma)}\right]^2 
     125    P(q_x, q_y) = \left[\frac{\sin(\tfrac{1}{2}qA\cos\alpha)}{(\tfrac{1} 
     126                   {2}qA\cos\alpha)}\right]^2 
     127                  \left[\frac{\sin(\tfrac{1}{2}qB\cos\beta)}{(\tfrac{1} 
     128                   {2}qB\cos\beta)}\right]^2 
     129                  \left[\frac{\sin(\tfrac{1}{2}qC\cos\gamma)}{(\tfrac{1} 
     130                   {2}qC\cos\gamma)}\right]^2 
    125131 
    126132with 
    … …  
    160166---------- 
    161167 
    162 P Mittelbach and G Porod, *Acta Physica Austriaca*, 14 (1961) 185-211 
    163  
    164 R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
     168.. [#Mittelbach] P Mittelbach and G Porod, *Acta Physica Austriaca*, 
     169   14 (1961) 185-211 
     170.. [#] R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
    165171 
    166172Authorship and Verification 
Note: See TracChangeset for help on using the changeset viewer.