Changes in / [4c1bfb3:bfea0c2] in sasmodels


Ignore:
Files:
8 added
7 deleted
6 edited

Legend:

Unmodified
Added
Removed

  • example/sesansfit.py

    r4c1bfb3 r4c1bfb3  
    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): 
    … …  
    3124        data = loader.load(file) 
    3225        if data is None: raise IOError("Could not load file %r"%(file)) 
    33         if HAS_CONVERTER == True: 
    34             default_unit = "A" 
    35             data_conv_q = Converter(data._xunit) 
    36             for x in data.x: 
    37                 print x 
    38             data.x = data_conv_q(data.x, units=default_unit) 
    39             for x in data.x: 
    40                 print x 
    41             data._xunit = default_unit 
    4226 
    4327    except: 

  • 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

    re98c1e0 r934f906  
    7575Iq = """ 
    7676      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 
    7786 
    7887      if(fabs(radius_effective) < 1.E-12) { 
    … …  
    97106      G=0.5*volfraction*A; 
    98107 
    99       if(X < 0.2) { 
    100       // 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 
    101110      // else no obvious way to rearrange the equations to avoid needing a very high number of significant figures.  
    102111      // Series expansion found using Mathematica software. Numerical test in .xls showed terms to X^2 are sufficient  

  • 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/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 
Note: See TracChangeset for help on using the changeset viewer.