Changeset d18582e in sasmodels

Timestamp:

Feb 4, 2016 10:44:23 AM (8 years ago)

Author:

Paul Kienzle <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:

5054e80

Parents:

bb6f0f3

Message:

default to double precision if single=False is set in model file

Location:

sasmodels

Files:

• core.py (modified) (3 diffs)
• data.py (modified) (4 diffs)
• direct_model.py (modified) (1 diff)
• kernelcl.py (modified) (11 diffs)
• models/pearl_necklace.py (modified) (1 diff)

sasmodels/core.py

@@ -73 +73 @@
     return True
 
-def load_model(model_definition, dtype="single", platform="ocl"):
+def load_model(model_definition, dtype=None, platform="ocl"):
     """
     Prepare the model for the default execution platform.
@@ -87 +87 @@
     for the calculation. Any valid numpy single or double precision identifier
     is valid, such as 'single', 'f', 'f32', or np.float32 for single, or
-    'double', 'd', 'f64'  and np.float64 for double.
+    'double', 'd', 'f64'  and np.float64 for double.  If *None*, then use
+    'single' unless the model defines single=False.
 
     *platform* should be "dll" to force the dll to be used for C models,
@@ -94 +95 @@
     if isstr(model_definition):
         model_definition = load_model_definition(model_definition)
+    if dtype is None:
+        dtype = 'single' if getattr(model_definition, 'single', True) else 'double'
     source, info = generate.make(model_definition)
     if callable(info.get('Iq', None)):

sasmodels/data.py

@@ -242 +242 @@
 
 
-def empty_data1D(q, resolution=0.05):
+def empty_data1D(q, resolution=0.0):
     """
     Create empty 1D data using the given *q* as the x value.
@@ -252 +252 @@
     #dIq = np.sqrt(Iq)
     Iq, dIq = None, None
+    q = np.asarray(q)
     data = Data1D(q, Iq, dx=resolution * q, dy=dIq)
     data.filename = "fake data"
@@ -257 +258 @@
 
 
-def empty_data2D(qx, qy=None, resolution=0.05):
+def empty_data2D(qx, qy=None, resolution=0.0):
     """
     Create empty 2D data using the given mesh.
@@ -267 +268 @@
     if qy is None:
         qy = qx
+    qx, qy = np.asarray(qx), np.asarray(qy)
     # 5% dQ/Q resolution
     Qx, Qy = np.meshgrid(qx, qy)

sasmodels/direct_model.py

@@ -234 +234 @@
 
     model_definition = load_model_definition(model_name)
-    model = load_model(model_definition, dtype='single')
+    model = load_model(model_definition)
     calculator = DirectModel(data, model)
     pars = dict((k, float(v))

sasmodels/kernelcl.py

@@ -172 +172 @@
         #self.data_boundary = max(d.min_data_type_align_size
         #                         for d in self.context.devices)
-        self.queues = [cl.CommandQueue(self.context, d)
-                       for d in self.context.devices]
+        self.queues = [cl.CommandQueue(context, context.devices[0])
+                       for context in self.context]
         self.compiled = {}
 
@@ -181 +181 @@
         """
         dtype = generate.F32 if dtype == 'fast' else np.dtype(dtype)
-        return all(has_type(d, dtype) for d in self.context.devices)
+        return any(has_type(d, dtype)
+                   for context in self.context
+                   for d in context.devices)
+
+    def get_queue(self, dtype):
+        """
+        Return a command queue for the kernels of type dtype.
+        """
+        for context, queue in zip(self.context, self.queues):
+            if all(has_type(d, dtype) for d in context.devices):
+                return queue
+
+    def get_context(self, dtype):
+        """
+        Return a OpenCL context for the kernels of type dtype.
+        """
+        for context, queue in zip(self.context, self.queues):
+            if all(has_type(d, dtype) for d in context.devices):
+                return context
 
     def _create_some_context(self):
@@ -190 +208 @@
         """
         try:
-            self.context = cl.create_some_context(interactive=False)
+            self.context = [cl.create_some_context(interactive=False)]
         except Exception as exc:
             warnings.warn(str(exc))
@@ -204 +222 @@
             #print("compiling",name)
             dtype = np.dtype(dtype)
-            program = compile_model(self.context, source, dtype, fast)
+            program = compile_model(self.get_context(dtype), source, dtype, fast)
             self.compiled[key] = program
         return self.compiled[key]
@@ -218 +236 @@
 def _get_default_context():
     """
