Changeset c1c37d5 in sasmodels


Ignore:
Timestamp:
Mar 18, 2016 12:06:08 PM (9 years ago)
Author:
krzywon
Branches:
master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
31c64d9, 3a45c2c
Parents:
89ab393 (diff), a629d8e (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
Message:

Merge remote-tracking branch 'origin/master'

Files:
4 added
9 edited
7 moved

Legend:

Unmodified
Added
Removed

  • sasmodels/core.py

    r02e70ff rd5ba841  
    170170    """ 
    171171    value, weight = zip(*pars) 
    172     if len(value) > 1: 
    173         value = [v.flatten() for v in np.meshgrid(*value)] 
    174         weight = np.vstack([v.flatten() for v in np.meshgrid(*weight)]) 
    175         weight = np.prod(weight, axis=0) 
     172    value = [v.flatten() for v in np.meshgrid(*value)] 
     173    weight = np.vstack([v.flatten() for v in np.meshgrid(*weight)]) 
     174    weight = np.prod(weight, axis=0) 
    176175    return value, weight 
    177176 

  • sasmodels/models/bessel.py

    rcbd37a7 ra5af4e1  
    7878 
    7979Iq = """ 
    80     return 2.0*j1(q)/q; 
     80    return J1(q); 
    8181    """ 
    8282 

  • sasmodels/models/hardsphere.py

    raa2edb2 r934f906  
    5858category = "structure-factor" 
    5959structure_factor = True 
     60single = False 
    6061 
    6162#             ["name", "units", default, [lower, upper], "type","description"], 
    … …  
    7475Iq = """ 
    7576      double D,A,B,G,X,X2,X4,S,C,FF,HARDSPH; 
     77      // these are c compiler instructions, can also put normal code inside the "if else" structure  
     78      #if FLOAT_SIZE > 4 
     79      // double precision    orig had 0.2, don't call the variable cutoff as PAK already has one called that! Must use UPPERCASE name please. 
     80      //  0.05 better, 0.1 OK 
     81      #define CUTOFFHS 0.05 
     82      #else 
     83      // 0.1 bad, 0.2 OK, 0.3 good, 0.4 better, 0.8 no good 
     84      #define CUTOFFHS 0.4   
     85      #endif 
    7686 
    7787      if(fabs(radius_effective) < 1.E-12) { 
    … …  
    96106      G=0.5*volfraction*A; 
    97107 
    98       if(X < 0.2) { 
    99       // RKH Feb 2016, use Taylor series expansion for small X, IT IS VERY PICKY ABOUT THE X CUT OFF VALUE, ought to be lower in double.  
     108      if(X < CUTOFFHS) { 
     109      // RKH Feb 2016, use Taylor series expansion for small X 
    100110      // else no obvious way to rearrange the equations to avoid needing a very high number of significant figures.  
    101111      // Series expansion found using Mathematica software. Numerical test in .xls showed terms to X^2 are sufficient  

  • sasmodels/models/lamellar.py

    raa2edb2 r7c391dd  
    55---------- 
    66 
    7 The scattering intensity $I(q)$ is 
     7The scattering intensity $I(q)$ for dilute, randomly oriented, "infinitely large" sheets or lamellae is 
    88 
    99.. math:: 
    1010 
    11     I(q) = \frac{2\pi P(q)}{\delta q^2} 
     11    I(q) = scale*\frac{2\pi P(q)}{q^2\delta } 
    1212 
    1313 
    … …  
    1616.. math:: 
    1717 
    18     P(q) = \frac{2\Delta\rho^2}{q^2}(1-cos(q\delta)) 
     18   P(q) = \frac{2\Delta\rho^2}{q^2}(1-cos(q\delta)) = \frac{4\Delta\rho^2}{q^2}sin^2(\frac{q\delta}{2}) 
    1919 
     20where $\delta$ is the total layer thickness and $\Delta\rho$ is the scattering length density difference. 
    2021 
    21 where $\delta$ is the bilayer thickness. 
     22This is the limiting form for a spherical shell of infinitely large radius. Note that the division by $\delta$ 
     23means that $scale$ in sasview is the volume fraction of sheet, $\phi = S\delta$ where $S$ is the area of  
     24sheet per unit volume. $S$ is half the Porod surface area per unit volume of a thicker layer (as that would  
     25include both faces of the sheet). 
    2226 
    2327The 2D scattering intensity is calculated in the same way as 1D, where 
    … …  
    5660 
    5761#             ["name", "units", default, [lower, upper], "type","description"], 
    58 parameters = [["sld", "1e-6/Ang^2", 1, [-inf, inf], "", 
    59                "Layer scattering length density" ], 
    60               ["solvent_sld", "1e-6/Ang^2", 6, [-inf, inf], "", 
    61                "Solvent scattering length density" ], 
    62               ["thickness", "Ang", 50, [0, inf], "volume","Bilayer thickness" ], 
     62parameters = [ ["thickness", "Ang", 50, [0, inf], "volume","total layer thickness" ], 
     63               ["sld", "1e-6/Ang^2", 1, [-inf, inf], "","Layer scattering length density" ], 
     64               ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "","Solvent scattering length density" ], 
    6365             ] 
    6466 
    65  
    6667# No volume normalization despite having a volume parameter 
    67 # This should perhaps be volume normalized? 
     68# This should perhaps be volume normalized? - it is! 
    6869form_volume = """ 
    6970    return 1.0; 
    … …  
    7172 
    7273Iq = """ 
    73     const double sub = sld - solvent_sld; 
     74    const double sub = sld - sld_solvent; 
    7475    const double qsq = q*q; 
    7576    // Original expression 
    … …  
    8990 
    9091demo = dict(scale=1, background=0, 
    91             sld=6, solvent_sld=1, 
     92            sld=6, sld_solvent=1, 
    9293            thickness=40, 
    9394            thickness_pd=0.2, thickness_pd_n=40) 
    9495oldname = 'LamellarModel' 
    95 oldpars = dict(sld='sld_bi', solvent_sld='sld_sol', thickness='bi_thick') 
     96oldpars = dict(sld='sld_bi', sld_solvent='sld_sol', thickness='bi_thick') 
    9697tests = [ 
    97         [ {'scale': 1.0, 'background' : 0.0, 'thickness' : 50.0, 'sld' : 1.0,'solvent_sld' : 6.3, 'thickness_pd' : 0.0,  
     98        [ {'scale': 1.0, 'background' : 0.0, 'thickness' : 50.0, 'sld' : 1.0,'sld_solvent' : 6.3, 'thickness_pd' : 0.0,  
    9899           }, [0.001], [882289.54309]] 
    99100        ] 

  • sasmodels/models/lamellar_hg.py

    raa2edb2 r52a3db3  
    2626 
    2727where $\delta_T$ is *tail_length*, $\delta_H$ is *head_length*, 
    28 $\Delta\rho_H$ is the head contrast (*head_sld* $-$ *solvent_sld*), 
    29 and $\Delta\rho_T$ is tail contrast (*sld* $-$ *solvent_sld*). 
     28$\Delta\rho_H$ is the head contrast (*sld_head* $-$ *sld_solvent*), 
     29and $\Delta\rho_T$ is tail contrast (*sld* $-$ *sld_solvent*). 
     30 
     31The total thickness of the lamellar sheet is $\delta_H$ + $\delta_T$ + $\delta_T$ + $\delta_H$. 
     32Note that in a non aqueous solvent the chemical "head" group may be the  
     33"Tail region" and vice-versa. 
    3034 
    3135The 2D scattering intensity is calculated in the same way as 1D, where 
    … …  
    4953from numpy import inf 
    5054 
    51 name = "lamellar_FFHG" 
    52 title = "Random lamellar phase with Head Groups" 
     55name = "lamellar_hg" 
     56title = "Random lamellar phase with Head and Tail Groups" 
    5357description = """\ 
    54     [Random lamellar phase with Head Groups] 
     58    [Random lamellar phase with Head and Tail Groups] 
    5559        I(q)= 2*pi*P(q)/(2(H+T)*q^(2)), where 
    5660        P(q)= see manual 
    5761        layer thickness =(H+T+T+H) = 2(Head+Tail) 
    5862        sld = Tail scattering length density 
    59         head_sld = Head scattering length density 
    60         solvent_sld = solvent scattering length density 
     63        sld_head = Head scattering length density 
     64        sld_solvent = solvent scattering length density 
    6165        background = incoherent background 
    6266        scale = scale factor 
    … …  
    6670# pylint: disable=bad-whitespace, line-too-long 
    6771#             ["name", "units", default, [lower, upper], "type","description"], 
    68 parameters = [["tail_length", "Ang",       15,   [0, inf],  "volume",  "Tail thickness"], 
    69               ["head_length", "Ang",       10,   [0, inf],  "volume",  "head thickness"], 
     72parameters = [["tail_length", "Ang",       15,   [0, inf],  "volume",  "Tail thickness ( total = H+T+T+H)"], 
     73              ["head_length", "Ang",       10,   [0, inf],  "volume",  "Head thickness"], 
    7074              ["sld",         "1e-6/Ang^2", 0.4, [-inf,inf], "",       "Tail scattering length density"], 
    71               ["head_sld",    "1e-6/Ang^2", 3.0, [-inf,inf], "",       "Head scattering length density"], 
    72               ["solvent_sld", "1e-6/Ang^2", 6,   [-inf,inf], "",       "Solvent scattering length density"]] 
     75              ["sld_head",    "1e-6/Ang^2", 3.0, [-inf,inf], "",       "Head scattering length density"], 
     76              ["sld_solvent", "1e-6/Ang^2", 6,   [-inf,inf], "",       "Solvent scattering length density"]] 
    7377# pylint: enable=bad-whitespace, line-too-long 
    7478 
    … …  
    8185Iq = """ 
    8286    const double qsq = q*q; 
    83     const double drh = head_sld - solvent_sld; 
    84     const double drt = sld - solvent_sld;    //correction 13FEB06 by L.Porcar 
     87    const double drh = sld_head - sld_solvent; 
     88    const double drt = sld - sld_solvent;    //correction 13FEB06 by L.Porcar 
    8589    const double qT = q*tail_length; 
    8690    double Pq, inten; 
    … …  
    106110demo = dict(scale=1, background=0, 
    107111            tail_length=15, head_length=10, 
    108             sld=0.4, head_sld=3.0, solvent_sld=6.0, 
     112            sld=0.4, sld_head=3.0, sld_solvent=6.0, 
    109113            tail_length_pd=0.2, tail_length_pd_n=40, 
    110114            head_length_pd=0.01, head_length_pd_n=40) 
    … …  
    112116oldname = 'LamellarFFHGModel' 
    113117oldpars = dict(head_length='h_thickness', tail_length='t_length', 
    114                sld='sld_tail', head_sld='sld_head', solvent_sld='sld_solvent') 
     118               sld='sld_tail', sld_head='sld_head', sld_solvent='sld_solvent') 
    115119# 
    116120tests = [[{'scale': 1.0, 'background': 0.0, 'tail_length': 15.0, 'head_length': 10.0, 
    117            'sld': 0.4, 'head_sld': 3.0, 'solvent_sld': 6.0}, [0.001], [653143.9209]]] 
     121           'sld': 0.4, 'sld_head': 3.0, 'sld_solvent': 6.0}, [0.001], [653143.9209]]] 
     122# ADDED by: RKH  ON: 18Mar2016  converted from sasview previously, now renaming everything & sorting the docs 

  • sasmodels/models/lamellar_hg_stack_caille.py

    raa2edb2 r6ab4ed8  
    99.. math:: 
    1010 
    11     I(q) = 2 \pi \frac{P(q)S(q)}{\delta q^2} 
     11    I(q) = 2 \pi \frac{P(q)S(q)}{q^2\delta } 
    1212 
    1313 
    … …  
    5656results for the next lower and higher values. 
    5757 
     58Be aware that the computations may be very slow. 
     59 
    5860The 2D scattering intensity is calculated in the same way as 1D, where 
    5961the $q$ vector is defined as 
    … …  
    7375from numpy import inf 
    7476 
    75 name = "lamellarCailleHG" 
    76 title = "Random lamellar sheet with Caille structure factor" 
     77name = "lamellar_hg_stack_caille" 
     78title = "Random lamellar head/tail/tail/head sheet with Caille structure factor" 
    7779description = """\ 
    7880    [Random lamellar phase with Caille  structure factor] 
    … …  
    104106    ["sld", "1e-6/Ang^2", 0.4, [-inf, inf], "", 
    105107     "Tail scattering length density"], 
    106     ["head_sld", "1e-6/Ang^2", 2.0, [-inf, inf], "", 
     108    ["sld_head", "1e-6/Ang^2", 2.0, [-inf, inf], "", 
    107109     "Head scattering length density"], 
    108     ["solvent_sld", "1e-6/Ang^2", 6, [-inf, inf], "", 
     110    ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "", 
    109111     "Solvent scattering length density"], 
    110112    ] 
    111113 
    112 source = ["lamellarCailleHG_kernel.c"] 
     114source = ["lamellar_hg_stack_caille_kernel.c"] 
    113115 
    114116# No volume normalization despite having a volume parameter 
    … …  
    129131    Nlayers=20, spacing=200., Caille_parameter=0.05, 
    130132    tail_length=15, head_length=10, 
    131     #sld=-1, head_sld=4.0, solvent_sld=6.0, 
    132     sld=-1, head_sld=4.1, solvent_sld=6.0, 
     133    #sld=-1, sld_head=4.0, sld_solvent=6.0, 
     134    sld=-1, sld_head=4.1, sld_solvent=6.0, 
    133135    tail_length_pd=0.1, tail_length_pd_n=20, 
    134136    head_length_pd=0.05, head_length_pd_n=30, 
    … …  
    140142oldpars = dict( 
    141143    tail_length='deltaT', head_length='deltaH', Nlayers='n_plates', 
    142     Caille_parameter='caille', sld='sld_tail', head_sld='sld_head', 
    143     solvent_sld='sld_solvent') 
     144    Caille_parameter='caille', sld='sld_tail', sld_head='sld_head', 
     145    sld_solvent='sld_solvent') 
    144146# 
    145147tests = [[{'scale': 1.0, 'background': 0.0, 'tail_length': 10.0, 'head_length': 2.0, 
    146148           'Nlayers': 30.0, 'spacing': 40., 'Caille_parameter': 0.001, 'sld': 0.4, 
    147            'head_sld': 2.0, 'solvent_sld': 6.0, 'tail_length_pd': 0.0, 
     149           'sld_head': 2.0, 'sld_solvent': 6.0, 'tail_length_pd': 0.0, 
    148150           'head_length_pd': 0.0, 'spacing_pd': 0.0}, [0.001], [6838238.571488]]] 
     151# ADDED by: RKH  ON: 18Mar2016  converted from sasview previously, now renaming everything & sorting the docs 

  • sasmodels/models/lamellar_stack_caille.py

    raa2edb2 r6ab4ed8  
    1111.. math:: 
    1212 
    13     I(q) = 2\pi \frac{P(q)S(q)}{\delta q^2} 
     13    I(q) = 2\pi \frac{P(q)S(q)}{q^2\delta } 
    1414 
    1515The form factor is 
    … …  
    7070from numpy import inf 
    7171 
    72 name = "lamellarPS" 
     72name = "lamellar_stack_caille" 
    7373title = "Random lamellar sheet with Caille structure factor" 
    7474description = """\ 
    … …  
    9292    ["Caille_parameter", "1/Ang^2",    0.1, [0.0,0.8],  "",       "Caille parameter"], 
    9393    ["sld",              "1e-6/Ang^2", 6.3, [-inf,inf], "",       "layer scattering length density"], 
    94     ["solvent_sld",      "1e-6/Ang^2", 1.0, [-inf,inf], "",       "Solvent scattering length density"], 
     94    ["sld_solvent",      "1e-6/Ang^2", 1.0, [-inf,inf], "",       "Solvent scattering length density"], 
    9595    ] 
    9696# pylint: enable=bad-whitespace, line-too-long 
    9797 
    98 source = ["lamellarCaille_kernel.c"] 
     98source = ["lamellar_stack_caille_kernel.c"] 
    9999 
    100100# No volume normalization despite having a volume parameter 
    … …  
    113113demo = dict(scale=1, background=0, 
    114114            thickness=67., Nlayers=3.75, spacing=200., 
    115             Caille_parameter=0.268, sld=1.0, solvent_sld=6.34, 
     115            Caille_parameter=0.268, sld=1.0, sld_solvent=6.34, 
    116116            thickness_pd=0.1, thickness_pd_n=100, 
    117117            spacing_pd=0.05, spacing_pd_n=40) 
    … …  
    120120oldpars = dict(thickness='delta', Nlayers='N_plates', 
    121121               Caille_parameter='caille', 
    122                sld='sld_bi', solvent_sld='sld_sol') 
     122               sld='sld_bi', sld_solvent='sld_sol') 
    123123# 
    124124tests = [ 
    125125    [{'scale': 1.0, 'background': 0.0, 'thickness': 30., 'Nlayers': 20.0, 
    126126      'spacing': 400., 'Caille_parameter': 0.1, 'sld': 6.3, 
    127       'solvent_sld': 1.0, 'thickness_pd': 0.0, 'spacing_pd': 0.0}, 
     127      'sld_solvent': 1.0, 'thickness_pd': 0.0, 'spacing_pd': 0.0}, 
    128128     [0.001], [28895.13397]] 
    129129    ] 
     130# ADDED by: RKH  ON: 18Mar2016  converted from sasview previously, now renaming everything & sorting the docs 

  • sasmodels/models/lamellar_stack_paracrystal.py

    raa2edb2 r6ab4ed8  
    1616  *not* the total excluded volume of the paracrystal), 
    1717 
    18 - *sld* $-$ *solvent_sld* is the contrast $\Delta \rho$, 
     18- *sld* $-$ *sld_solvent* is the contrast $\Delta \rho$, 
    1919 
    2020- *thickness* is the layer thickness $t$, 
    … …  
    3636 
    3737The form factor of the bilayer is approximated as the cross section of an 
    38 infinite, planar bilayer of thickness $t$ 
     38infinite, planar bilayer of thickness $t$ (compare the equations for the 
     39lamellar model). 
    3940 
    4041.. math:: 
    … …  
    9495from numpy import inf 
    9596 
    96 name = "lamellarPC" 
     97name = "lamellar_stack_paracrystal" 
    9798title = "Random lamellar sheet with paracrystal structure factor" 
    9899description = """\ 
    … …  
    120121              ["sld", "1e-6/Ang^2", 1.0, [-inf, inf], "", 
    121122               "layer scattering length density"], 
    122               ["solvent_sld", "1e-6/Ang^2", 6.34, [-inf, inf], "", 
     123              ["sld_solvent", "1e-6/Ang^2", 6.34, [-inf, inf], "", 
    123124               "Solvent scattering length density"], 
    124125             ] 
    125126 
    126127 
    127 source = ["lamellarPC_kernel.c"] 
     128source = ["lamellar_stack_paracrystal_kernel.c"] 
    128129 
    129130form_volume = """ 
    … …  
    140141demo = dict(scale=1, background=0, 
    141142            thickness=33, Nlayers=20, spacing=250, spacing_polydisp=0.2, 
    142             sld=1.0, solvent_sld=6.34, 
     143            sld=1.0, sld_solvent=6.34, 
    143144            thickness_pd=0.2, thickness_pd_n=40) 
    144145 
    145146oldname = 'LamellarPCrystalModel' 
    146147oldpars = dict(spacing_polydisp='pd_spacing', sld='sld_layer', 
    147                solvent_sld='sld_solvent') 
     148               sld_solvent='sld_solvent') 
    148149# 
    149150tests = [ 
    150151        [ {'scale': 1.0, 'background' : 0.0, 'thickness' : 33.,'Nlayers' : 20.0, 'spacing' : 250., 'spacing_polydisp' : 0.2, 
    151             'sld' : 1.0, 'solvent_sld' : 6.34, 
     152            'sld' : 1.0, 'sld_solvent' : 6.34, 
    152153            'thickness_pd' : 0.0, 'thickness_pd_n' : 40 }, [0.001, 0.215268], [21829.3, 0.00487686]] 
    153154        ] 
     155# ADDED by: RKH  ON: 18Mar2016  converted from sasview previously, now renaming everything & sorting the docs 

  • sasmodels/models/lib/j1_cephes.c

    re2af2a9 rfad5dc1  
    3939Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier 
    4040*/ 
    41 double j1(double ); 
    42  
    43 double j1(double x) { 
    44  
    45 //Cephes double pression function 
    46 #if FLOAT_SIZE>4 
    47  
    48     double w, z, p, q, xn; 
    49  
    50     const double DR1 = 5.78318596294678452118E0; 
    51     const double DR2 = 3.04712623436620863991E1; 
    52     const double Z1 = 1.46819706421238932572E1; 
    53     const double Z2 = 4.92184563216946036703E1; 
    54     const double THPIO4 =  2.35619449019234492885; 
    55     const double SQ2OPI = 0.79788456080286535588; 
    56  
    57     double RP[8] = { 
     41 
     42constant double RP[8] = { 
    5843    -8.99971225705559398224E8, 
    5944    4.52228297998194034323E11, 
    … …  
    6348    0.0, 
    6449    0.0, 
    65     0.0 
    66     }; 
    67  
    68     double RQ[8] = { 
    69     /* 1.00000000000000000000E0,*/ 
    70     6.20836478118054335476E2, 
    71     2.56987256757748830383E5, 
    72     8.35146791431949253037E7, 
    73     2.21511595479792499675E10, 
    74     4.74914122079991414898E12, 
    75     7.84369607876235854894E14, 
    76     8.95222336184627338078E16, 
    77     5.32278620332680085395E18, 
    78     }; 
    79  
    80     double PP[8] = { 
     50    0.0 }; 
     51 
     52constant double RQ[8] = { 
     53        6.20836478118054335476E2, 
     54        2.56987256757748830383E5, 
     55        8.35146791431949253037E7, 
     56        2.21511595479792499675E10, 
     57        4.74914122079991414898E12, 
     58        7.84369607876235854894E14, 
     59        8.95222336184627338078E16, 
     60        5.32278620332680085395E18 
     61    }; 
     62 
     63constant double PP[8] = { 
    8164    7.62125616208173112003E-4, 
    8265    7.31397056940917570436E-2, 
    … …  
    8669    5.21451598682361504063E0, 
    8770    1.00000000000000000254E0, 
    88     0.0, 
    89     }; 
    90     double PQ[8] = { 
     71    0.0} ; 
     72 
     73 
     74constant double PQ[8] = { 
    9175    5.71323128072548699714E-4, 
    9276    6.88455908754495404082E-2, 
    … …  
    9680    5.20982848682361821619E0, 
    9781    9.99999999999999997461E-1, 
    98     0.0, 
    99     }; 
    100  
    101     double QP[8] = { 
     82    0.0 }; 
     83 
     84constant double QP[8] = { 
    10285    5.10862594750176621635E-2, 
    10386    4.98213872951233449420E0, 
    … …  
    10790    5.97489612400613639965E2, 
    10891    2.11688757100572135698E2, 
    109     2.52070205858023719784E1, 
    110     }; 
    111  
    112     double QQ[8] = { 
    113     /* 1.00000000000000000000E0,*/ 
     92    2.52070205858023719784E1 }; 
     93 
     94constant double QQ[8] = { 
    11495    7.42373277035675149943E1, 
    11596    1.05644886038262816351E3, 
    … …  
    119100    2.82619278517639096600E3, 
    120101    3.36093607810698293419E2, 
    121     0.0, 
    122     }; 
    123  
    124     w = x; 
    125     if( x < 0 ) 
    126             w = -x; 
    127  
    128     if( w <= 5.0 ) 
    129         { 
    130             z = x * x; 
    131             w = polevl( z, RP, 3 ) / p1evl( z, RQ, 8 ); 
    132             w = w * x * (z - Z1) * (z - Z2); 
    133             return( w ); 
    134         } 
    135  
    136     w = 5.0/x; 
    137     z = w * w; 
    138     p = polevl( z, PP, 6)/polevl( z, PQ, 6 ); 
    139     q = polevl( z, QP, 7)/p1evl( z, QQ, 7 ); 
    140     xn = x - THPIO4; 
    141  
    142     double sn, cn; 
    143     SINCOS(xn, sn, cn); 
    144     p = p * cn - w * q * sn; 
    145  
    146     return( p * SQ2OPI / sqrt(x) ); 
    147  
    148  
    149 //Single precission version of cephes 
    150 #else 
    151     double xx, w, z, p, q, xn; 
    152  
    153     const double Z1 = 1.46819706421238932572E1; 
    154     const double THPIO4F =  2.35619449019234492885;    /* 3*pi/4 */ 
    155  
    156  
    157     double JP[8] = { 
     102    0.0 }; 
     103 
     104constant double JP[8] = { 
    158105        -4.878788132172128E-009, 
    159106        6.009061827883699E-007, 
    … …  
    166113    }; 
    167114 
    168     double MO1[8] = { 
     115constant double MO1[8] = { 
    169116        6.913942741265801E-002, 
    170117        -2.284801500053359E-001, 
    … …  
    177124    }; 
    178125 
    179     double PH1[8] = { 
     126constant double PH1[8] = { 
    180127        -4.497014141919556E+001, 
    181128        5.073465654089319E+001, 
    … …  
    188135    }; 
    189136 
     137double J1(double x) { 
     138 
     139//Cephes double pression function 
     140#if FLOAT_SIZE>4 
     141 
     142    double w, z, p, q, xn; 
     143 
     144    const double Z1 = 1.46819706421238932572E1; 
     145    const double Z2 = 4.92184563216946036703E1; 
     146    const double THPIO4 =  2.35619449019234492885; 
     147    const double SQ2OPI = 0.79788456080286535588; 
     148 
     149    w = x; 
     150    if( x < 0 ) 
     151            w = -x; 
     152 
     153    if( w <= 5.0 ) 
     154        { 
     155            z = x * x; 
     156            w = polevl( z, RP, 3 ) / p1evl( z, RQ, 8 ); 
     157            w = w * x * (z - Z1) * (z - Z2); 
     158            return( w ); 
     159        } 
     160 
     161    w = 5.0/x; 
     162    z = w * w; 
     163 
     164    p = polevl( z, PP, 6)/polevl( z, PQ, 6 ); 
     165    q = polevl( z, QP, 7)/p1evl( z, QQ, 7 ); 
     166 
     167    xn = x - THPIO4; 
     168 
     169    double sn, cn; 
     170    SINCOS(xn, sn, cn); 
     171    p = p * cn - w * q * sn; 
     172 
     173    return( p * SQ2OPI / sqrt(x) ); 
     174 
     175 
     176//Single precission version of cephes 
     177#else 
     178    double xx, w, z, p, q, xn; 
     179 
     180    const double Z1 = 1.46819706421238932572E1; 
     181    const double THPIO4F =  2.35619449019234492885;    /* 3*pi/4 */ 
     182 
     183 
    190184    xx = x; 
    191185    if( xx < 0 ) 
    … …  
    211205} 
    212206 
     207 
     208/*double J1c(double x) { 
     209    return 
     210}*/ 

  • sasmodels/models/lib/polevl.c

    re2af2a9 rfad5dc1  
    5151*/ 
    5252 
    53 double polevl( double x, double coef[8], int N ); 
    54 double p1evl( double x, double coef[8], int N ); 
    5553 
    56 double polevl( double x, double coef[8], int N ) { 
     54double polevl( double x, constant double *coef, int N ) { 
    5755 
    5856    int i = 0; 
    … …  
    7472 */ 
    7573 
    76 double p1evl( double x, double coef[8], int N ) { 
     74double p1evl( double x, constant double *coef, int N ) { 
    7775    int i=0; 
    7876    double ans = x+coef[i]; 
    … …  
    8482 
    8583    return( ans ); 
    86  
    8784} 

  • sasmodels/product.py

    r3bcd03d rd5ba841  
    124124        # so borrow values from end of p_fixed.  This makes volfraction the 
    125125        # first S parameter. 
    126         start += num_p_fixed - 2 
    127         par_map['s_fixed'] = np.arange(start, start+num_s_fixed) 
     126        start += num_p_fixed 
     127        par_map['s_fixed'] = np.hstack(([start,start], 
     128                                        np.arange(start, start+num_s_fixed-2))) 
    128129        par_map['volfraction'] = num_p_fixed 
    129         start += num_s_fixed 
     130        start += num_s_fixed-2 
    130131        # vol pars offset from the start of pd pars 
    131132        par_map['vol_pars'] = [start+k for k in vol_par_idx] 
    132133        par_map['p_pd'] = np.arange(start, start+num_p_pd) 
    133         start += num_p_pd 
    134         par_map['s_pd'] = np.arange(start-1, start+num_s_pd)  # should be empty... 
     134        start += num_p_pd-1 
     135        par_map['s_pd'] = np.hstack((start, 
     136                                     np.arange(start, start+num_s_pd-1))) 
    135137 
    136138        self.fixed_pars = model_info['partype']['fixed-' + dim] 

  • sasmodels/sasview_model.py

    r2622b3f r08376e7  
    2323from . import core 
    2424 
    25 def make_class(model_info, dtype='single', namestyle='name'): 
     25def standard_models(): 
     26    return [make_class(model_name) for model_name in core.list_models()] 
     27 
     28def make_class(model_name, namestyle='name'): 
    2629    """ 
    2730    Load the sasview model defined in *kernel_module*. 
    2831 
    29     Returns a class that can be used directly as a sasview model. 
     32    Returns a class that can be used directly as a sasview model.t 
    3033 
    3134    Defaults to using the new name for a model.  Setting 
    … …  
    3336    compatible with SasView. 
    3437    """ 
    35     model = core.build_model(model_info, dtype=dtype) 
     38    model_info = core.load_model_info(model_name) 
    3639    def __init__(self, multfactor=1): 
    37         SasviewModel.__init__(self, model) 
    38     attrs = dict(__init__=__init__) 
    39     ConstructedModel = type(model.info[namestyle], (SasviewModel,), attrs) 
     40        SasviewModel.__init__(self) 
     41    attrs = dict(__init__=__init__, _model_info=model_info) 
     42    ConstructedModel = type(model_info[namestyle], (SasviewModel,), attrs) 
    4043    return ConstructedModel 
    4144 
    … …  
    4447    Sasview wrapper for opencl/ctypes model. 
    4548    """ 
    46     def __init__(self, model): 
    47         """Initialization""" 
    48         self._model = model 
    49  
    50         self.name = model.info['name'] 
    51         self.oldname = model.info['oldname'] 
    52         self.description = model.info['description'] 
     49    def __init__(self): 
     50        self._kernel = None 
     51        model_info = self._model_info 
     52 
     53        self.name = model_info['name'] 
     54        self.oldname = model_info['oldname'] 
     55        self.description = model_info['description'] 
    5356        self.category = None 
    5457        self.multiplicity_info = None 
    … …  
    6063        self.params = collections.OrderedDict() 
    6164        self.dispersion = dict() 
    62         partype = model.info['partype'] 
    63  
    64         for p in model.info['parameters']: 
     65        partype = model_info['partype'] 
     66 
     67        for p in model_info['parameters']: 
    6568            self.params[p.name] = p.default 
    6669            self.details[p.name] = [p.units] + p.limits 
    … …  
    8386 
    8487        ## independent parameter name and unit [string] 
    85         self.input_name = model.info.get("input_name", "Q") 
    86         self.input_unit = model.info.get("input_unit", "A^{-1}") 
    87         self.output_name = model.info.get("output_name", "Intensity") 
    88         self.output_unit = model.info.get("output_unit", "cm^{-1}") 
     88        self.input_name = model_info.get("input_name", "Q") 
     89        self.input_unit = model_info.get("input_unit", "A^{-1}") 
     90        self.output_name = model_info.get("output_name", "Intensity") 
     91        self.output_unit = model_info.get("output_unit", "cm^{-1}") 
    8992 
    9093        ## _persistency_dict is used by sas.perspectives.fitting.basepage 
    … …  
    9598        ## New fields introduced for opencl rewrite 
    9699        self.cutoff = 1e-5 
     100 
     101    def __get_state__(self): 
     102        state = self.__dict__.copy() 
     103        model_id = self._model_info['id'] 
     104        state.pop('_kernel') 
     105        # May need to reload model info on set state since it has pointers 
     106        # to python implementations of Iq, etc. 
     107        #state.pop('_model_info') 
     108        return state 
     109 
     110    def __set_state__(self, state): 
     111        self.__dict__ = state 
     112        self._kernel = None 
    97113 
    98114    def __str__(self): 
    … …  
    187203        # TODO: fix test so that parameter order doesn't matter 
    188204        ret = ['%s.%s' % (d.lower(), p) 
    189                for d in self._model.info['partype']['pd-2d'] 
     205               for d in self._model_info['partype']['pd-2d'] 
    190206               for p in ('npts', 'nsigmas', 'width')] 
    191207        #print(ret) 
    … …  
    261277            # Check whether we have a list of ndarrays [qx,qy] 
    262278            qx, qy = qdist 
    263             partype = self._model.info['partype'] 
     279            partype = self._model_info['partype'] 
    264280            if not partype['orientation'] and not partype['magnetic']: 
    265281                return self.calculate_Iq(np.sqrt(qx ** 2 + qy ** 2)) 
    … …  
    284300        to the card for each evaluation. 
    285301        """ 
     302        if self._kernel is None: 
     303            self._kernel = core.build_model(self._model_info) 
    286304        q_vectors = [np.asarray(q) for q in args] 
    287         fn = self._model(q_vectors) 
     305        fn = self._kernel(q_vectors) 
    288306        pars = [self.params[v] for v in fn.fixed_pars] 
    289307        pd_pars = [self._get_weights(p) for p in fn.pd_pars] 
    … …  
    299317        :return: the value of the effective radius 
    300318        """ 
    301         ER = self._model.info.get('ER', None) 
     319        ER = self._model_info.get('ER', None) 
    302320        if ER is None: 
    303321            return 1.0 
    … …  
    314332        :return: the value of the volf ratio 
    315333        """ 
    316         VR = self._model.info.get('VR', None) 
     334        VR = self._model_info.get('VR', None) 
    317335        if VR is None: 
    318336            return 1.0 
    … …  
    358376        parameter set in the vector. 
    359377        """ 
    360         pars = self._model.info['partype']['volume'] 
     378        pars = self._model_info['partype']['volume'] 
    361379        return core.dispersion_mesh([self._get_weights(p) for p in pars]) 
    362380 
    … …  
    368386        from . import weights 
    369387 
    370         relative = self._model.info['partype']['pd-rel'] 
    371         limits = self._model.info['limits'] 
     388        relative = self._model_info['partype']['pd-rel'] 
     389        limits = self._model_info['limits'] 
    372390        dis = self.dispersion[par] 
    373391        value, weight = weights.get_weights( 
    … …  
    375393            self.params[par], limits[par], par in relative) 
    376394        return value, weight / np.sum(weight) 
     395 

  • example/sesansfit.py

    r170ea69 r89ab393  
    22from sasmodels import core, bumps_model, sesans 
    33from sas.sascalc.dataloader.loader import Loader 
    4  
    5 HAS_CONVERTER = True 
    6 try: 
    7     from sas.sascalc.data_util.nxsunit import Converter 
    8 except ImportError: 
    9     HAS_CONVERTER = False 
    10  
    114 
    125def get_bumps_model(model_name): 
    … …  
    3023        data = loader.load(file) 
    3124        if data is None: raise IOError("Could not load file %r"%(file)) 
    32         if HAS_CONVERTER == True: 
    33             default_unit = "A" 
    34             data_conv_q = Converter(data._xunit) 
    35             for x in data.x: 
    36                 print x 
    37             data.x = data_conv_q(data.x, units=default_unit) 
    38             for x in data.x: 
    39                 print x 
    40             data._xunit = default_unit 
    4125 
    4226    except: 
Note: See TracChangeset for help on using the changeset viewer.