Changeset b3f6bc3 in sasmodels


Ignore:
Timestamp:
Feb 15, 2015 4:07:00 PM (8 years ago)
Author:
pkienzle
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:
f786ff3
Parents:
cb6ecf4
Message:

support pure python model Iq/Iqxy?

Location:
sasmodels
Files:
1 added
4 edited

Legend:

Unmodified
Added
Removed
  • sasmodels/dll.py

    rcb6ecf4 rb3f6bc3  
    1111 
    1212from . import generate 
     13from .pykernel import PyInput, PyKernel 
    1314 
    1415from .generate import F32, F64 
     
    116117 
    117118    def __call__(self, input): 
     119        # Support pure python kernel call 
     120        if input.is_2D and callable(self.info['Iqxy']): 
     121            return PyKernel(self.info['Iqxy'], self.info, input) 
     122        elif not input.is_2D and callable(self.info['Iq']): 
     123            return PyKernel(self.info['Iq'], self.info, input) 
     124 
    118125        if self.dll is None: self._load_dll() 
    119  
    120126        kernel = self.Iqxy if input.is_2D else self.Iq 
    121127        return DllKernel(kernel, self.info, input) 
     
    125131        Make q input vectors available to the model. 
    126132 
    127         This only needs to be done once for all models that operate on the 
    128         same input.  So for example, if you are adding two different models 
    129         together to compare to a data set, then only one model needs to 
    130         needs to call make_input, so long as the models have the same dtype. 
     133        Note that each model needs its own q vector even if the case of 
     134        mixture models because some models may be OpenCL, some may be 
     135        ctypes and some may be pure python. 
    131136        """ 
    132         return DllInput(q_vectors) 
     137        return PyInput(q_vectors, dtype=F64) 
    133138 
    134139    def release(self): 
     
    136141 
    137142 
    138 class DllInput(object): 
    139     """ 
    140     Make q data available to the gpu. 
    141  
    142     *q_vectors* is a list of q vectors, which will be *[q]* for 1-D data, 
    143     and *[qx, qy]* for 2-D data.  Internally, the vectors will be reallocated 
    144     to get the best performance on OpenCL, which may involve shifting and 
    145     stretching the array to better match the memory architecture.  Additional 
    146     points will be evaluated with *q=1e-3*. 
    147  
    148     *dtype* is the data type for the q vectors. The data type should be 
    149     set to match that of the kernel, which is an attribute of 
    150     :class:`GpuProgram`.  Note that not all kernels support double 
    151     precision, so even if the program was created for double precision, 
    152     the *GpuProgram.dtype* may be single precision. 
    153  
    154     Call :meth:`release` when complete.  Even if not called directly, the 
    155     buffer will be released when the data object is freed. 
    156     """ 
    157     def __init__(self, q_vectors): 
    158         self.nq = q_vectors[0].size 
    159         self.dtype = np.dtype('double') 
    160         self.is_2D = (len(q_vectors) == 2) 
    161         self.q_vectors = [np.ascontiguousarray(q,self.dtype) for q in q_vectors] 
    162         self.q_pointers = [q.ctypes.data for q in q_vectors] 
    163  
    164     def release(self): 
    165         self.q_vectors = [] 
    166  
    167143class DllKernel(object): 
    168144    """ 
    169145    Callable SAS kernel. 
    170146 
    171     *kernel* is the DllKernel object to call. 
     147    *kernel* is the c function to call. 
    172148 
    173149    *info* is the module information 
  • sasmodels/generate.py

    rcb6ecf4 rb3f6bc3  
    354354""" 
    355355 
     356 
     357 
     358########################################################## 
     359#                                                        # 
     360#   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!   # 
     361#   !!                                              !!   # 
     362#   !!  KEEP THIS CODE CONSISTENT WITH PYKERNEL.PY  !!   # 
     363#   !!                                              !!   # 
     364#   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!   # 
     365#                                                        # 
     366########################################################## 
     367 
     368 
    356369# Polydispersity loop body. 
    357370# This computes the weight, and if it is sufficient, calls the scattering 
     
    387400# to call the form_volume function from the user supplied kernel, and accumulate 
    388401# a normalized weight. 
    389 VOLUME_NORM="""const double vol_weight = %(weight)s; 
    390     vol += vol_weight*form_volume(%(pars)s); 
     402VOLUME_NORM="""const double vol_weight = %(vol_weight)s; 
     403    vol += vol_weight*form_volume(%(vol_pars)s); 
    391404    norm_vol += vol_weight;\ 
    392405""" 
     
    490503    if vol_pars: 
    491504        subst = { 
    492             'weight': "*".join(p+"_w" for p in vol_pars), 
    493             'pars': ", ".join(vol_pars), 
     505            'vol_weight': "*".join(p+"_w" for p in vol_pars), 
     506            'vol_pars': ", ".join(vol_pars), 
    494507            } 
    495508        volume_norm = VOLUME_NORM%subst 
     
    637650            } 
    638651        source.append(WORK_FUNCTION%subst) 
    639     kernel_Iq = make_kernel(info, is_2D=False) 
    640     kernel_Iqxy = make_kernel(info, is_2D=True) 
     652    kernel_Iq = make_kernel(info, is_2D=False) if not callable(info['Iq']) else "" 
     653    kernel_Iqxy = make_kernel(info, is_2D=True) if not callable(info['Iqxy']) else "" 
    641654    kernel = "\n\n".join([KERNEL_HEADER]+source+[kernel_Iq, kernel_Iqxy]) 
    642655    return kernel 
     
    718731def demo_time(): 
    719732    import datetime 
     733    from .models import cylinder 
     734    toc = lambda: (datetime.datetime.now()-tic).total_seconds() 
    720735    tic = datetime.datetime.now() 
    721     toc = lambda: (datetime.datetime.now()-tic).total_seconds() 
    722     path = os.path.dirname("__file__") 
    723     doc, c = make_model(os.path.join(path, "models", "cylinder.c")) 
     736    source, info = make(cylinder) 
    724737    print "time:",toc() 
    725738 
    726739def demo(): 
    727     from os.path import join as joinpath, dirname 
    728     c, info, doc = make_model(joinpath(dirname(__file__), "models", "cylinder.c")) 
    729     #print doc 
    730     #print c 
     740    from .models import cylinder 
     741    source, info = make(cylinder) 
     742    #print doc(cylinder) 
     743    print source 
    731744 
    732745if __name__ == "__main__": 
  • sasmodels/gpu.py

    rcb6ecf4 rb3f6bc3  
    2424""" 
    2525import numpy as np 
     26 
    2627import pyopencl as cl 
    2728from pyopencl import mem_flags as mf 
    2829 
    2930from . import generate 
     31from .pykernel import PyInput, PyKernel 
    3032 
    3133F64_DEFS = """\ 
     
    185187 
    186188    def __call__(self, input): 
     189        # Support pure python kernel call 
     190        if input.is_2D and callable(self.info['Iqxy']): 
     191            return PyKernel(self.info['Iqxy'], self.info, input) 
     192        elif not input.is_2D and callable(self.info['Iq']): 
     193            return PyKernel(self.info['Iq'], self.info, input) 
     194 
    187195        if self.dtype != input.dtype: 
    188196            raise TypeError("data and kernel have different types") 
     
    202210        Make q input vectors available to the model. 
    203211 
    204         This only needs to be done once for all models that operate on the 
    205         same input.  So for example, if you are adding two different models 
    206         together to compare to a data set, then only one model needs to 
    207         needs to call make_input, so long as the models have the same dtype. 
     212        Note that each model needs its own q vector even if the case of 
     213        mixture models because some models may be OpenCL, some may be 
     214        ctypes and some may be pure python. 
    208215        """ 
    209         # Note: the weird interface, where make_input doesn't care which 
    210         # model calls it, allows us to ask the model to define the data 
    211         # and the caller does not need to know if it is opencl or ctypes. 
    212         # The returned data object is opaque. 
    213         return GpuInput(q_vectors, dtype=self.dtype) 
     216        # Support pure python kernel call 
     217        if len(q_vectors) == 1 and callable(self.info['Iq']): 
     218            return PyInput(q_vectors, dtype=self.dtype) 
     219        elif callable(self.info['Iqxy']): 
     220            return PyInput(q_vectors, dtype=self.dtype) 
     221        else: 
     222            return GpuInput(q_vectors, dtype=self.dtype) 
    214223 
    215224# TODO: check that we don't need a destructor for buffers which go out of scope 
  • sasmodels/models/spherepy.py

    rbe802cb rb3f6bc3  
    8383 
    8484 
    85 # No volume normalization despite having a volume parameter 
    86 # This should perhaps be volume normalized? 
    8785def form_volume(radius): 
    88     return 1.333333333333333*pi*radius*radius*radius 
     86    return 1.333333333333333*pi*radius**3 
    8987 
    9088def Iq(q, sld, solvent_sld, radius): 
     89    #print "q",q 
     90    #print "sld,r",sld,solvent_sld,radius 
    9191    qr = q*radius 
    9292    sn, cn = sin(qr), cos(qr) 
    93     bes = 3 * (sn-qr*cn)/qr**3 
     93    # FOR VECTORIZED VERSION, UNCOMMENT THE NEXT TWO LINES 
     94    bes = 3 * (sn-qr*cn)/qr**3 # may be 0/0 but we fix that next line 
    9495    bes[qr==0] = 1 
     96    # FOR NON VECTORIZED VERSION, UNCOMMENT THE NEXT LINE 
     97    #bes = 3 * (sn-qr*cn)/qr**3 if qr>0 else 1 
    9598    fq = bes * (sld - solvent_sld) * form_volume(radius) 
    96     return 1.0e-4*fq*fq 
     99    return 1.0e-4*fq**2 
     100# FOR VECTORIZED VERSION, UNCOMMENT THE NEXT LINE 
     101Iq.vectorized = True 
    97102 
    98103def Iqxy(qx, qy, sld, solvent_sld, radius): 
    99104    return Iq(sqrt(qx**2 + qy**2), sld, solvent_sld, radius) 
     105Iqxy.vectorized = True 
    100106 
    101107def ER(radius): 
Note: See TracChangeset for help on using the changeset viewer.