Changes in / [26a6608:a70959a] in sasmodels


Ignore:
Location:
sasmodels
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/core_shell_parallelepiped.c

    rc69d6d6 r92dfe0c  
    1 double form_volume(double length_a, double length_b, double length_c, 
     1double form_volume(double length_a, double length_b, double length_c,  
    22                   double thick_rim_a, double thick_rim_b, double thick_rim_c); 
    33double Iq(double q, double core_sld, double arim_sld, double brim_sld, double crim_sld, 
    … …  
    99            double thick_rim_c, double theta, double phi, double psi); 
    1010 
    11 double form_volume(double length_a, double length_b, double length_c, 
     11double form_volume(double length_a, double length_b, double length_c,  
    1212                   double thick_rim_a, double thick_rim_b, double thick_rim_c) 
    1313{ 
    1414    //return length_a * length_b * length_c; 
    15     return length_a * length_b * length_c + 
    16            2.0 * thick_rim_a * length_b * length_c + 
     15    return length_a * length_b * length_c +  
     16           2.0 * thick_rim_a * length_b * length_c +  
    1717           2.0 * thick_rim_b * length_a * length_c + 
    1818           2.0 * thick_rim_c * length_a * length_b; 
    … …  
    3434    // Code converted from functions CSPPKernel and CSParallelepiped in libCylinder.c_scaled 
    3535    // Did not understand the code completely, it should be rechecked (Miguel Gonzalez) 
    36     //Code is rewritten,the code is compliant with Diva Singhs thesis now (Dirk Honecker) 
    37  
     36     
    3837    const double mu = 0.5 * q * length_b; 
    39  
     38     
    4039    //calculate volume before rescaling (in original code, but not used) 
    41     //double vol = form_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c); 
    42     //double vol = length_a * length_b * length_c + 
    43     //       2.0 * thick_rim_a * length_b * length_c + 
     40    //double vol = form_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c);         
     41    //double vol = length_a * length_b * length_c +  
     42    //       2.0 * thick_rim_a * length_b * length_c +  
    4443    //       2.0 * thick_rim_b * length_a * length_c + 
    4544    //       2.0 * thick_rim_c * length_a * length_b; 
    46  
     45     
    4746    // Scale sides by B 
    4847    const double a_scaled = length_a / length_b; 
    4948    const double c_scaled = length_c / length_b; 
    5049 
    51     double ta = a_scaled + 2.0*thick_rim_a/length_b; // incorrect ta = (a_scaled + 2.0*thick_rim_a)/length_b; 
    52     double tb = 1+ 2.0*thick_rim_b/length_b; // incorrect tb = (a_scaled + 2.0*thick_rim_b)/length_b; 
    53     double tc = c_scaled + 2.0*thick_rim_c/length_b; //not present 
     50    // ta and tb correspond to the definitions in CSPPKernel, but they don't make sense to me (MG) 
     51    // the a_scaled in the definition of tb was present in CSPPKernel in libCylinder.c, 
     52    // while in cspkernel in csparallelepiped.cpp (used for the 2D), all the definitions 
     53    // for ta, tb, tc use also A + 2*rim_thickness (but not scaled by B!!!) 
     54    double ta = (a_scaled + 2.0*thick_rim_a)/length_b; 
     55    double tb = (a_scaled + 2.0*thick_rim_b)/length_b; 
    5456 
    5557    double Vin = length_a * length_b * length_c; 
    … …  
    6062    double V1 = (2.0 * thick_rim_a * length_b * length_c);    // incorrect V1 (aa*bb*cc+2*ta*bb*cc) 
    6163    double V2 = (2.0 * length_a * thick_rim_b * length_c);    // incorrect V2(aa*bb*cc+2*aa*tb*cc) 
    62     double V3 = (2.0 * length_a * length_b * thick_rim_c);    //not present 
    63     double Vot = Vin + V1 + V2 + V3; 
    6464 
    6565    // Scale factors (note that drC is not used later) 
    … …  
    6767    const double drhoA = (arim_sld-solvent_sld); 
    6868    const double drhoB = (brim_sld-solvent_sld); 
    69     const double drhoC = (crim_sld-solvent_sld);  // incorrect const double drC_Vot = (crim_sld-solvent_sld)*Vot; 
    70  
     69    //const double drC_Vot = (crim_sld-solvent_sld)*Vot; 
    7170 
    7271    // Precompute scale factors for combining cross terms from the shape 
    7372    const double scale23 = drhoA*V1; 
    7473    const double scale14 = drhoB*V2; 
    75     const double scale24 = drhoC*V3; 
    76     const double scale11 = drho0*Vin; 
    77     const double scale12 = drho0*Vin - scale23 - scale14 - scale24; 
     74    const double scale12 = drho0*Vin - scale23 - scale14; 
    7875 
    7976    // outer integral (with gauss points), integration limits = 0, 1 
    … …  
    8683        // inner integral (with gauss points), integration limits = 0, 1 
    8784        double inner_total = 0.0; 
    88         double inner_total_crim = 0.0; 
    8985        for(int j=0; j<76; j++) { 
    9086            const double uu = 0.5 * ( Gauss76Z[j] + 1.0 ); 
    … …  
    9288            SINCOS(M_PI_2*uu, sin_uu, cos_uu); 
    9389            const double si1 = sas_sinx_x(mu_proj * sin_uu * a_scaled); 
    94             const double si2 = sas_sinx_x(mu_proj * cos_uu ); 
     90            const double si2 = sas_sinx_x(mu_proj * cos_uu); 
    9591            const double si3 = sas_sinx_x(mu_proj * sin_uu * ta); 
    9692            const double si4 = sas_sinx_x(mu_proj * cos_uu * tb); 
    … …  
    9894            // Expression in libCylinder.c (neither drC nor Vot are used) 
    9995            const double form = scale12*si1*si2 + scale23*si2*si3 + scale14*si1*si4; 
    100             const double form_crim = scale11*si1*si2; 
    10196 
    102  
     97            // To note also that in csparallelepiped.cpp, there is a function called 
     98            // cspkernel, which seems to make more sense and has the following comment: 
     99            //   This expression is different from NIST/IGOR package (I strongly believe the IGOR is wrong!!!). 10/15/2010. 
     100            //   tmp =( dr0*tmp1*tmp2*tmp3*Vin + drA*(tmpt1-tmp1)*tmp2*tmp3*V1+ drB*tmp1*(tmpt2-tmp2)*tmp3*V2 + drC*tmp1*tmp2*(tmpt3-tmp3)*V3)* 
     101            //   ( dr0*tmp1*tmp2*tmp3*Vin + drA*(tmpt1-tmp1)*tmp2*tmp3*V1+ drB*tmp1*(tmpt2-tmp2)*tmp3*V2 + drC*tmp1*tmp2*(tmpt3-tmp3)*V3);   //  correct FF : square of sum of phase factors 
     102            // This is the function called by csparallelepiped_analytical_2D_scaled, 
     103            // while CSParallelepipedModel calls CSParallelepiped in libCylinder.c         
     104             
    103105            //  correct FF : sum of square of phase factors 
    104106            inner_total += Gauss76Wt[j] * form * form; 
    105             inner_total_crim += Gauss76Wt[j] * form_crim * form_crim; 
    106107        } 
    107108        inner_total *= 0.5; 
    108         inner_total_crim *= 0.5; 
     109 
    109110        // now sum up the outer integral 
    110111        const double si = sas_sinx_x(mu * c_scaled * sigma); 
    111         const double si_crim = sas_sinx_x(mu * tc * sigma); 
    112         outer_total += Gauss76Wt[i] * (inner_total * si * si + inner_total_crim * si_crim * si_crim); 
     112        outer_total += Gauss76Wt[i] * inner_total * si * si; 
    113113    } 
    114114    outer_total *= 0.5; 
    … …  
    154154 
    155155    // The definitions of ta, tb, tc are not the same as in the 1D case because there is no 
    156     // the scaling by B. 
     156    // the scaling by B. The use of length_a for the 3 of them seems clearly a mistake to me, 
     157    // but for the moment I let it like this until understanding better the code. 
    157158    double ta = length_a + 2.0*thick_rim_a; 
    158     double tb = length_b + 2.0*thick_rim_b; 
    159     double tc = length_c + 2.0*thick_rim_c; 
     159    double tb = length_a + 2.0*thick_rim_b; 
     160    double tc = length_a + 2.0*thick_rim_c; 
    160161    //handle arg=0 separately, as sin(t)/t -> 1 as t->0 
    161162    double siA = sas_sinx_x(0.5*q*length_a*xhat); 
    … …  
    165166    double siBt = sas_sinx_x(0.5*q*tb*yhat); 
    166167    double siCt = sas_sinx_x(0.5*q*tc*zhat); 
    167  
     168     
    168169 
    169170    // f uses Vin, V1, V2, and V3 and it seems to have more sense than the value computed 
    … …  
    172173               + drA*(siAt-siA)*siB*siC*V1 
    173174               + drB*siA*(siBt-siB)*siC*V2 
    174                + drC*siA*siB*(siCt-siC)*V3); 
    175  
     175               + drC*siA*siB*(siCt*siCt-siC)*V3); 
     176    
    176177    return 1.0e-4 * f * f; 
    177178} 

  • sasmodels/models/core_shell_parallelepiped.py

    rfa70e04 r8f04da4  
    211211 
    212212# rkh 7/4/17 add random unit test for 2d, note make all params different, 2d values not tested against other codes or models 
    213 if 0:  # pak: model rewrite; need to update tests 
    214     qx, qy = 0.2 * cos(pi/6.), 0.2 * sin(pi/6.) 
    215     tests = [[{}, 0.2, 0.533149288477], 
    216             [{}, [0.2], [0.533149288477]], 
    217             [{'theta':10.0, 'phi':20.0}, (qx, qy), 0.0853299803222], 
    218             [{'theta':10.0, 'phi':20.0}, [(qx, qy)], [0.0853299803222]], 
    219             ] 
    220     del qx, qy  # not necessary to delete, but cleaner 
     213qx, qy = 0.2 * cos(pi/6.), 0.2 * sin(pi/6.) 
     214tests = [[{}, 0.2, 0.533149288477], 
     215         [{}, [0.2], [0.533149288477]], 
     216         [{'theta':10.0, 'phi':20.0}, (qx, qy), 0.0853299803222], 
     217         [{'theta':10.0, 'phi':20.0}, [(qx, qy)], [0.0853299803222]], 
     218        ] 
     219del qx, qy  # not necessary to delete, but cleaner 

  • sasmodels/product.py

    r146793b r058460c  
    100100    # Remember the component info blocks so we can build the model 
    101101    model_info.composition = ('product', [p_info, s_info]) 
    102     model_info.control = p_info.control 
    103     model_info.hidden = p_info.hidden 
    104     if getattr(p_info, 'profile', None) is not None: 
    105         profile_pars = set(p.id for p in p_info.parameters.kernel_parameters) 
    106         def profile(**kwargs): 
    107             # extract the profile args 
    108             kwargs = dict((k, v) for k, v in kwargs.items() if k in profile_pars) 
    109             return p_info.profile(**kwargs) 
    110     else: 
    111         profile = None 
    112     model_info.profile = profile 
    113     model_info.profile_axes = p_info.profile_axes 
    114  
    115102    # TODO: delegate random to p_info, s_info 
    116103    #model_info.random = lambda: {} 
    … …  
    142129    def __init__(self, model_info, P, S): 
    143130        # type: (ModelInfo, KernelModel, KernelModel) -> None 
    144         #: Combined info plock for the product model 
    145131        self.info = model_info 
    146         #: Form factor modelling individual particles. 
    147132        self.P = P 
    148         #: Structure factor modelling interaction between particles. 
    149133        self.S = S 
    150         #: Model precision. This is not really relevant, since it is the 
    151         #: individual P and S models that control the effective dtype, 
    152         #: converting the q-vectors to the correct type when the kernels 
    153         #: for each are created. Ideally this should be set to the more 
    154         #: precise type to avoid loss of precision, but precision in q is 
    155         #: not critical (single is good enough for our purposes), so it just 
    156         #: uses the precision of the form factor. 
    157         self.dtype = P.dtype  # type: np.dtype 
     134        self.dtype = P.dtype 
    158135 
    159136    def make_kernel(self, q_vectors): 

  • sasmodels/sasview_model.py

    redb0f85 r9f8ade1  
    205205                                           structure_factor._model_info) 
    206206    ConstructedModel = make_model_from_info(model_info) 
    207     return ConstructedModel(form_factor.multiplicity) 
     207    return ConstructedModel() 
    208208 
    209209 
    … …  
    323323    #: True if model has multiplicity 
    324324    is_multiplicity_model = False 
    325     #: Multiplicity information 
     325    #: Mulitplicity information 
    326326    multiplicity_info = None # type: MultiplicityInfoType 
    327327 
    … …  
    354354        # and lines to plot. 
    355355 
    356         # Get the list of hidden parameters given the multiplicity 
     356        # Get the list of hidden parameters given the mulitplicity 
    357357        # Don't include multiplicity in the list of parameters 
    358358        self.multiplicity = multiplicity 
Note: See TracChangeset for help on using the changeset viewer.