-    Get an OpenCL context, preferring GPU over CPU.
-    """
-    default = None
+    Get an OpenCL context, preferring GPU over CPU, and preferring Intel
+    drivers over AMD drivers.
+    """
+    # Note: on mobile devices there is automatic clock scaling if either the
+    # CPU or the GPU is underutilized; probably doesn't affect us, but we if
+    # it did, it would mean that putting a busy loop on the CPU while the GPU
+    # is running may increase throughput.
+    #
+    # Macbook pro, base install:
+    #     {'Apple': [Intel CPU, NVIDIA GPU]}
+    # Macbook pro, base install:
+    #     {'Apple': [Intel CPU, Intel GPU]}
+    # 2 x nvidia 295 with Intel and NVIDIA opencl drivers installed
+    #     {'Intel': [CPU], 'NVIDIA': [GPU, GPU, GPU, GPU]}
+    gpu, cpu = None, None
     for platform in cl.get_platforms():
         for device in platform.get_devices():
             if device.type == cl.device_type.GPU:
-                return cl.Context([device])
-            if default is None:
-                default = device
-
-    if not default:
-        raise RuntimeError("OpenCL device not found")
-
-    return cl.Context([default])
+                gpu = device
+            else:
+                cpu = device
+    single = gpu if gpu is not None else cpu
+    double = gpu if has_type(gpu, np.dtype('double')) else cpu
+    single_context = cl.Context([single])
+    double_context = cl.Context([double]) if single != double else single_context
+    return single_context, double_context
 
 
@@ -314 +344 @@
         # architectures tested so far.
         self.q_vectors = [_stretch_input(q, self.dtype, 32) for q in q_vectors]
+        context = env.get_context(self.dtype)
         self.q_buffers = [
-            cl.Buffer(env.context, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=q)
+            cl.Buffer(context, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=q)
             for q in self.q_vectors
         ]
@@ -363 +394 @@
         # Note: res may be shorter than res_b if global_size != nq
         env = environment()
-        self.loops_b = [cl.Buffer(env.context, mf.READ_WRITE,
-                                  2 * MAX_LOOPS * q_input.dtype.itemsize)
-                        for _ in env.queues]
-        self.res_b = [cl.Buffer(env.context, mf.READ_WRITE,
-                                q_input.global_size[0] * q_input.dtype.itemsize)
-                      for _ in env.queues]
+        self.queue = env.get_queue(dtype)
+        self.loops_b = cl.Buffer(self.queue.context, mf.READ_WRITE,
+                                 2 * MAX_LOOPS * q_input.dtype.itemsize)
+        self.res_b = cl.Buffer(self.queue.context, mf.READ_WRITE,
+                               q_input.global_size[0] * q_input.dtype.itemsize)
         self.q_input = q_input
+
+        self._need_release = [self.loops_b, self.res_b, self.q_input]
 
     def __call__(self, fixed_pars, pd_pars, cutoff=1e-5):
@@ -377 +409 @@
                 else np.float32)  # will never get here, so use np.float32
 
-        device_num = 0
-        queuei = environment().queues[device_num]
-        res_bi = self.res_b[device_num]
+        res_bi = self.res_b
         nq = np.uint32(self.q_input.nq)
         if pd_pars:
@@ -394 +424 @@
                 raise ValueError("too many polydispersity points")
 
-            loops_bi = self.loops_b[device_num]
-            cl.enqueue_copy(queuei, loops_bi, loops)
+            loops_bi = self.loops_b
+            cl.enqueue_copy(self.queue, loops_bi, loops)
             loops_l = cl.LocalMemory(len(loops.data))
             #ctx = environment().context
@@ -404 +434 @@
         fixed = [real(p) for p in fixed_pars]
         args = self.q_input.q_buffers + [res_bi, nq] + dispersed + fixed
-        self.kernel(queuei, self.q_input.global_size, None, *args)
-        cl.enqueue_copy(queuei, self.res, res_bi)
+        self.kernel(self.queue, self.q_input.global_size, None, *args)
+        cl.enqueue_copy(self.queue, self.res, res_bi)
 
         return self.res
@@ -413 +443 @@
         Release resources associated with the kernel.
         """
-        for b in self.loops_b:
-            b.release()
-        self.loops_b = []
-        for b in self.res_b:
-            b.release()
-        self.res_b = []
-        self.q_input.release()
+        for v in self._need_release:
+            v.release()
+        self._need_release = []
 
     def __del__(self):

sasmodels/models/pearl_necklace.py

@@ -96 +96 @@
 
 source = ["lib/Si.c", "pearl_necklace.c"]
-# new flag to let the compiler know to never use single precision
-single = False
+single = False  # use double precision unless told otherwise
 
 def volume(radius, edge_separation, string_thickness, number_of_pearls):