Changes in / [6196487:d24e390] in sasmodels


Ignore:
Location:
sasmodels
Files:
27 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/compare.py

    rd504bcd rfe25eda  
    340340        if pars['radius'] < pars['thick_string']: 
    341341            pars['radius'], pars['thick_string'] = pars['thick_string'], pars['radius'] 
     342        pars['num_pearls'] = math.ceil(pars['num_pearls']) 
    342343        pass 
    343344 
    … …  
    352353        for c in '1234': 
    353354            pars['Phi'+c] /= total 
     355 
     356    elif name == 'stacked_disks': 
     357        pars['n_stacking'] = math.ceil(pars['n_stacking']) 
    354358 
    355359def parlist(model_info, pars, is2d): 

  • sasmodels/models/core_multi_shell.c

    r925ad6e r925ad6e  
    88} 
    99 
    10 static double 
    11 form_volume(double core_radius, double fp_n, double thickness[]) 
     10double 
     11form_volume(double core_radius, double n, double thickness[]); 
     12double 
     13form_volume(double core_radius, double n, double thickness[]) 
    1214{ 
    1315  double r = core_radius; 
    14   int n = (int)(fp_n+0.5); 
    1516  for (int i=0; i < n; i++) { 
    1617    r += thickness[i]; 
    … …  
    1920} 
    2021 
    21 static double 
     22double 
    2223Iq(double q, double core_sld, double core_radius, 
    23    double solvent_sld, double fp_n, double sld[], double thickness[]) 
     24   double solvent_sld, double num_shells, double sld[], double thickness[]); 
     25double 
     26Iq(double q, double core_sld, double core_radius, 
     27   double solvent_sld, double num_shells, double sld[], double thickness[]) 
    2428{ 
    25   const int n = (int)(fp_n+0.5); 
     29  const int n = (int)ceil(num_shells); 
    2630  double f, r, last_sld; 
    2731  r = core_radius; 

  • sasmodels/models/core_multi_shell.py

    r5a0b3d7 r925ad6e  
    107107    Returns the SLD profile *r* (Ang), and *rho* (1e-6/Ang^2). 
    108108    """ 
    109     n = int(n+0.5) 
    110109    z = [] 
    111110    rho = [] 
    … …  
    134133def ER(radius, n, thickness): 
    135134    """Effective radius""" 
    136     n = int(n[0]+0.5)  # n is a control parameter and is not polydisperse 
     135    n = int(n[0])  # n cannot be polydisperse 
    137136    return np.sum(thickness[:n], axis=0) + radius 
    138137 

  • sasmodels/models/flexible_cylinder.c

    r592343f r592343f  
    11static double 
    2 form_volume(double length, double kuhn_length, double radius) 
     2form_volume(length, kuhn_length, radius) 
    33{ 
    44    return 1.0; 

  • sasmodels/models/lamellar_hg_stack_caille.c

    r1c6286d ra807206  
    33*/ 
    44 
    5 static double 
    6 Iq(double qval, 
    7    double length_tail, 
    8    double length_head, 
    9    double fp_Nlayers, 
    10    double dd, 
    11    double Cp, 
    12    double tail_sld, 
    13    double head_sld, 
    14    double solvent_sld) 
     5double Iq(double qval, 
     6      double length_tail, 
     7      double length_head, 
     8      double Nlayers,  
     9      double dd, 
     10      double Cp, 
     11      double tail_sld, 
     12      double head_sld, 
     13      double solvent_sld); 
     14 
     15double Iq(double qval, 
     16      double length_tail, 
     17      double length_head, 
     18      double Nlayers,  
     19      double dd, 
     20      double Cp, 
     21      double tail_sld, 
     22      double head_sld, 
     23      double solvent_sld) 
    1524{ 
    16   int Nlayers = (int)(fp_Nlayers+0.5);    //cast to an integer for the loop 
     25  double NN;   //local variables of coefficient wave 
    1726  double inten,Pq,Sq,alpha,temp,t2; 
    1827  //double dQ, dQDefault, t1, t3; 
     28  int ii,NNint; 
    1929  // from wikipedia 0.577215664901532860606512090082402431042159335 
    2030  const double Euler = 0.577215664901533;   // Euler's constant, increased sig figs for new models Feb 2015 
    … …  
    2232  //dQ = dQDefault; // REMOVED UNUSED dQ calculations for new models Feb 2015 
    2333 
    24   Pq = (head_sld-solvent_sld)*(sin(qval*(length_head+length_tail))-sin(qval*length_tail)) 
    25        + (tail_sld-solvent_sld)*sin(qval*length_tail); 
     34  NN = trunc(Nlayers);    //be sure that NN is an integer 
     35   
     36  Pq = (head_sld-solvent_sld)*(sin(qval*(length_head+length_tail))-sin(qval*length_tail)) + 
     37              (tail_sld-solvent_sld)*sin(qval*length_tail); 
    2638  Pq *= Pq; 
    2739  Pq *= 4.0/(qval*qval); 
    2840 
     41  NNint = (int)NN;    //cast to an integer for the loop 
     42  ii=0; 
    2943  Sq = 0.0; 
    30   for(int ii=1; ii < Nlayers; ii++) { 
     44  for(ii=1;ii<=(NNint-1);ii+=1) { 
     45 
     46    //fii = (double)ii;   //do I really need to do this? - unused variable, removed 18Feb2015 
     47 
    3148    temp = 0.0; 
    3249    alpha = Cp/4.0/M_PI/M_PI*(log(M_PI*ii) + Euler); 
    … …  
    3552    //t3 = dQ*dQ*dd*dd*ii*ii; 
    3653 
    37     temp = 1.0-(double)ii/(double)Nlayers; 
     54    temp = 1.0-ii/NN; 
    3855    //temp *= cos(dd*qval*ii/(1.0+t1)); 
    3956    temp *= cos(dd*qval*ii); 
    4057    //if (temp < 0) printf("q=%g: ii=%d, cos(dd*q*ii)=cos(%g) < 0\n",qval,ii,dd*qval*ii); 
    4158    //temp *= exp(-1.0*(t2 + t3)/(2.0*(1.0+t1)) ); 
    42     temp *= exp(-t2/2.0); 
     59    temp *= exp(-t2/2.0 ); 
    4360    //temp /= sqrt(1.0+t1); 
    4461 
    … …  
    5471 
    5572  inten *= 1.0e-04;   // 1/A to 1/cm 
    56   return inten; 
     73  return(inten); 
    5774} 
    5875 

  • sasmodels/models/lamellar_hg_stack_caille.py

    ra57b31d r7c57861  
    9898    ["length_head", "Ang", 2, [0, inf], "volume", 
    9999     "head thickness"], 
    100     ["Nlayers", "", 30, [1, inf], "", 
     100    ["Nlayers", "", 30, [0, inf], "", 
    101101     "Number of layers"], 
    102102    ["d_spacing", "Ang", 40., [0.0, inf], "volume", 

  • sasmodels/models/lamellar_stack_caille.c

    r1c6286d r0bef47b  
    33*/ 
    44 
    5 static double 
    6 Iq(double qval, 
    7    double del, 
    8    double fp_Nlayers, 
    9    double dd, 
    10    double Cp, 
    11    double sld, 
    12    double solvent_sld) 
     5double Iq(double qval, 
     6      double del, 
     7      double Nlayers,  
     8      double dd, 
     9          double Cp,  
     10      double sld, 
     11      double solvent_sld); 
     12 
     13double Iq(double qval, 
     14      double del, 
     15      double Nlayers,  
     16      double dd, 
     17          double Cp,  
     18      double sld, 
     19      double solvent_sld) 
    1320{ 
    14   int Nlayers = (int)(fp_Nlayers+0.5);    //cast to an integer for the loop 
    15   double contr;   //local variables of coefficient wave 
     21  double contr,NN;   //local variables of coefficient wave 
    1622  double inten,Pq,Sq,alpha,temp,t2; 
    1723  //double dQ, dQDefault, t1, t3; 
     24  int ii,NNint; 
    1825  // from wikipedia 0.577215664901532860606512090082402431042159335 
    1926  const double Euler = 0.577215664901533;   // Euler's constant, increased sig figs for new models Feb 2015 
    … …  
    2128  //dQ = dQDefault; // REMOVED UNUSED dQ calculations for new models Feb 2015 
    2229 
     30  NN = trunc(Nlayers);    //be sure that NN is an integer 
     31   
    2332  contr = sld - solvent_sld; 
    2433 
    2534  Pq = 2.0*contr*contr/qval/qval*(1.0-cos(qval*del)); 
    2635 
     36  NNint = (int)NN;    //cast to an integer for the loop 
     37  ii=0; 
    2738  Sq = 0.0; 
    28   for (int ii=1; ii < Nlayers; ii++) { 
     39  // the vital "=" in ii<=  added March 2015 
     40  for(ii=1;ii<=(NNint-1);ii+=1) { 
     41 
     42    //fii = (double)ii;   //do I really need to do this? - unused variable, removed 18Feb2015 
     43 
    2944    temp = 0.0; 
    3045    alpha = Cp/4.0/M_PI/M_PI*(log(M_PI*ii) + Euler); 
    … …  
    3348    //t3 = dQ*dQ*dd*dd*ii*ii; 
    3449 
    35     temp = 1.0 - (double)ii / (double)Nlayers; 
     50    temp = 1.0-ii/NN; 
    3651    //temp *= cos(dd*qval*ii/(1.0+t1)); 
    3752    temp *= cos(dd*qval*ii); 
    3853    //temp *= exp(-1.0*(t2 + t3)/(2.0*(1.0+t1)) ); 
    39     temp *= exp(-t2/2.0); 
     54    temp *= exp(-t2/2.0 ); 
    4055    //temp /= sqrt(1.0+t1); 
    4156 

  • sasmodels/models/lamellar_stack_caille.py

    ra57b31d r7c57861  
    8888parameters = [ 
    8989    ["thickness",        "Ang",      30.0,  [0, inf],   "volume", "sheet thickness"], 
    90     ["Nlayers",          "",          20,   [1, inf],   "",       "Number of layers"], 
     90    ["Nlayers",          "",          20,   [0, inf],   "",       "Number of layers"], 
    9191    ["d_spacing",        "Ang",      400.,  [0.0,inf],  "volume", "lamellar d-spacing of Caille S(Q)"], 
    9292    ["Caille_parameter", "1/Ang^2",    0.1, [0.0,0.8],  "",       "Caille parameter"], 

  • sasmodels/models/lamellar_stack_paracrystal.c

    r5467cd8 r4962519  
    22 
    33*/ 
    4 double paraCryst_sn(double ww, double qval, double davg, int Nlayers, double an); 
    5 double paraCryst_an(double ww, double qval, double davg, int Nlayers); 
     4double Iq(double qval, 
     5      double th, 
     6      double Nlayers,  
     7          double davg,  
     8          double pd, 
     9      double sld, 
     10      double solvent_sld); 
     11double paraCryst_sn(double ww, double qval, double davg, long Nlayers, double an); 
     12double paraCryst_an(double ww, double qval, double davg, long Nlayers); 
    613 
    7 static double 
    8 Iq(double qval, 
    9    double th, 
    10    double fp_Nlayers, 
    11    double davg, 
    12    double pd, 
    13    double sld, 
    14    double solvent_sld) 
     14double Iq(double qval, 
     15      double th, 
     16      double Nlayers,  
     17          double davg,  
     18          double pd, 
     19      double sld, 
     20      double solvent_sld) 
    1521{ 
     22     
     23        double inten,contr,xn; 
     24        double xi,ww,Pbil,Znq,Snq,an; 
     25        long n1,n2; 
     26         
     27        contr = sld - solvent_sld; 
    1628        //get the fractional part of Nlayers, to determine the "mixing" of N's 
    17         int n1 = (int)(fp_Nlayers);             //truncate towards zero 
    18         int n2 = n1 + 1; 
    19         const double xn = (double)n2 - fp_Nlayers;      //fractional contribution of n1 
    2029         
    21         const double ww = exp(-0.5*square(qval*pd*davg)); 
     30        n1 = (long)trunc(Nlayers);              //rounds towards zero 
     31        n2 = n1 + 1; 
     32        xn = (double)n2 - Nlayers;                      //fractional contribution of n1 
     33         
     34        ww = exp(-qval*qval*pd*pd*davg*davg/2.0); 
    2235 
    2336        //calculate the n1 contribution 
    24         double Znq,Snq,an; 
    2537        an = paraCryst_an(ww,qval,davg,n1); 
    2638        Snq = paraCryst_sn(ww,qval,davg,n1,an); 
    27  
     39         
    2840        Znq = xn*Snq; 
    2941         
    … …  
    4052//      Zq = (1-ww^2)/(1+ww^2-2*ww*cos(qval*davg)) 
    4153         
    42         const double xi = th/2.0;               //use 1/2 the bilayer thickness 
    43         const double Pbil = square(sas_sinx_x(qval*xi)); 
     54        xi = th/2.0;            //use 1/2 the bilayer thickness 
     55        Pbil = (sin(qval*xi)/(qval*xi))*(sin(qval*xi)/(qval*xi)); 
    4456         
    45         const double contr = sld - solvent_sld; 
    46         const double inten = 2.0*M_PI*contr*contr*Pbil*Znq/(qval*qval); 
     57        inten = 2.0*M_PI*contr*contr*Pbil*Znq/(qval*qval); 
     58        inten *= 1.0e-04; 
    4759//printf("q=%.7e wwm1=%g ww=%.5e an=% 12.5e Snq=% 12.5e Znq=% 12.5e Pbil=% 12.5e\n",qval,wwm1,ww,an,Snq,Znq,Pbil); 
    48         return 1.0e-4*inten; 
     60        return(inten); 
    4961} 
    5062 
    5163// functions for the lamellar paracrystal model 
    5264double 
    53 paraCryst_sn(double ww, double qval, double davg, int Nlayers, double an) { 
     65paraCryst_sn(double ww, double qval, double davg, long Nlayers, double an) { 
    5466         
    5567        double Snq; 
    … …  
    5769        Snq = an/( (double)Nlayers*square(1.0+ww*ww-2.0*ww*cos(qval*davg)) ); 
    5870         
    59         return Snq; 
     71        return(Snq); 
    6072} 
    6173 
    6274double 
    63 paraCryst_an(double ww, double qval, double davg, int Nlayers) { 
     75paraCryst_an(double ww, double qval, double davg, long Nlayers) { 
     76         
    6477        double an; 
    6578         
    … …  
    6982        an += 2.0*pow(ww,(Nlayers+1))*cos((double)(Nlayers+1)*qval*davg); 
    7083         
    71         return an; 
     84        return(an); 
    7285} 
    7386 

  • sasmodels/models/lamellar_stack_paracrystal.py

    ra0168e8 r7c57861  
    113113parameters = [["thickness", "Ang", 33.0, [0, inf], "volume", 
    114114               "sheet thickness"], 
    115               ["Nlayers", "", 20, [1, inf], "", 
     115              ["Nlayers", "", 20, [0, inf], "", 
    116116               "Number of layers"], 
    117117              ["d_spacing", "Ang", 250., [0.0, inf], "", 

  • sasmodels/models/linear_pearls.c

    r925ad6e r925ad6e  
    44            double radius, 
    55            double edge_sep, 
    6             double fp_num_pearls, 
     6            double num_pearls, 
    77            double pearl_sld, 
    88            double solvent_sld); 
    … …  
    1111            double radius, 
    1212            double edge_sep, 
    13             int num_pearls, 
     13            double num_pearls, 
    1414            double pearl_sld, 
    1515            double solvent_sld); 
    1616 
    1717 
    18 double form_volume(double radius, double fp_num_pearls) 
     18double form_volume(double radius, double num_pearls) 
    1919{ 
    20     int num_pearls = (int)(fp_num_pearls + 0.5); 
    2120    // Pearl volume 
    2221    double pearl_vol = M_4PI_3 * cube(radius); 
    … …  
    2827            double radius, 
    2928            double edge_sep, 
    30             int num_pearls, 
     29            double num_pearls, 
    3130            double pearl_sld, 
    3231            double solvent_sld) 
    3332{ 
     33    double n_contrib; 
    3434    //relative sld 
    3535    double contrast_pearl = pearl_sld - solvent_sld; 
    … …  
    4646    double psi = sas_3j1x_x(q * radius); 
    4747 
    48     // N pearls interaction terms  
    49     double structure_factor = (double)num_pearls; 
    50     for(int num=1; num<num_pearls; num++) { 
    51         structure_factor += 2.0*(num_pearls-num)*sas_sinx_x(q*separation*num); 
     48    // N pearls contribution 
     49    int n_max = num_pearls - 1; 
     50    n_contrib = num_pearls; 
     51    for(int num=1; num<=n_max; num++) { 
     52        n_contrib += (2.0*(num_pearls-num)*sas_sinx_x(q*separation*num)); 
    5253    } 
    5354    // form factor for num_pearls 
    54     double form_factor = 1.0e-4 * structure_factor * square(m_s*psi) / tot_vol; 
     55    double form_factor = 1.0e-4 * n_contrib * square(m_s*psi) / tot_vol; 
    5556 
    5657    return form_factor; 
    … …  
    6061            double radius, 
    6162            double edge_sep, 
    62             double fp_num_pearls, 
     63            double num_pearls, 
    6364            double pearl_sld, 
    6465            double solvent_sld) 
    6566{ 
    6667 
    67     int num_pearls = (int)(fp_num_pearls + 0.5); 
    6868        double result = linear_pearls_kernel(q, 
    6969                    radius, 

  • sasmodels/models/linear_pearls.py

    r925ad6e r925ad6e  
    1616.. math:: 
    1717 
    18     P(Q) = \frac{\text{scale}}{V}\left[ m_{p}^2 
    19     \left(N+2\sum_{n-1}^{N-1}(N-n)\frac{\sin(qnl)}{qnl}\right) 
    20     \left( 3\frac{\sin(qR)-qR\cos(qR)}{(qr)^3}\right)^2\right] 
     18    P(Q) = \frac{scale}{V}\left[ m_{p}^2 
     19    \left(N+2\sum_{n-1}^{N-1}(N-n)\frac{sin(qnl)}{qnl}\right) 
     20    \left( 3\frac{sin(qR)-qRcos(qR)}{(qr)^3}\right)^2\right] 
    2121 
    2222where the mass $m_p$ is $(SLD_{pearl}-SLD_{solvent})*(volume\ of\ N\ pearls)$. 
    … …  
    5656    ["radius",      "Ang",       80.0, [0, inf],     "", "Radius of the pearls"], 
    5757    ["edge_sep",    "Ang",      350.0, [0, inf],     "", "Length of the string segment - surface to surface"], 
    58     ["num_pearls",  "",           3.0, [1, inf],     "", "Number of the pearls"], 
     58    ["num_pearls",  "",           3.0, [0, inf],     "", "Number of the pearls"], 
    5959    ["sld",   "1e-6/Ang^2", 1.0, [-inf, inf],  "sld", "SLD of the pearl spheres"], 
    6060    ["sld_solvent", "1e-6/Ang^2", 6.3, [-inf, inf],  "sld", "SLD of the solvent"], 

  • sasmodels/models/multilayer_vesicle.c

    r925ad6e r925ad6e  
    77          double sld_solvent, 
    88          double sld, 
    9           int n_pairs) 
     9          double n_pairs) 
    1010{ 
    1111    //calculate with a loop, two shells at a time 
    … …  
    4747          double sld_solvent, 
    4848          double sld, 
    49           double fp_n_pairs) 
     49          double n_pairs) 
    5050{ 
    51     int n_pairs = (int)(fp_n_pairs + 0.5); 
    5251    return multilayer_vesicle_kernel(q, 
    5352           volfraction, 

  • sasmodels/models/multilayer_vesicle.py

    r925ad6e r925ad6e  
    1919 
    2020.. math:: 
    21     P(q) = \text{scale} \cdot \frac{V_f}{V_t} F^2(q) + \text{background} 
    2221 
    23 for 
     22    P(q) = \frac{\text{scale.volfraction}}{V_t} F^2(q) + \text{background} 
     23 
     24where 
    2425 
    2526.. math:: 
    26     F(q) = (\rho_\text{shell}-\rho_\text{solv}) \sum_{i=1}^{n_\text{pairs}} 
    27         \left[ 
    28           3V(R_i)\frac{\sin(qR_i)-qR_i\cos(qR_i)}{(qR_i)^3} \\ 
    29           - 3V(R_i+t_s)\frac{\sin(q(R_i+t_s))-q(R_i+t_s)\cos(q(R_i+t_s))}{(q(R_i+t_s))^3} 
    30         \right] 
    3127 
    32 and 
     28     F(q) = (\rho_{shell}-\rho_{solv}) \sum_{i=1}^{n\_pairs} \left[ 
     29     3V(R_i)\frac{\sin(qR_i)-qR_i\cos(qR_i)}{(qR_i)^3} \\ 
     30      - 3V(R_i+t_s)\frac{\sin(q(R_i+t_s))-q(R_i+t_s)\cos(q(R_i+t_s))}{(q(R_i+t_s))^3} 
     31     \right] 
    3332 
    34 .. math:: 
    35      R_i = r_c + (i-1)(t_s + t_w) 
    3633 
    37 where $V_f$ is the volume fraction of particles, $V_t$ is the volume of the 
    38 whole particle, $V(r)$ is the volume of a sphere of radius $r$, $r_c$ is the 
    39 radius of the core, $\rho_\text{shell}$ is the scattering length density of a 
    40 shell, $\rho_\text{solv}$ is the scattering length density of the solvent. 
     34where $R_i = r_c + (i-1)(t_s + t_w)$ 
     35    
     36where $V_t$ is the volume of the whole particle, $V(R)$ is the volume of a sphere 
     37of radius $R$, $r_c$ is the radius of the core, $\rho_{shell}$ is the scattering length  
     38density of a shell, $\rho_{solv}$ is the scattering length density of the solvent. 
    4139 
    42 The outer most radius, $r_o = R_n + t_s$, is used for both the volume fraction 
    43 normalization and for the effective radius for *S(Q)* when $P(Q) * S(Q)$ 
    44 is applied. 
    4540 
    4641The 2D scattering intensity is the same as 1D, regardless of the orientation 
    … …  
    5045 
    5146    q = \sqrt{q_x^2 + q_y^2} 
     47 
     48 
     49The outer most radius 
     50 
     51$radius + n\_pairs * thick\_shell + (n\_pairs- 1) * thick\_solvent$ 
     52 
     53is used for both the volume fraction normalization and for the  
     54effective radius for *S(Q)* when $P(Q) * S(Q)$ is applied. 
    5255 
    5356For information about polarised and magnetic scattering, see 
    … …  
    6770**Author:** NIST IGOR/DANSE **on:** pre 2010 
    6871 
    69 **Last Modified by:** Piotr Rozyczko **on:** Feb 24, 2016 
     72**Last Modified by:** Piotr Rozyczko**on:** Feb 24, 2016 
    7073 
    7174**Last Reviewed by:** Paul Butler **on:** March 20, 2016 
    … …  
    106109source = ["lib/sas_3j1x_x.c", "multilayer_vesicle.c"] 
    107110 
    108 # TODO: the following line does nothing 
    109111polydispersity = ["radius", "n_pairs"] 
    110112 

  • sasmodels/models/onion.c

    r925ad6e r925ad6e  
    3030 
    3131static double 
    32 form_volume(double radius_core, double n_shells, double thickness[]) 
     32form_volume(double radius_core, double n, double thickness[]) 
    3333{ 
    34   int n = (int)(n_shells+0.5); 
     34  int i; 
    3535  double r = radius_core; 
    36   for (int i=0; i < n; i++) { 
     36  for (i=0; i < n; i++) { 
    3737    r += thickness[i]; 
    3838  } 

  • sasmodels/models/onion.py

    r925ad6e r925ad6e  
    323323    Returns shape profile with x=radius, y=SLD. 
    324324    """ 
    325     n_shells = int(n_shells+0.5) 
     325 
    326326    total_radius = 1.25*(sum(thickness[:n_shells]) + radius_core + 1) 
    327327    dz = total_radius/400  # 400 points for a smooth plot 
    … …  
    366366    return np.asarray(z), np.asarray(rho) 
    367367 
    368 def ER(radius_core, n_shells, thickness): 
     368def ER(radius_core, n, thickness): 
    369369    """Effective radius""" 
    370     n = int(n_shells[0]+0.5) 
    371     return np.sum(thickness[:n], axis=0) + radius_core 
     370    return np.sum(thickness[:int(n[0])], axis=0) + radius_core 
    372371 
    373372demo = { 

  • sasmodels/models/pearl_necklace.c

    r3f853beb r4b541ac  
    11double form_volume(double radius, double edge_sep, 
    2     double thick_string, double fp_num_pearls); 
     2    double thick_string, double num_pearls); 
    33double Iq(double q, double radius, double edge_sep, 
    4     double thick_string, double fp_num_pearls, double sld, 
     4    double thick_string, double num_pearls, double sld,  
    55    double string_sld, double solvent_sld); 
    66 
    … …  
    99// From Igor library 
    1010static double 
    11 pearl_necklace_kernel(double q, double radius, double edge_sep, double thick_string, 
    12     int num_pearls, double sld_pearl, double sld_string, double sld_solv) 
     11_pearl_necklace_kernel(double q, double radius, double edge_sep, double thick_string, 
     12    double num_pearls, double sld_pearl, double sld_string, double sld_solv) 
    1313{ 
    1414    // number of string segments 
    15     const int num_strings = num_pearls - 1; 
     15    num_pearls = floor(num_pearls + 0.5); //Force integer number of pearls 
     16    const double num_strings = num_pearls - 1.0; 
    1617 
    1718    //each masses: contrast * volume 
    … …  
    6869 
    6970double form_volume(double radius, double edge_sep, 
    70     double thick_string, double fp_num_pearls) 
     71    double thick_string, double num_pearls) 
    7172{ 
    72     const int num_pearls = (int)(fp_num_pearls + 0.5); //Force integer number of pearls 
    73     const int num_strings = num_pearls - 1; 
     73    num_pearls = floor(num_pearls + 0.5); //Force integer number of pearls 
     74 
     75    const double num_strings = num_pearls - 1.0; 
    7476    const double string_vol = edge_sep * M_PI_4 * thick_string * thick_string; 
    7577    const double pearl_vol = M_4PI_3 * radius * radius * radius; 
    … …  
    8082 
    8183double Iq(double q, double radius, double edge_sep, 
    82     double thick_string, double fp_num_pearls, double sld, 
     84    double thick_string, double num_pearls, double sld,  
    8385    double string_sld, double solvent_sld) 
    8486{ 
    85     const int num_pearls = (int)(fp_num_pearls + 0.5); //Force integer number of pearls 
    86     const double form = pearl_necklace_kernel(q, radius, edge_sep, 
     87    const double form = _pearl_necklace_kernel(q, radius, edge_sep, 
    8788        thick_string, num_pearls, sld, string_sld, solvent_sld); 
    8889 

  • sasmodels/models/pearl_necklace.py

    r1bd1ea2 r4b541ac  
    8282              ["thick_string", "Ang", 2.5, [0, inf], "volume", 
    8383               "Thickness of the chain linkage"], 
    84               ["num_pearls", "none", 3, [1, inf], "volume", 
     84              ["num_pearls", "none", 3, [0, inf], "volume", 
    8585               "Number of pearls in the necklace (must be integer)"], 
    8686              ["sld", "1e-6/Ang^2", 1.0, [-inf, inf], "sld", 
    … …  
    100100    Redundant with form_volume. 
    101101    """ 
    102     num_pearls = int(num_pearls + 0.5) 
    103102    number_of_strings = num_pearls - 1.0 
    104103    string_vol = edge_sep * pi * pow((thick_string / 2.0), 2.0) 
    … …  
    112111    Calculation for effective radius. 
    113112    """ 
    114     num_pearls = int(num_pearls + 0.5) 
    115113    tot_vol = volume(radius, edge_sep, thick_string, num_pearls) 
    116114    rad_out = (tot_vol/(4.0/3.0*pi)) ** (1./3.) 

  • sasmodels/models/raspberry.py

    r8e68ea0 r8e68ea0  
    1010    Schematic of the raspberry model 
    1111 
    12 In order to calculate the form factor of the entire complex, the 
    13 self-correlation of the large droplet, the self-correlation of the particles, 
    14 the correlation terms between different particles and the cross terms between 
    15 large droplet and small particles all need to be calculated. 
     12In order to calculate the form factor of the entire complex, the self- 
     13correlation of the large droplet, the self-correlation of the particles, the 
     14correlation terms between different particles and the cross terms between large 
     15droplet and small particles all need to be calculated. 
    1616 
    17 Consider two infinitely thin shells of radii $R_1$ and $R_2$ separated by 
    18 distance $r$. The general structure of the equation is then the form factor 
    19 of the two shells multiplied by the phase factor that accounts for the 
    20 separation of their centers. 
     17Consider two infinitely thin shells of radii R1 and R2 separated by distance r. 
     18The general structure of the equation is then the form factor of the two shells 
     19multiplied by the phase factor that accounts for the separation of their 
     20centers. 
    2121 
    2222.. math:: 

  • sasmodels/models/rpa.c

    r20c856a racfb094  
    1 double Iq(double q, double fp_case_num, 
     1double Iq(double q, double case_num, 
    22    double N[], double Phi[], double v[], double L[], double b[], 
    33    double Kab, double Kac, double Kad, 
    … …  
    55    ); 
    66 
    7 double Iq(double q, double fp_case_num, 
     7double Iq(double q, double case_num, 
    88    double N[],    // DEGREE OF POLYMERIZATION 
    99    double Phi[],  // VOL FRACTION 
    … …  
    1515    ) 
    1616{ 
    17   int icase = (int)(fp_case_num+0.5); 
     17  int icase = (int)case_num; 
    1818 
    1919  double Nab,Nac,Nad,Nbc,Nbd,Ncd; 
    … …  
    309309 
    310310  //calculate contrast where L[i] is the scattering length of i and D is the matrix 
    311   //Note that should multiply by Nav to get units of SLD which will become 
    312   // Nav*2 in the next line (SLD^2) but then normalization in that line would 
    313   //need to divide by Nav leaving only Nav or sqrt(Nav) before squaring.  
     311  //need to verify why the sqrt of Nav rather than just Nav (assuming v is molar volume) 
    314312  Nav=6.022045e+23; 
    315313  Lad=(L[0]/v[0]-L[3]/v[3])*sqrt(Nav); 
    … …  
    319317  Intg=Lad*Lad*S11+Lbd*Lbd*S22+Lcd*Lcd*S33+2.0*Lad*Lbd*S12+2.0*Lbd*Lcd*S23+2.0*Lad*Lcd*S13; 
    320318 
    321   //rescale for units of Lij^2 (fm^2 to cm^2) 
    322   Intg *= 1.0e-26;     
     319  //rescale for units of Lij^2 (in 10e-12 m^2 to m^2 ?) 
     320  Intg *= 1.0e-24;     
    323321 
    324322  return Intg; 

  • sasmodels/models/rpa.py

    r20c856a rbb73096  
    11r""" 
    2 Definition 
    3 ---------- 
    4  
    52Calculates the macroscopic scattering intensity for a multi-component 
    63homogeneous mixture of polymers using the Random Phase Approximation. 
    … …  
    2724Case 9: A-B-C-D tetra-block copolymer 
    2825 
    29 .. note:: 
    30     These case numbers are different from those in the NIST SANS package! 
     26**NB: these case numbers are different from those in the NIST SANS package!** 
    3127 
    32 USAGE NOTES: 
     28Only one case can be used at any one time. 
    3329 
    34 * Only one case can be used at any one time. 
    35 * The RPA (mean field) formalism only applies only when the multicomponent 
    36   polymer mixture is in the homogeneous mixed-phase region. 
    37 * **Component D is assumed to be the "background" component (ie, all contrasts 
    38   are calculated with respect to component D).** So the scattering contrast 
    39   for a C/D blend = [SLD(component C) - SLD(component D)]\ :sup:`2`. 
    40 * Depending on which case is being used, the number of fitting parameters can  
    41   vary.  Note that in general the degrees of polymerization, the volume 
    42   fractions, the molar volumes, and the neutron scattering lengths for each 
    43   component are obtained from other methods and held fixed while the segment 
    44   lengths (b\ :sub:`a`, b\ :sub:`b`, etc) and $\chi$ parameters (K\ :sub:`ab`, 
    45   K\ :sub:`ac`, etc). The *scale* parameter should be held equal to unity. 
     30The RPA (mean field) formalism only applies only when the multicomponent 
     31polymer mixture is in the homogeneous mixed-phase region. 
     32 
     33**Component D is assumed to be the "background" component (ie, all contrasts 
     34are calculated with respect to component D).** So the scattering contrast 
     35for a C/D blend = [SLD(component C) - SLD(component D)]\ :sup:`2`. 
     36 
     37Depending on which case is being used, the number of fitting parameters - the 
     38segment lengths (ba, bb, etc) and $\chi$ parameters (Kab, Kac, etc) - vary. 
     39The *scale* parameter should be held equal to unity. 
     40 
     41The input parameters are the degrees of polymerization, the volume fractions, 
     42the specific volumes, and the neutron scattering length densities for each 
     43component. 
    4644 
    4745 
    … …  
    4947---------- 
    5048 
    51 .. [#] A Z Akcasu, R Klein and B Hammouda, *Macromolecules*, 26 (1993) 4136 
     49A Z Akcasu, R Klein and B Hammouda, *Macromolecules*, 26 (1993) 4136 
    5250""" 
    5351 
    … …  
    5553 
    5654name = "rpa" 
    57 title = "Random Phase Approximation" 
     55title = "Random Phase Approximation - unfinished work in progress" 
    5856description = """ 
    5957This formalism applies to multicomponent polymer mixtures in the 
    … …  
    9290    ["N[4]", "", 1000.0, [1, inf], "", "Degree of polymerization"], 
    9391    ["Phi[4]", "", 0.25, [0, 1], "", "volume fraction"], 
    94     ["v[4]", "mL/mol", 100.0, [0, inf], "", "molar volume"], 
     92    ["v[4]", "mL/mol", 100.0, [0, inf], "", "specific volume"], 
    9593    ["L[4]", "fm", 10.0, [-inf, inf], "", "scattering length"], 
    9694    ["b[4]", "Ang", 5.0, [0, inf], "", "segment length"], 
    … …  
    116114    Return a list of parameters to hide depending on the multiplicity parameter. 
    117115    """ 
    118     case_num = int(case_num+0.5) 
    119116    if case_num < 2: 
    120117        return HIDE_AB 

  • sasmodels/models/spherical_sld.c

    r925ad6e r925ad6e  
    11static double form_volume( 
    2     double fp_n_shells, 
     2    int n_shells, 
    33    double thickness[], 
    44    double interface[]) 
    55{ 
    6     int n_shells= (int)(fp_n_shells + 0.5); 
    76    double r = 0.0; 
    87    for (int i=0; i < n_shells; i++) { 
    … …  
    2120        return pow(z, nu); 
    2221    } else if (shape==2) { 
    23         return 1.0 - pow(1.0 - z, nu); 
     22        return 1.0 - pow(1. - z, nu); 
    2423    } else if (shape==3) { 
    2524        return expm1(-nu*z)/expm1(-nu); 
    … …  
    4544static double Iq( 
    4645    double q, 
    47     double fp_n_shells, 
     46    int n_shells, 
    4847    double sld_solvent, 
    4948    double sld[], 
    … …  
    5251    double shape[], 
    5352    double nu[], 
    54     double fp_n_steps) 
     53    int n_steps) 
    5554{ 
    5655    // iteration for # of shells + core + solvent 
    57     int n_shells = (int)(fp_n_shells + 0.5); 
    58     int n_steps = (int)(fp_n_steps + 0.5); 
    5956    double f=0.0; 
    6057    double r=0.0; 

  • sasmodels/models/spherical_sld.py

    r925ad6e r925ad6e  
    233233    """ 
    234234 
    235     n_shells = int(n_shells + 0.5) 
    236     n_steps = int(n_steps + 0.5) 
    237235    z = [] 
    238236    rho = [] 
    … …  
    242240    rho.append(sld[0]) 
    243241 
    244     for i in range(0, n_shells): 
     242    for i in range(0, int(n_shells)): 
    245243        z_next += thickness[i] 
    246244        z.append(z_next) 
    … …  
    263261def ER(n_shells, thickness, interface): 
    264262    """Effective radius""" 
    265     n_shells = int(n_shells + 0.5) 
     263    n_shells = int(n_shells) 
    266264    total = (np.sum(thickness[:n_shells], axis=1) 
    267265             + np.sum(interface[:n_shells], axis=1)) 

  • sasmodels/models/stacked_disks.c

    r19f996b r6c3e266  
    1 static double stacked_disks_kernel( 
    2     double q, 
    3     double halfheight, 
    4     double thick_layer, 
    5     double radius, 
    6     int n_stacking, 
    7     double sigma_dnn, 
    8     double core_sld, 
    9     double layer_sld, 
    10     double solvent_sld, 
    11     double sin_alpha, 
    12     double cos_alpha, 
    13     double d) 
     1double form_volume(double thick_core, 
     2                   double thick_layer, 
     3                   double radius, 
     4                   double n_stacking); 
     5 
     6double Iq(double q, 
     7          double thick_core, 
     8          double thick_layer, 
     9          double radius, 
     10          double n_stacking, 
     11          double sigma_dnn, 
     12          double core_sld, 
     13          double layer_sld, 
     14          double solvent_sld); 
     15 
     16double Iqxy(double qx, double qy, 
     17          double thick_core, 
     18          double thick_layer, 
     19          double radius, 
     20          double n_stacking, 
     21          double sigma_dnn, 
     22          double core_sld, 
     23          double layer_sld, 
     24          double solvent_sld, 
     25          double theta, 
     26          double phi); 
     27 
     28static 
     29double _kernel(double q, 
     30               double radius, 
     31               double core_sld, 
     32               double layer_sld, 
     33               double solvent_sld, 
     34               double halfheight, 
     35               double thick_layer, 
     36               double sin_alpha, 
     37               double cos_alpha, 
     38               double sigma_dnn, 
     39               double d, 
     40               double n_stacking) 
    1441 
    1542{ 
    … …  
    6188 
    6289 
    63 static double stacked_disks_1d( 
    64     double q, 
    65     double thick_core, 
    66     double thick_layer, 
    67     double radius, 
    68     int n_stacking, 
    69     double sigma_dnn, 
    70     double core_sld, 
    71     double layer_sld, 
    72     double solvent_sld) 
     90static 
     91double stacked_disks_kernel(double q, 
     92                            double thick_core, 
     93                            double thick_layer, 
     94                            double radius, 
     95                            double n_stacking, 
     96                            double sigma_dnn, 
     97                            double core_sld, 
     98                            double layer_sld, 
     99                            double solvent_sld) 
    73100{ 
    74101/*    StackedDiscsX  :  calculates the form factor of a stacked "tactoid" of core shell disks 
    … …  
    84111        double sin_alpha, cos_alpha; // slots to hold sincos function output 
    85112        SINCOS(zi, sin_alpha, cos_alpha); 
    86         double yyy = stacked_disks_kernel(q, 
    87                            halfheight, 
    88                            thick_layer, 
     113        double yyy = _kernel(q, 
    89114                           radius, 
    90                            n_stacking, 
    91                            sigma_dnn, 
    92115                           core_sld, 
    93116                           layer_sld, 
    94117                           solvent_sld, 
     118                           halfheight, 
     119                           thick_layer, 
    95120                           sin_alpha, 
    96121                           cos_alpha, 
    97                            d); 
     122                           sigma_dnn, 
     123                           d, 
     124                           n_stacking); 
    98125        summ += Gauss76Wt[i] * yyy * sin_alpha; 
    99126    } 
    … …  
    105132} 
    106133 
    107 static double form_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); 
     134double form_volume(double thick_core, 
     135                   double thick_layer, 
     136                   double radius, 
     137                   double n_stacking){ 
    114138    double d = 2.0 * thick_layer + thick_core; 
    115139    return M_PI * radius * radius * d * n_stacking; 
    116140} 
    117141 
    118 static double Iq( 
    119     double q, 
    120     double thick_core, 
    121     double thick_layer, 
    122     double radius, 
    123     double fp_n_stacking, 
    124     double sigma_dnn, 
    125     double core_sld, 
    126     double layer_sld, 
    127     double solvent_sld) 
    128 { 
    129     int n_stacking = (int)(fp_n_stacking + 0.5); 
    130     return stacked_disks_1d(q, 
     142double Iq(double q, 
     143          double thick_core, 
     144          double thick_layer, 
     145          double radius, 
     146          double n_stacking, 
     147          double sigma_dnn, 
     148          double core_sld, 
     149          double layer_sld, 
     150          double solvent_sld) 
     151{ 
     152    return stacked_disks_kernel(q, 
    131153                    thick_core, 
    132154                    thick_layer, 
    … …  
    140162 
    141163 
    142 static double Iqxy(double qx, double qy, 
    143     double thick_core, 
    144     double thick_layer, 
    145     double radius, 
    146     double fp_n_stacking, 
    147     double sigma_dnn, 
    148     double core_sld, 
    149     double layer_sld, 
    150     double solvent_sld, 
    151     double theta, 
    152     double phi) 
    153 { 
    154     int n_stacking = (int)(fp_n_stacking + 0.5); 
     164double 
     165Iqxy(double qx, double qy, 
     166     double thick_core, 
     167     double thick_layer, 
     168     double radius, 
     169     double n_stacking, 
     170     double sigma_dnn, 
     171     double core_sld, 
     172     double layer_sld, 
     173     double solvent_sld, 
     174     double theta, 
     175     double phi) 
     176{ 
    155177    double q, sin_alpha, cos_alpha; 
    156178    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); 
    … …  
    158180    double d = 2.0 * thick_layer + thick_core; 
    159181    double halfheight = 0.5*thick_core; 
    160     double answer = stacked_disks_kernel(q, 
    161                      halfheight, 
    162                      thick_layer, 
     182    double answer = _kernel(q, 
    163183                     radius, 
    164                      n_stacking, 
    165                      sigma_dnn, 
    166184                     core_sld, 
    167185                     layer_sld, 
    168186                     solvent_sld, 
     187                     halfheight, 
     188                     thick_layer, 
    169189                     sin_alpha, 
    170190                     cos_alpha, 
    171                      d); 
     191                     sigma_dnn, 
     192                     d, 
     193                     n_stacking); 
    172194 
    173195    //convert to [cm-1] 

  • sasmodels/models/stacked_disks.py

    rb7e8b94 rb7e8b94  
    126126    ["thick_layer", "Ang",        10.0, [0, inf],    "volume",      "Thickness of layer each side of core"], 
    127127    ["radius",      "Ang",        15.0, [0, inf],    "volume",      "Radius of the stacked disk"], 
    128     ["n_stacking",  "",            1.0, [1, inf],    "volume",      "Number of stacked layer/core/layer disks"], 
     128    ["n_stacking",  "",            1.0, [0, inf],    "volume",      "Number of stacked layer/core/layer disks"], 
    129129    ["sigma_d",     "Ang",         0,   [0, inf],    "",            "Sigma of nearest neighbor spacing"], 
    130130    ["sld_core",    "1e-6/Ang^2",  4,   [-inf, inf], "sld",         "Core scattering length density"], 

  • sasmodels/models/star_polymer.c

    r2586093f r3a48772  
    33double Iq(double q, double radius2, double arms); 
    44 
    5 static double star_polymer_kernel(double q, double radius2, double arms) 
     5static double _mass_fractal_kernel(double q, double radius2, double arms) 
    66{ 
    77 
    … …  
    2323double Iq(double q, double radius2, double arms) 
    2424{ 
    25     return star_polymer_kernel(q, radius2, arms); 
     25    return _mass_fractal_kernel(q, radius2, arms); 
    2626} 

  • sasmodels/models/unified_power_Rg.py

    r66ca2a6 r66ca2a6  
    9797 
    9898def Iq(q, level, rg, power, B, G): 
    99     level = int(level + 0.5) 
    100     if level == 0: 
     99    ilevel = int(level) 
     100    if ilevel == 0: 
    101101        with errstate(divide='ignore'): 
    102102            return 1./q 
    … …  
    104104    with errstate(divide='ignore', invalid='ignore'): 
    105105        result = np.zeros(q.shape, 'd') 
    106         for i in range(level): 
     106        for i in range(ilevel): 
    107107            exp_now = exp(-(q*rg[i])**2/3.) 
    108108            pow_now = (erf(q*rg[i]/sqrt(6.))**3/q)**power[i] 
    109             if i < level-1: 
     109            if i < ilevel-1: 
    110110                exp_next = exp(-(q*rg[i+1])**2/3.) 
    111111            else: 
    … …  
    113113            result += G[i]*exp_now + B[i]*exp_next*pow_now 
    114114 
    115     result[q == 0] = np.sum(G[:level]) 
     115    result[q == 0] = np.sum(G[:ilevel]) 
    116116    return result 
    117117 
Note: See TracChangeset for help on using the changeset viewer.