Changeset 2a0b2b1 in sasmodels for sasmodels/models/hollow_cylinder.c


Ignore:
Timestamp:
Apr 14, 2017 8:30:29 AM (7 years ago)
Author:
Paul Kienzle <pkienzle@…>
Branches:
master, core_shell_microgels, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
fdd56a1
Parents:
9901384
Message:

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

File:
1 edited

Legend:

Unmodified
Added
Removed
  • sasmodels/models/hollow_cylinder.c

    r592343f r2a0b2b1  
    11double form_volume(double radius, double thickness, double length); 
    22double Iq(double q, double radius, double thickness, double length, double sld, 
    3         double solvent_sld); 
     3    double solvent_sld); 
    44double Iqxy(double qx, double qy, double radius, double thickness, double length, double sld, 
    5         double solvent_sld, double theta, double phi); 
     5    double solvent_sld, double theta, double phi); 
    66 
    77//#define INVALID(v) (v.radius_core >= v.radius) 
     
    1414} 
    1515 
    16  
    1716static double 
    18 _hollow_cylinder_kernel(double q, 
    19     double radius, double thickness, double length, double sin_val, double cos_val) 
     17_fq(double qab, double qc, 
     18    double radius, double thickness, double length) 
    2019{ 
    21     const double qs = q*sin_val; 
    22     const double lam1 = sas_2J1x_x((radius+thickness)*qs); 
    23     const double lam2 = sas_2J1x_x(radius*qs); 
     20    const double lam1 = sas_2J1x_x((radius+thickness)*qab); 
     21    const double lam2 = sas_2J1x_x(radius*qab); 
    2422    const double gamma_sq = square(radius/(radius+thickness)); 
    25     //Note: lim_{thickness -> 0} psi = sas_J0(radius*qs) 
    26     //Note: lim_{radius -> 0} psi = sas_2J1x_x(thickness*qs) 
    27     const double psi = (lam1 - gamma_sq*lam2)/(1.0 - gamma_sq); //SRK 10/19/00 
    28     const double t2 = sas_sinx_x(0.5*q*length*cos_val); 
     23    //Note: lim_{thickness -> 0} psi = sas_J0(radius*qab) 
     24    //Note: lim_{radius -> 0} psi = sas_2J1x_x(thickness*qab) 
     25    const double psi = (lam1 - gamma_sq*lam2)/(1.0 - gamma_sq);    //SRK 10/19/00 
     26    const double t2 = sas_sinx_x(0.5*length*qc); 
    2927    return psi*t2; 
    3028} 
     
    4341{ 
    4442    const double lower = 0.0; 
    45     const double upper = 1.0;           //limits of numerical integral 
     43    const double upper = 1.0;        //limits of numerical integral 
    4644 
    47     double summ = 0.0;                  //initialize intergral 
     45    double summ = 0.0;            //initialize intergral 
    4846    for (int i=0;i<76;i++) { 
    49         const double cos_val = 0.5*( Gauss76Z[i] * (upper-lower) + lower + upper ); 
    50         const double sin_val = sqrt(1.0 - cos_val*cos_val); 
    51         const double inter = _hollow_cylinder_kernel(q, radius, thickness, length, 
    52                                                      sin_val, cos_val); 
    53         summ += Gauss76Wt[i] * inter * inter; 
     47        const double cos_theta = 0.5*( Gauss76Z[i] * (upper-lower) + lower + upper ); 
     48        const double sin_theta = sqrt(1.0 - cos_theta*cos_theta); 
     49        const double form = _fq(q*sin_theta, q*cos_theta, 
     50                                radius, thickness, length); 
     51        summ += Gauss76Wt[i] * form * form; 
    5452    } 
    5553 
     
    6664    double q, sin_alpha, cos_alpha; 
    6765    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); 
    68     const double Aq = _hollow_cylinder_kernel(q, radius, thickness, length, 
    69         sin_alpha, cos_alpha); 
     66    const double qab = q*sin_alpha; 
     67    const double qc = q*cos_alpha; 
     68 
     69    const double form = _fq(qab, qc, radius, thickness, length); 
    7070 
    7171    const double vol = form_volume(radius, thickness, length); 
    72     return _hollow_cylinder_scaling(Aq*Aq, solvent_sld-sld, vol); 
     72    return _hollow_cylinder_scaling(form*form, solvent_sld-sld, vol); 
    7373} 
    74  
Note: See TracChangeset for help on using the changeset viewer.