source: sasmodels/sasmodels/core.py @ 5d316e9

core_shell_microgelscostrafo411magnetic_modelrelease_v0.94release_v0.95ticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 5d316e9 was 5d316e9, checked in by Paul Kienzle <pkienzle@…>, 7 years ago

support fast and loose single precision and half precision

  • Property mode set to 100644
File size: 6.9 KB
Line 
1"""
2Core model handling routines.
3"""
4__all__ = ["list_models", "load_model_definition",  "precompile_dll",
5           "load_model", "make_kernel", "call_kernel", "call_ER", "call_VR" ]
6
7from os.path import basename, dirname, join as joinpath
8from glob import glob
9
10import numpy as np
11
12from . import models
13from . import weights
14from . import generate
15
16from . import kernelpy
17from . import kerneldll
18try:
19    from . import kernelcl
20    HAVE_OPENCL = True
21except:
22    HAVE_OPENCL = False
23
24
25def list_models():
26    """
27    Return the list of available models on the model path.
28    """
29    root = dirname(__file__)
30    files = sorted(glob(joinpath(root, 'models', "[a-zA-Z]*.py")))
31    available_models = [basename(f)[:-3] for f in files]
32    return available_models
33
34
35def load_model_definition(model_name):
36    """
37    Load a model definition given the model name.
38    """
39    __import__('sasmodels.models.'+model_name)
40    model_definition = getattr(models, model_name, None)
41    return model_definition
42
43
44def precompile_dll(model_name, dtype="double"):
45    """
46    Precompile the dll for a model.
47
48    Returns the path to the compiled model.
49
50    This can be used when build the windows distribution of sasmodels
51    (which may be missing the OpenCL driver and the dll compiler), or
52    otherwise sharing models with windows users who do not have a compiler.
53
54    See :func:`sasmodels.kerneldll.make_dll` for details on controlling the
55    dll path and the allowed floating point precision.
56    """
57    model_definition = load_model_definition(model_name)
58    source, info = generate.make(model_definition)
59    return kerneldll.make_dll(source, info, dtype=dtype)
60
61
62def isstr(s):
63    try: s + ''
64    except: return False
65    return True
66
67def load_model(model_definition, dtype="single", platform="ocl", fast=False):
68    """
69    Prepare the model for the default execution platform.
70
71    This will return an OpenCL model, a DLL model or a python model depending
72    on the model and the computing platform.
73
74    *model_definition* is the python module which defines the model.  If the
75    model name is given instead, then :func:`load_model_definition` will be
76    called with the model name.
77
78    *dtype* indicates whether the model should use single or double precision
79    for the calculation. Any valid numpy single or double precision identifier
80    is valid, such as 'single', 'f', 'f32', or np.float32 for single, or
81    'double', 'd', 'f64'  and np.float64 for double.
82
83    *platform* should be "dll" to force the dll to be used for C models,
84    otherwise it uses the default "ocl".
85
86    *fast* is True if fast inaccurate math is acceptable (40% speed increase).
87    """
88    if isstr(model_definition):
89        model_definition = load_model_definition(model_definition)
90    source, info = generate.make(model_definition)
91    if callable(info.get('Iq', None)):
92        return kernelpy.PyModel(info)
93
94    ## for debugging:
95    ##  1. uncomment open().write so that the source will be saved next time
96    ##  2. run "python -m sasmodels.direct_model $MODELNAME" to save the source
97    ##  3. recomment the open.write() and uncomment open().read()
98    ##  4. rerun "python -m sasmodels.direct_model $MODELNAME"
99    ##  5. uncomment open().read() so that source will be regenerated from model
100    # open(info['name']+'.c','w').write(source)
101    # source = open(info['name']+'.cl','r').read()
102
103    dtype = np.dtype(dtype)
104    if (platform=="dll"
105            or not HAVE_OPENCL
106            or not kernelcl.environment().has_type(dtype)):
107        return kerneldll.load_dll(source, info, dtype)
108    else:
109        return kernelcl.GpuModel(source, info, dtype, fast)
110
111def make_kernel(model, q_vectors):
112    """
113    Return a computation kernel from the model definition and the q input.
114    """
115    model_input = model.make_input(q_vectors)
116    return model(model_input)
117
118def get_weights(info, pars, name):
119    """
120    Generate the distribution for parameter *name* given the parameter values
121    in *pars*.
122
123    Uses "name", "name_pd", "name_pd_type", "name_pd_n", "name_pd_sigma"
124    from the *pars* dictionary for parameter value and parameter dispersion.
125    """
126    relative = name in info['partype']['pd-rel']
127    limits = info['limits'][name]
128    disperser = pars.get(name+'_pd_type', 'gaussian')
129    value = pars.get(name, info['defaults'][name])
130    npts = pars.get(name+'_pd_n', 0)
131    width = pars.get(name+'_pd', 0.0)
132    nsigma = pars.get(name+'_pd_nsigma', 3.0)
133    value,weight = weights.get_weights(
134        disperser, npts, width, nsigma, value, limits, relative)
135    return value, weight / np.sum(weight)
136
137def dispersion_mesh(pars):
138    """
139    Create a mesh grid of dispersion parameters and weights.
140
141    Returns [p1,p2,...],w where pj is a vector of values for parameter j
142    and w is a vector containing the products for weights for each
143    parameter set in the vector.
144    """
145    value, weight = zip(*pars)
146    if len(value) > 1:
147        value = [v.flatten() for v in np.meshgrid(*value)]
148        weight = np.vstack([v.flatten() for v in np.meshgrid(*weight)])
149        weight = np.prod(weight, axis=0)
150    return value, weight
151
152def call_kernel(kernel, pars, cutoff=0):
153    """
154    Call *kernel* returned from :func:`make_kernel` with parameters *pars*.
155
156    *cutoff* is the limiting value for the product of dispersion weights used
157    to perform the multidimensional dispersion calculation more quickly at a
158    slight cost to accuracy. The default value of *cutoff=0* integrates over
159    the entire dispersion cube.  Using *cutoff=1e-5* can be 50% faster, but
160    with an error of about 1%, which is usually less than the measurement
161    uncertainty.
162    """
163    fixed_pars = [pars.get(name, kernel.info['defaults'][name])
164                  for name in kernel.fixed_pars]
165    pd_pars = [get_weights(kernel.info, pars, name) for name in kernel.pd_pars]
166    return kernel(fixed_pars, pd_pars, cutoff=cutoff)
167
168def call_ER(info, pars):
169    """
170    Call the model ER function using *pars*.
171
172    *info* is either *model.info* if you have a loaded model, or *kernel.info*
173    if you have a model kernel prepared for evaluation.
174    """
175    ER = info.get('ER', None)
176    if ER is None:
177        return 1.0
178    else:
179        vol_pars = [get_weights(info, pars, name)
180                    for name in info['partype']['volume']]
181        value, weight = dispersion_mesh(vol_pars)
182        individual_radii = ER(*value)
183        #print(values[0].shape, weights.shape, fv.shape)
184        return np.sum(weight*individual_radii) / np.sum(weight)
185
186def call_VR(info, pars):
187    """
188    Call the model VR function using *pars*.
189
190    *info* is either *model.info* if you have a loaded model, or *kernel.info*
191    if you have a model kernel prepared for evaluation.
192    """
193    VR = info.get('VR', None)
194    if VR is None:
195        return 1.0
196    else:
197        vol_pars = [get_weights(info, pars, name)
198                    for name in info['partype']['volume']]
199        value, weight = dispersion_mesh(vol_pars)
200        whole,part = VR(*value)
201        return np.sum(weight*part)/np.sum(weight*whole)
202
Note: See TracBrowser for help on using the repository browser.