Changes in / [a0632d9:6abd46f] in sasmodels


Ignore:
Location:
sasmodels/models
Files:
19 added
8 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/HayterMSAsq.py

    r13ed84c rec2ca99  
    5252    [Hayter-Penfold MSA charged sphere interparticle S(Q) structure factor] 
    5353        Interparticle structure factor S(Q)for a charged hard spheres. 
    54         Routine takes absolute value of charge, use HardSphere if charge goes to zero. 
     54        Routine takes absolute value of charge, use HardSphere if charge 
     55        goes to zero. 
    5556        In sasview the effective radius will be calculated from the 
    5657        parameters used in P(Q). 
    5758""" 
    5859single = False  # double precision only for now 
     60 
     61# pylint: disable=bad-whitespace, line-too-long 
    5962#             [ "name", "units", default, [lower, upper], "type", "description" ], 
    60 parameters = [["effect_radius", "Ang", 20.75, [0, inf], "volume", 
    61                "effective radius of hard sphere"], 
    62               ["charge", "e", 19.0, [0, inf], "", 
    63                "charge on sphere (in electrons)"], 
    64               ["volfraction", "", 0.0192, [0, 0.74], "", 
    65                "volume fraction of spheres"], 
    66               ["temperature", "K", 318.16, [0, inf], "", 
    67                "temperature, in Kelvin, for Debye length calculation"], 
    68               ["saltconc", "M", 0.0, [-inf, inf], "", 
    69                "conc of salt, 1:1 electolyte, for Debye length"], 
    70               ["dielectconst", "", 71.08, [-inf, inf], "", 
    71                "dielectric constant of solvent (default water), for Debye length"], 
    72              ] 
     63parameters = [ 
     64    ["effect_radius", "Ang", 20.75,   [0, inf],    "volume", "effective radius of hard sphere"], 
     65    ["charge",        "e",   19.0,    [0, inf],    "", "charge on sphere (in electrons)"], 
     66    ["volfraction",   "",     0.0192, [0, 0.74],   "", "volume fraction of spheres"], 
     67    ["temperature",   "K",  318.16,   [0, inf],    "", "temperature, in Kelvin, for Debye length calculation"], 
     68    ["saltconc",      "M",    0.0,    [-inf, inf], "", "conc of salt, 1:1 electolyte, for Debye length"], 
     69    ["dielectconst",  "",    71.08,   [-inf, inf], "", "dielectric constant of solvent (default water), for Debye length"], 
     70    ] 
     71# pylint: enable=bad-whitespace, line-too-long 
    7372category = "structure-factor" 
    7473 
    … …  
    8887oldpars = dict() 
    8988# default parameter set,  use  compare.sh -midQ -linear 
    90 # note the calculation varies in different limiting cases so a wide range of parameters will be required for a thorough test! 
     89# note the calculation varies in different limiting cases so a wide range of 
     90# parameters will be required for a thorough test! 
    9191# odd that the default st has saltconc zero 
    92 demo = dict(effect_radius = 20.75,charge=19.0,volfraction = 0.0192,temperature=318.16,saltconc=0.05,dielectconst=71.08,effect_radius_pd = 0.1,effect_radius_pd_n = 40) 
     92demo = dict(effect_radius=20.75, 
     93            charge=19.0, 
     94            volfraction=0.0192, 
     95            temperature=318.16, 
     96            saltconc=0.05, 
     97            dielectconst=71.08, 
     98            effect_radius_pd=0.1, 
     99            effect_radius_pd_n=40) 
    93100# 
    94101# attempt to use same values as old sasview unit test 
    95102tests = [ 
    96         [ {'scale': 1.0, 'background' : 0.0, 'effect_radius' : 20.75, 'charge' : 19.0, 'volfraction' : 0.0192, 'temperature' : 298.0, 
    97           'saltconc' : 0,'dielectconst' : 78.0, 'effect_radius_pd' : 0}, [0.0010], [0.0712928]] 
    98         ] 
     103    [{'scale': 1.0, 
     104      'background': 0.0, 
     105      'effect_radius': 20.75, 
     106      'charge': 19.0, 
     107      'volfraction': 0.0192, 
     108      'temperature': 298.0, 
     109      'saltconc': 0, 
     110      'dielectconst': 78.0, 
     111      'effect_radius_pd': 0}, 
     112     [0.0010], [0.0712928]] 
     113    ] 
    99114 

  • sasmodels/models/bcc.c

    r7ed702f r9aac25d  
    4040} 
    4141 
    42 double _sphereform(double q, double radius, double sld, double solvent_sld){ 
    43     const double qr = q*radius; 
    44     double sn, cn; 
    45     SINCOS(qr, sn, cn); 
    46     const double bes = (qr == 0.0 ? 1.0 : 3.0*(sn-qr*cn)/(qr*qr*qr)); 
    47     const double fq = bes * (sld - solvent_sld)*form_volume(radius); 
    48     return 1.0e-4*fq*fq; 
    49 } 
    5042 
    5143double form_volume(double radius){ 
    … …  
    8779 
    8880        answer = (vb-va)/2.0*summ; 
    89         answer = answer*_sphereform(q,radius,sld,solvent_sld)*latticescale; 
     81        answer = answer*sphere_form(q,radius,sld,solvent_sld)*latticescale; 
    9082 
    9183    return answer; 
    … …  
    174166 
    175167  // Use SphereForm directly from libigor 
    176   answer = _sphereform(q,radius,sld,solvent_sld)*Zq*latticescale; 
     168  answer = sphere_form(q,radius,sld,solvent_sld)*Zq*latticescale; 
    177169 
    178170  return answer; 

  • sasmodels/models/bcc.py

    r13ed84c r9aac25d  
    132132# pylint: enable=bad-whitespace, line-too-long 
    133133 
    134 source = ["lib/J1.c", "lib/gauss150.c", "bcc.c"] 
     134source = ["lib/J1.c", "lib/gauss150.c", "lib/sphere_form.c", "bcc.c"] 
    135135 
    136136# parameters for demo 

  • sasmodels/models/fcc.c

    reeb8bac r9aac25d  
    77double _FCC_Integrand(double q, double dnn, double d_factor, double theta, double phi); 
    88double _FCCeval(double Theta, double Phi, double temp1, double temp3); 
    9 double _sphereform(double q, double radius, double sld, double solvent_sld); 
    109 
    1110 
    … …  
    3938 
    4039        return (result); 
    41 } 
    42  
    43 double _sphereform(double q, double radius, double sld, double solvent_sld){ 
    44     const double qr = q*radius; 
    45     double sn, cn; 
    46     SINCOS(qr, sn, cn); 
    47     const double bes = (qr == 0.0 ? 1.0 : 3.0*(sn-qr*cn)/(qr*qr*qr)); 
    48     const double fq = bes * (sld - solvent_sld)*form_volume(radius); 
    49     return 1.0e-4*fq*fq; 
    5040} 
    5141 
    … …  
    8878 
    8979        answer = (vb-va)/2.0*summ; 
    90         answer = answer*_sphereform(q,radius,sld,solvent_sld)*latticescale; 
     80        answer = answer*sphere_form(q,radius,sld,solvent_sld)*latticescale; 
    9181 
    9282    return answer; 
    … …  
    175165 
    176166  // Use SphereForm directly from libigor 
    177   answer = _sphereform(q,radius,sld,solvent_sld)*Zq*latticescale; 
     167  answer = sphere_form(q,radius,sld,solvent_sld)*Zq*latticescale; 
    178168 
    179169  return answer; 

  • sasmodels/models/fcc.py

    r13ed84c r9aac25d  
    125125             ] 
    126126 
    127 source = ["lib/J1.c", "lib/gauss150.c", "fcc.c"] 
     127source = ["lib/J1.c", "lib/gauss150.c", "lib/sphere_form.c", "fcc.c"] 
    128128 
    129129# parameters for demo 

  • sasmodels/models/hardsphere.py

    r13ed84c r093f754  
    3535from numpy import inf 
    3636 
    37 name = "hardsphere" 
    38 title = "Hard sphere structure factor, with Percus-Yevick closure" 
     37name = "hardsphere_fish" 
     38title = "Hard sphere structure factor from FISH, with Percus-Yevick closure" 
    3939description = """\ 
    4040    [Hard sphere structure factor, with Percus-Yevick closure] 
    … …  
    5555               "volume fraction of hard spheres"], 
    5656             ] 
    57 single = False 
    5857 
    5958# No volume normalization despite having a volume parameter 
    … …  
    6463 
    6564Iq = """ 
    66     double denom,dnum,alpha,beta,gamm,a,asq,ath,afor,rca,rsa; 
    67     double calp,cbeta,cgam,prefac,c,vstruc; 
    68     double struc; 
     65      double D,A,B,G,X,X2,X4,S,C,FF,HARDSPH; 
    6966 
    70     //  compute constants 
    71     denom = pow((1.0-volfraction),4); 
    72     dnum = pow((1.0 + 2.0*volfraction),2); 
    73     alpha = dnum/denom; 
    74     beta = -6.0*volfraction*pow((1.0 + volfraction/2.0),2)/denom; 
    75     gamm = 0.50*volfraction*dnum/denom; 
    76     // 
    77     //  calculate the structure factor 
    78     // 
    79     a = 2.0*q*effect_radius; 
    80     asq = a*a; 
    81     ath = asq*a; 
    82     afor = ath*a; 
    83     SINCOS(a,rsa,rca); 
    84     //rca = cos(a); 
    85     //rsa = sin(a); 
    86     calp = alpha*(rsa/asq - rca/a); 
    87     cbeta = beta*(2.0*rsa/asq - (asq - 2.0)*rca/ath - 2.0/ath); 
    88     cgam = gamm*(-rca/a + (4.0/a)*((3.0*asq - 6.0)*rca/afor + (asq - 6.0)*rsa/ath + 6.0/afor)); 
    89     prefac = -24.0*volfraction/a; 
    90     c = prefac*(calp + cbeta + cgam); 
    91     vstruc = 1.0/(1.0-c); 
    92     struc = vstruc; 
     67      if(fabs(effect_radius) < 1.E-12) { 
     68               HARDSPH=1.0; 
     69               return(HARDSPH); 
     70      } 
     71      D=pow((1.-volfraction),2); 
     72      A=pow((1.+2.*volfraction)/D, 2); 
     73      X=fabs(q*effect_radius*2.0); 
    9374 
    94     return(struc); 
     75      if(X < 5.E-06) { 
     76                 HARDSPH=1./A; 
     77                 return(HARDSPH); 
     78      } 
     79      X2=pow(X,2); 
     80      X4=pow(X2,2); 
     81      B=-6.*volfraction* pow((1.+0.5*volfraction)/D ,2); 
     82      G=0.5*volfraction*A; 
     83 
     84      if(X < 0.2) { 
     85      // use Taylor series expansion for small X, IT IS VERY PICKY ABOUT THE X CUT OFF VALUE, ought to be lower in double.  
     86      // No obvious way to rearrange the equations to avoid needing a very high number of significant figures.  
     87      // Series expansion found using Mathematica software. Numerical test in .xls showed terms to X^2 are sufficient  
     88      // for 5 or 6 significant figures but I put the X^4 one in anyway  
     89            FF = 8*A +6*B + 4*G - (0.8*A +2.0*B/3.0 +0.5*G)*X2 +(A/35. +B/40. +G/50.)*X4; 
     90            // combining the terms makes things worse at smallest Q in single precision 
     91            //FF = (8-0.8*X2)*A +(3.0-X2/3.)*2*B + (4+0.5*X2)*G +(A/35. +B/40. +G/50.)*X4; 
     92            // note that G = -volfraction*A/2, combining this makes no further difference at smallest Q 
     93            //FF = (8 +2.*volfraction + ( volfraction/4. -0.8 +(volfraction/100. -1./35.)*X2 )*X2 )*A  + (3.0 -X2/3. +X4/40)*2*B; 
     94            HARDSPH= 1./(1. + volfraction*FF ); 
     95            return(HARDSPH); 
     96      } 
     97      SINCOS(X,S,C); 
     98 
     99// RKH Feb 2016, use version from FISH code as it is better than original sasview one at small Q in single precision 
     100      FF=A*(S-X*C)/X + B*(2.*X*S -(X2-2.)*C -2.)/X2 + G*( (4.*X2*X -24.*X)*S -(X4 -12.*X2 +24.)*C +24. )/X4; 
     101      HARDSPH= 1./(1. + 24.*volfraction*FF/X2 ); 
     102 
     103// rearrange the terms, is now about same as sasmodels 
     104//     FF=A*(S/X-C) + B*(2.*S/X - C +2.0*(C-1.0)/X2) + G*( (4./X -24./X3)*S -(1.0 -12./X2 +24./X4)*C +24./X4 ); 
     105//     HARDSPH= 1./(1. + 24.*volfraction*FF/X2 ); 
     106// remove 1/X2 from final line, take more powers of X inside the brackets, stil bad 
     107//      FF=A*(S/X3-C/X2) + B*(2.*S/X3 - C/X2 +2.0*(C-1.0)/X4) + G*( (4./X -24./X3)*S -(1.0 -12./X2 +24./X4)*C +24./X4 )/X2; 
     108//      HARDSPH= 1./(1. + 24.*volfraction*FF ); 
     109      return(HARDSPH); 
    95110   """ 
    96111 
    … …  
    106121oldname = 'HardsphereStructure' 
    107122oldpars = dict() 
    108  
     123# Q=0.001 is in the Taylor series, low Q part, so add Q=0.1, assuming double precision sasview is correct 
    109124tests = [ 
    110125        [ {'scale': 1.0, 'background' : 0.0, 'effect_radius' : 50.0, 'volfraction' : 0.2, 
    111            'effect_radius_pd' : 0}, [0.001], [0.209128]] 
     126           'effect_radius_pd' : 0}, [0.001,0.1], [0.209128,0.930587]] 
    112127        ] 
    113128 

  • sasmodels/models/hollow_cylinder.py

    r0420af7 rec2ca99  
    7878""" 
    7979category = "shape:cylinder" 
    80  
    81 #             ["name", "units", default, [lower, upper], "type","description"], 
     80# pylint: disable=bad-whitespace, line-too-long 
     81#   ["name", "units", default, [lower, upper], "type","description"], 
    8282parameters = [ 
    83               ["radius", "Ang", 30.0, [0, inf], "volume", "Cylinder radius"], 
    84               ["core_radius", "Ang", 20.0, [0, inf], "volume", "Hollow core radius"], 
    85               ["length", "Ang", 400.0, [0, inf], "volume", "Cylinder length"], 
    86               ["sld", "1/Ang^2", 6.3, [-inf, inf], "", "Cylinder sld"], 
    87               ["solvent_sld", "1/Ang^2", 1, [-inf, inf], "", "Solvent sld"], 
    88               ["theta", "degrees", 90, [-360, 360], "orientation", "Theta angle"], 
    89               ["phi", "degrees", 0, [-360, 360], "orientation", "Phi angle"], 
    90               ] 
     83    ["radius",      "Ang",     30.0, [0, inf],    "volume",      "Cylinder radius"], 
     84    ["core_radius", "Ang",     20.0, [0, inf],    "volume",      "Hollow core radius"], 
     85    ["length",      "Ang",    400.0, [0, inf],    "volume",      "Cylinder length"], 
     86    ["sld",         "1/Ang^2",  6.3, [-inf, inf], "",            "Cylinder sld"], 
     87    ["solvent_sld", "1/Ang^2",  1,   [-inf, inf], "",            "Solvent sld"], 
     88    ["theta",       "degrees", 90,   [-360, 360], "orientation", "Theta angle"], 
     89    ["phi",         "degrees",  0,   [-360, 360], "orientation", "Phi angle"], 
     90    ] 
     91# pylint: enable=bad-whitespace, line-too-long 
    9192 
    9293source = ["lib/J1.c", "lib/gauss76.c", "hollow_cylinder.c"] 
    9394 
    9495def ER(radius, core_radius, length): 
     96    """ 
     97    :param radius:      Cylinder radius 
     98    :param core_radius: Hollow core radius, UNUSED 
     99    :param length:      Cylinder length 
     100    :return:            Effective radius 
     101    """ 
    95102    if radius == 0 or length == 0: 
    96103        return 0.0 
    … …  
    104111 
    105112def VR(radius, core_radius, length): 
     113    """ 
     114    :param radius:      Cylinder radius 
     115    :param core_radius: Hollow core radius 
     116    :param length:      Cylinder length 
     117    :return:            Volf ratio for P(q)*S(q) 
     118    """ 
    106119    vol_core = pi*core_radius*core_radius*length 
    107120    vol_total = pi*radius*radius*length 
    … …  
    110123 
    111124# parameters for demo 
    112 demo = dict(scale=1.0,background=0.0,length=400.0,radius=30.0,core_radius=20.0, 
    113             sld=6.3,solvent_sld=1,theta=90,phi=0, 
     125demo = dict(scale=1.0, background=0.0, length=400.0, radius=30.0, 
     126            core_radius=20.0, sld=6.3, solvent_sld=1, theta=90, phi=0, 
    114127            radius_pd=.2, radius_pd_n=9, 
    115128            length_pd=.2, length_pd_n=10, 
    116129            core_radius_pd=.2, core_radius_pd_n=9, 
    117130            theta_pd=10, theta_pd_n=5, 
    118             ) 
     131           ) 
    119132 
    120133# For testing against the old sasview models, include the converted parameter 
    121134# names and the target sasview model name. 
    122135oldname = 'HollowCylinderModel' 
    123 oldpars = dict(scale='scale',background='background',radius='radius', 
    124                core_radius='core_radius',sld='sldCyl',length='length', 
    125                solvent_sld='sldSolv',phi='axis_phi',theta='axis_theta') 
     136oldpars = dict(scale='scale', background='background', radius='radius', 
     137               core_radius='core_radius', sld='sldCyl', length='length', 
     138               solvent_sld='sldSolv', phi='axis_phi', theta='axis_theta') 
    126139 
    127140# Parameters for unit tests 
    128141tests = [ 
    129          [{"radius" : 30.0},0.00005,1764.926], 
    130          [{},'VR',1.8], 
    131          [{},0.001,1756.76] 
    132          ] 
     142    [{"radius": 30.0}, 0.00005, 1764.926], 
     143    [{}, 'VR', 1.8], 
     144    [{}, 0.001, 1756.76] 
     145    ] 

  • sasmodels/models/lamellarCaille.py

    r13ed84c rec2ca99  
    1 # Note: model title and parameter table are inserted automatically 
    21r""" 
    32This model provides the scattering intensity, $I(q) = P(q) S(q)$, for a 
    … …  
    7776title = "Random lamellar sheet with Caille structure factor" 
    7877description = """\ 
    79         [Random lamellar phase with Caille  structure factor] 
    80         randomly oriented stacks of infinite sheets 
    81                 with Caille S(Q), having polydisperse spacing. 
    82         sld = sheet scattering length density 
    83                 sld_solvent = solvent scattering length density 
    84                 background = incoherent background 
    85                 scale = scale factor 
     78    [Random lamellar phase with Caille  structure factor] 
     79    randomly oriented stacks of infinite sheets 
     80    with Caille S(Q), having polydisperse spacing. 
     81    sld = sheet scattering length density 
     82    sld_solvent = solvent scattering length density 
     83    background = incoherent background 
     84    scale = scale factor 
    8685""" 
    8786category = "shape:lamellae" 
    8887 
    8988single = False 
    90  
     89# pylint: disable=bad-whitespace, line-too-long 
    9190#             ["name", "units", default, [lower, upper], "type","description"], 
    92 parameters = [["thickness", "Ang",  30.0, [0, inf], "volume", "sheet thickness"], 
    93               ["Nlayers", "",  20, [0, inf], "", "Number of layers"], 
    94               ["spacing", "Ang", 400., [0.0,inf], "volume", "d-spacing of Caille S(Q)"], 
    95               ["Caille_parameter", "1/Ang^2", 0.1, [0.0,0.8], "", "Caille parameter"], 
    96               ["sld", "1e-6/Ang^2", 6.3, [-inf,inf], "", 
    97                "layer scattering length density"], 
    98               ["solvent_sld", "1e-6/Ang^2", 1.0, [-inf,inf], "", 
    99                "Solvent scattering length density"], 
    100              ] 
     91parameters = [ 
     92    ["thickness",        "Ang",      30.0,  [0, inf],   "volume", "sheet thickness"], 
     93    ["Nlayers",          "",          20,   [0, inf],   "",       "Number of layers"], 
     94    ["spacing",          "Ang",      400.,  [0.0,inf],  "volume", "d-spacing of Caille S(Q)"], 
     95    ["Caille_parameter", "1/Ang^2",    0.1, [0.0,0.8],  "",       "Caille parameter"], 
     96    ["sld",              "1e-6/Ang^2", 6.3, [-inf,inf], "",       "layer scattering length density"], 
     97    ["solvent_sld",      "1e-6/Ang^2", 1.0, [-inf,inf], "",       "Solvent scattering length density"], 
     98    ] 
     99# pylint: enable=bad-whitespace, line-too-long 
    101100 
    102101source = ["lamellarCaille_kernel.c"] 
    … …  
    116115 
    117116demo = dict(scale=1, background=0, 
    118             thickness=67.,Nlayers=3.75,spacing=200., 
    119             Caille_parameter=0.268,sld=1.0, solvent_sld=6.34, 
    120             thickness_pd= 0.1, thickness_pd_n=100, 
    121             spacing_pd= 0.05, spacing_pd_n=40) 
     117            thickness=67., Nlayers=3.75, spacing=200., 
     118            Caille_parameter=0.268, sld=1.0, solvent_sld=6.34, 
     119            thickness_pd=0.1, thickness_pd_n=100, 
     120            spacing_pd=0.05, spacing_pd_n=40) 
    122121 
    123122oldname = 'LamellarPSModel' 
    124 oldpars = dict(thickness='delta', Nlayers='N_plates', Caille_parameter='caille', 
    125                sld='sld_bi',solvent_sld='sld_sol') 
     123oldpars = dict(thickness='delta', Nlayers='N_plates', 
     124               Caille_parameter='caille', 
     125               sld='sld_bi', solvent_sld='sld_sol') 
    126126# 
    127127tests = [ 
    128         [ {'scale': 1.0, 'background' : 0.0, 'thickness' : 30.,'Nlayers' : 20.0, 'spacing' : 400., 
    129             'Caille_parameter' : 0.1, 'sld' : 6.3, 'solvent_sld' : 1.0, 
    130             'thickness_pd' : 0.0, 'spacing_pd' : 0.0 }, [0.001], [28895.13397]] 
    131         ] 
     128    [{'scale': 1.0, 'background': 0.0, 'thickness': 30., 'Nlayers': 20.0, 
     129      'spacing': 400., 'Caille_parameter': 0.1, 'sld': 6.3, 
     130      'solvent_sld': 1.0, 'thickness_pd': 0.0, 'spacing_pd': 0.0}, 
     131     [0.001], [28895.13397]] 
     132    ] 
Note: See TracChangeset for help on using the changeset viewer.