diff --git a/python/cuda_cccl/benchmarks/compute/bench_reduce.py b/python/cuda_cccl/benchmarks/compute/bench_reduce.py
index 46952a29e78..8bf483702a3 100644
--- a/python/cuda_cccl/benchmarks/compute/bench_reduce.py
+++ b/python/cuda_cccl/benchmarks/compute/bench_reduce.py
@@ -91,10 +91,7 @@ def run():
         reduce_pointer(input_array, build_only=(bench_fixture == "compile_benchmark"))
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -109,10 +106,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -127,10 +121,7 @@ def run():
         reduce_struct(input_array, build_only=(bench_fixture == "compile_benchmark"))
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -141,10 +132,7 @@ def run():
         reduce_pointer_custom_op(input_array, build_only=False)
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 def bench_reduce_pointer_single_phase(benchmark, size):
diff --git a/python/cuda_cccl/benchmarks/compute/bench_scan.py b/python/cuda_cccl/benchmarks/compute/bench_scan.py
index 6cd6b7a62bf..95f4e892c06 100644
--- a/python/cuda_cccl/benchmarks/compute/bench_scan.py
+++ b/python/cuda_cccl/benchmarks/compute/bench_scan.py
@@ -108,13 +108,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        if scan_type == "exclusive":
-            fixture(cuda.compute.make_exclusive_scan, run)
-        else:
-            fixture(cuda.compute.make_inclusive_scan, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("scan_type", ["exclusive", "inclusive"])
@@ -145,13 +139,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        if scan_type == "exclusive":
-            fixture(cuda.compute.make_exclusive_scan, run)
-        else:
-            fixture(cuda.compute.make_inclusive_scan, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("scan_type", ["exclusive", "inclusive"])
@@ -171,13 +159,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        if scan_type == "exclusive":
-            fixture(cuda.compute.make_exclusive_scan, run)
-        else:
-            fixture(cuda.compute.make_inclusive_scan, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 def scan_pointer_single_phase(input_array, build_only, scan_type):
diff --git a/python/cuda_cccl/benchmarks/compute/bench_select.py b/python/cuda_cccl/benchmarks/compute/bench_select.py
new file mode 100644
index 00000000000..8d284c63a74
--- /dev/null
+++ b/python/cuda_cccl/benchmarks/compute/bench_select.py
@@ -0,0 +1,148 @@
+import cupy as cp
+import numpy as np
+import pytest
+
+import cuda.compute
+from cuda.compute import (
+    CacheModifiedInputIterator,
+    gpu_struct,
+)
+
+
+def select_pointer(inp, out, num_selected, build_only):
+    size = len(inp)
+
+    def even_op(x):
+        return x % 2 == 0
+
+    selector = cuda.compute.make_select(inp, out, num_selected, even_op)
+    if not build_only:
+        temp_bytes = selector(None, inp, out, num_selected, size)
+        temp_storage = cp.empty(temp_bytes, dtype=np.uint8)
+        selector(temp_storage, inp, out, num_selected, size)
+
+    cp.cuda.runtime.deviceSynchronize()
+
+
+def select_iterator(size, d_in, out, num_selected, build_only):
+    d_in_iter = CacheModifiedInputIterator(d_in, modifier="stream")
+
+    def less_than_50(x):
+        return x < 50
+
+    selector = cuda.compute.make_select(d_in_iter, out, num_selected, less_than_50)
+    if not build_only:
+        temp_bytes = selector(None, d_in_iter, out, num_selected, size)
+        temp_storage = cp.empty(temp_bytes, dtype=np.uint8)
+        selector(temp_storage, d_in_iter, out, num_selected, size)
+
+    cp.cuda.runtime.deviceSynchronize()
+
+
+@gpu_struct
+class Point:
+    x: np.int32
+    y: np.int32
+
+
+def select_struct(inp, out, num_selected, build_only):
+    size = len(inp)
+
+    def in_first_quadrant(p: Point) -> np.uint8:
+        return (p.x > 50) and (p.y > 50)
+
+    selector = cuda.compute.make_select(inp, out, num_selected, in_first_quadrant)
+    if not build_only:
+        temp_bytes = selector(None, inp, out, num_selected, size)
+        temp_storage = cp.empty(temp_bytes, dtype=np.uint8)
+        selector(temp_storage, inp, out, num_selected, size)
+
+    cp.cuda.runtime.deviceSynchronize()
+
+
+def select_stateful(inp, out, num_selected, threshold_state, build_only):
+    size = len(inp)
+
+    def threshold_select(x):
+        return x > threshold_state[0]
+
+    selector = cuda.compute.make_select(inp, out, num_selected, threshold_select)
+    if not build_only:
+        temp_bytes = selector(None, inp, out, num_selected, size)
+        temp_storage = cp.empty(temp_bytes, dtype=np.uint8)
+        selector(temp_storage, inp, out, num_selected, size)
+
+    cp.cuda.runtime.deviceSynchronize()
+
+
+@pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
+def bench_select_pointer(bench_fixture, request, size):
+    actual_size = 100 if bench_fixture == "compile_benchmark" else size
+    inp = cp.random.randint(0, 100, actual_size, dtype=np.int32)
+    out = cp.empty_like(inp)
+    num_selected = cp.empty(2, dtype=np.uint64)
+
+    def run():
+        select_pointer(
+            inp, out, num_selected, build_only=(bench_fixture == "compile_benchmark")
+        )
+
+    fixture = request.getfixturevalue(bench_fixture)
+    fixture(run)
+
+
+@pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
+def bench_select_iterator(bench_fixture, request, size):
+    actual_size = 100 if bench_fixture == "compile_benchmark" else size
+    d_in = cp.random.randint(0, 100, actual_size, dtype=np.int32)
+    out = cp.empty(actual_size, dtype=np.int32)
+    num_selected = cp.empty(2, dtype=np.uint64)
+
+    def run():
+        select_iterator(
+            actual_size,
+            d_in,
+            out,
+            num_selected,
+            build_only=(bench_fixture == "compile_benchmark"),
+        )
+
+    fixture = request.getfixturevalue(bench_fixture)
+    fixture(run)
+
+
+@pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
+def bench_select_struct(bench_fixture, request, size):
+    actual_size = 100 if bench_fixture == "compile_benchmark" else size
+    inp = cp.random.randint(0, 100, (actual_size, 2), dtype=np.int32).view(Point.dtype)
+    out = cp.empty_like(inp)
+    num_selected = cp.empty(2, dtype=np.uint64)
+
+    def run():
+        select_struct(
+            inp, out, num_selected, build_only=(bench_fixture == "compile_benchmark")
+        )
+
+    fixture = request.getfixturevalue(bench_fixture)
+    fixture(run)
+
+
+@pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
+def bench_select_stateful(bench_fixture, request, size):
+    actual_size = 100 if bench_fixture == "compile_benchmark" else size
+    inp = cp.random.randint(0, 100, actual_size, dtype=np.int32)
+    out = cp.empty_like(inp)
+    num_selected = cp.empty(2, dtype=np.uint64)
+    threshold_state = cp.array([50], dtype=np.int32)
+
+    def run():
+        select_stateful(
+            inp,
+            out,
+            num_selected,
+            threshold_state,
+            build_only=(bench_fixture == "compile_benchmark"),
+        )
+
+    fixture = request.getfixturevalue(bench_fixture)
+    fixture(run)
diff --git a/python/cuda_cccl/benchmarks/compute/bench_three_way_partition.py b/python/cuda_cccl/benchmarks/compute/bench_three_way_partition.py
index e88b377f138..a9bef4e27a7 100644
--- a/python/cuda_cccl/benchmarks/compute/bench_three_way_partition.py
+++ b/python/cuda_cccl/benchmarks/compute/bench_three_way_partition.py
@@ -167,10 +167,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_three_way_partition, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -192,10 +189,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_three_way_partition, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -218,10 +212,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_three_way_partition, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 def three_way_partition_pointer_single_phase(inp):
diff --git a/python/cuda_cccl/benchmarks/compute/bench_transform.py b/python/cuda_cccl/benchmarks/compute/bench_transform.py
index e17044afd00..79c1cce4e77 100644
--- a/python/cuda_cccl/benchmarks/compute/bench_transform.py
+++ b/python/cuda_cccl/benchmarks/compute/bench_transform.py
@@ -117,10 +117,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_unary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -135,10 +132,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_unary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -154,10 +148,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_unary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -174,10 +165,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_binary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -192,10 +180,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_binary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -212,10 +197,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_binary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -231,7 +213,4 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_binary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
diff --git a/python/cuda_cccl/benchmarks/compute/bench_zip_iterator.py b/python/cuda_cccl/benchmarks/compute/bench_zip_iterator.py
index 41387086a9e..c9e27dcd5ae 100644
--- a/python/cuda_cccl/benchmarks/compute/bench_zip_iterator.py
+++ b/python/cuda_cccl/benchmarks/compute/bench_zip_iterator.py
@@ -106,10 +106,7 @@ def run():
         reduce_zip_array(input_array, build_only=(bench_fixture == "compile_benchmark"))
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -122,10 +119,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -138,10 +132,7 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_reduce_into, run)
-    else:
-        fixture(run)
+    fixture(run)
 
 
 @pytest.mark.parametrize("bench_fixture", ["compile_benchmark", "benchmark"])
@@ -163,7 +154,4 @@ def run():
         )
 
     fixture = request.getfixturevalue(bench_fixture)
-    if bench_fixture == "compile_benchmark":
-        fixture(cuda.compute.make_binary_transform, run)
-    else:
-        fixture(run)
+    fixture(run)
diff --git a/python/cuda_cccl/benchmarks/compute/conftest.py b/python/cuda_cccl/benchmarks/compute/conftest.py
index 0e4f9c73829..981b203b7b4 100644
--- a/python/cuda_cccl/benchmarks/compute/conftest.py
+++ b/python/cuda_cccl/benchmarks/compute/conftest.py
@@ -1,5 +1,7 @@
 import pytest
 
+import cuda.compute
+
 
 @pytest.fixture(params=[True, False])
 def build_only(request):
@@ -13,11 +15,11 @@ def size(request):
 
 @pytest.fixture
 def compile_benchmark(benchmark):
-    def run_compile_benchmark(algorithm, function):
+    def run_compile_benchmark(function):
         def setup():
             # This function is called once before the benchmark runs
             # to set up the environment.
-            algorithm.cache_clear()
+            cuda.compute.clear_all_caches()
 
         benchmark.pedantic(
             function,
diff --git a/python/cuda_cccl/cuda/compute/__init__.py b/python/cuda_cccl/cuda/compute/__init__.py
index 71f0ad70f4b..0ed7da29e34 100644
--- a/python/cuda_cccl/cuda/compute/__init__.py
+++ b/python/cuda_cccl/cuda/compute/__init__.py
@@ -2,6 +2,7 @@
 #
 # SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
 
+from ._caching import clear_all_caches
 from .algorithms import (
     DoubleBuffer,
     SortOrder,
@@ -49,13 +50,13 @@
 
 __all__ = [
     "binary_transform",
+    "clear_all_caches",
     "CacheModifiedInputIterator",
     "ConstantIterator",
     "CountingIterator",
     "DiscardIterator",
     "DoubleBuffer",
     "exclusive_scan",
-    "select",
     "gpu_struct",
     "histogram_even",
     "inclusive_scan",
@@ -81,6 +82,7 @@
     "ReverseIterator",
     "segmented_reduce",
     "segmented_sort",
+    "select",
     "SortOrder",
     "TransformIterator",
     "TransformOutputIterator",
@@ -89,3 +91,7 @@
     "unique_by_key",
     "ZipIterator",
 ]
+
+# import extensions to numba CUDA enabling global capture of
+# device arrays.
+from . import _device_array_capture  # noqa: F403, F401
diff --git a/python/cuda_cccl/cuda/compute/_caching.py b/python/cuda_cccl/cuda/compute/_caching.py
index 394125017a9..dfc55d1bdac 100644
--- a/python/cuda_cccl/cuda/compute/_caching.py
+++ b/python/cuda_cccl/cuda/compute/_caching.py
@@ -7,6 +7,9 @@
 
 from cuda.core.experimental import Device
 
+# Central registry of all algorithm caches
+_cache_registry: dict[str, object] = {}
+
 
 def cache_with_key(key):
     """
@@ -18,6 +21,9 @@ def cache_with_key(key):
     -----
     The CUDA compute capability of the current device is appended to
     the cache key returned by `key`.
+
+    The decorated function is automatically registered in the central
+    cache registry for easy cache management.
     """
 
     def deco(func):
@@ -36,11 +42,60 @@ def cache_clear():
             cache.clear()
 
         inner.cache_clear = cache_clear
+
+        # Register the cache in the central registry
+        cache_name = func.__qualname__
+        _cache_registry[cache_name] = inner
+
         return inner
 
     return deco
 
 
+def clear_all_caches():
+    """
+    Clear all algorithm caches.
+
+    This function clears all cached algorithm build results, forcing
+    recompilation on the next invocation. Useful for benchmarking
+    compilation time.
+
+    Example
+    -------
+    >>> import cuda.compute
+    >>> cuda.compute.clear_all_caches()
+    """
+    for cached_func in _cache_registry.values():
+        cached_func.cache_clear()
+
+
+def _hash_device_array_like(value):
+    # hash based on pointer, shape, and dtype
+    ptr = value.__cuda_array_interface__["data"][0]
+    shape = value.__cuda_array_interface__["shape"]
+    dtype = value.__cuda_array_interface__["typestr"]
+    return hash((ptr, shape, dtype))
+
+
+def _make_hashable(value):
+    import numba.cuda.dispatcher
+
+    from .typing import DeviceArrayLike
+
+    if isinstance(value, numba.cuda.dispatcher.CUDADispatcher):
+        return CachableFunction(value.py_func)
+    elif isinstance(value, DeviceArrayLike):
+        return _hash_device_array_like(value)
+    elif isinstance(value, (list, tuple)):
+        return tuple(_make_hashable(v) for v in value)
+    elif isinstance(value, dict):
+        return tuple(
+            sorted((_make_hashable(k), _make_hashable(v)) for k, v in value.items())
+        )
+    else:
+        return id(value)
+
+
 class CachableFunction:
     """
     A type that wraps a function and provides custom comparison
@@ -52,8 +107,6 @@ class CachableFunction:
     """
 
     def __init__(self, func):
-        import numba.cuda.dispatcher
-
         self._func = func
 
         closure = func.__closure__ if func.__closure__ is not None else []
@@ -61,16 +114,16 @@ def __init__(self, func):
         # if any of the contents is a numba.cuda.dispatcher.CUDADispatcher
         # use the function for caching purposes:
         for cell in closure:
-            if isinstance(cell.cell_contents, numba.cuda.dispatcher.CUDADispatcher):
-                contents.append(CachableFunction(cell.cell_contents.py_func))
-            else:
-                contents.append(cell.cell_contents)
+            contents.append(_make_hashable(cell.cell_contents))
         self._identity = (
             func.__name__,
             func.__code__.co_code,
             func.__code__.co_consts,
             tuple(contents),
-            tuple(func.__globals__.get(name, None) for name in func.__code__.co_names),
+            tuple(
+                _make_hashable(func.__globals__.get(name, None))
+                for name in func.__code__.co_names
+            ),
         )
 
     def __eq__(self, other):
diff --git a/python/cuda_cccl/cuda/compute/_device_array_capture.py b/python/cuda_cccl/cuda/compute/_device_array_capture.py
new file mode 100644
index 00000000000..321b92e8f70
--- /dev/null
+++ b/python/cuda_cccl/cuda/compute/_device_array_capture.py
@@ -0,0 +1,198 @@
+"""
+Extends Numba CUDA to enable global capture of objects implementing
+__cuda_array_interface__.
+
+HOW IT WORKS:
+1. Patches typeof_impl to check for DeviceArrayLike objects
+2. Defines a DeviceArrayType (subclass of Array) to represent DeviceArrayLike objects
+3. Registers data model and lower_constant for the type
+"""
+
+import numpy as np
+from numba.cuda import types
+from numba.cuda.core.imputils import builtin_registry
+from numba.cuda.datamodel.models import ArrayModel
+from numba.cuda.datamodel.registry import register_default
+from numba.cuda.np import numpy_support
+from numba.cuda.typing import typeof as typeof_module
+
+from .typing import DeviceArrayLike
+
+# =============================================================================
+# Step 1: Define a custom Array subclass for DeviceArrayLike objects
+# =============================================================================
+
+
+class DeviceArrayType(types.Array):
+    """
+    Type for any DeviceArrayLike object.
+
+    By subclassing Array and overriding copy() to ignore readonly,
+    we ensure these arrays remain mutable when captured from globals.
+    """
+
+    def __init__(self, dtype, ndim, layout, readonly=False, aligned=False):
+        type_name = "device_array"
+        if readonly:
+            type_name = "readonly " + type_name
+        if not aligned:
+            type_name = "unaligned " + type_name
+        name = f"{type_name}({dtype}, {ndim}d, {layout})"
+        super().__init__(
+            dtype, ndim, layout, readonly=readonly, name=name, aligned=aligned
+        )
+
+    def copy(self, dtype=None, ndim=None, layout=None, readonly=None):
+        # IGNORE readonly parameter - CAI objects are device memory references
+        # that should remain mutable
+        return DeviceArrayType(
+            dtype=dtype if dtype is not None else self.dtype,
+            ndim=ndim if ndim is not None else self.ndim,
+            layout=layout if layout is not None else self.layout,
+            readonly=not self.mutable,  # Preserve original mutability
+            aligned=self.aligned,
+        )
+
+
+# =============================================================================
+# Step 2: Define _typeof_device_array to type DeviceArrayLike objects
+# =============================================================================
+
+
+def _typeof_device_array(val, c):
+    """Type of a DeviceArrayLike object when captured as a constant."""
+
+    # don't affect typing of arguments:
+    if c.purpose == typeof_module.Purpose.argument:
+        return None
+
+    interface = val.__cuda_array_interface__
+
+    dtype = numpy_support.from_dtype(np.dtype(interface["typestr"]))
+    shape = interface["shape"]
+    ndim = len(shape)
+    strides = interface.get("strides")
+
+    # Determine layout
+    if ndim == 0:
+        layout = "C"
+    elif strides is None:
+        layout = "C"
+    else:
+        itemsize = np.dtype(interface["typestr"]).itemsize
+        # Check C-contiguous
+        c_strides = []
+        stride = itemsize
+        for i in range(ndim - 1, -1, -1):
+            c_strides.insert(0, stride)
+            stride *= shape[i]
+
+        if tuple(strides) == tuple(c_strides):
+            layout = "C"
+        else:
+            # Check F-contiguous
+            f_strides = []
+            stride = itemsize
+            for i in range(ndim):
+                f_strides.append(stride)
+                stride *= shape[i]
+            if tuple(strides) == tuple(f_strides):
+                layout = "F"
+            else:
+                layout = "A"
+
+    readonly = interface["data"][1]
+    return DeviceArrayType(dtype, ndim, layout, readonly=readonly)
+
+
+# =============================================================================
+# Step 3: Patch typeof_impl to handle CAI objects
+# =============================================================================
+
+# Get the original generic fallback
+_original_typeof_impl = typeof_module.typeof_impl
+_original_generic = _original_typeof_impl.dispatch(object)
+
+
+def _patched_generic(val, c):
+    """
+    Patched generic typeof handler that checks for DeviceArrayLike objects.
+    """
+    # Check for CAI objects that we should handle (constants only)
+    if isinstance(val, DeviceArrayLike):
+        result = _typeof_device_array(val, c)
+        if result is not None:
+            return result
+        # Fall through to original handler for arguments
+
+    # Fall back to original behavior
+    return _original_generic(val, c)
+
+
+# Register our patched handler for the generic object case
+_original_typeof_impl.register(object)(_patched_generic)
+
+
+# =============================================================================
+# Step 4: Register data model for our type
+# =============================================================================
+
+
+register_default(DeviceArrayType)(ArrayModel)
+# =============================================================================
+# Step 5: Register lower_constant for our type
+# =============================================================================
+
+
+@builtin_registry.lower_constant(DeviceArrayType)
+def lower_device_array(context, builder, ty, pyval):
+    """
+    Lower DeviceArrayLike objects by embedding the device pointer as a constant.
+    """
+
+    interface = pyval.__cuda_array_interface__
+
+    # hold on to the device-array-like object
+    lib = context.active_code_library
+    referenced_objects = getattr(lib, "referenced_objects", None)
+    if referenced_objects is None:
+        lib.referenced_objects = referenced_objects = {}
+    referenced_objects[id(pyval)] = pyval
+
+    shape = interface["shape"]
+    strides = interface.get("strides")
+    data_ptr = interface["data"][0]
+    typestr = interface["typestr"]
+
+    # Calculate strides if not provided (C-contiguous)
+    if strides is None:
+        itemsize = np.dtype(typestr).itemsize
+        ndim = len(shape)
+        strides = []
+        stride = itemsize
+        for i in range(ndim - 1, -1, -1):
+            strides.insert(0, stride)
+            stride *= shape[i]
+        strides = tuple(strides)
+
+    # Embed device pointer as constant
+    llvoidptr = context.get_value_type(types.voidptr)
+    data = context.get_constant(types.uintp, data_ptr).inttoptr(llvoidptr)
+
+    # Build array structure
+    ary = context.make_array(ty)(context, builder)
+    kshape = [context.get_constant(types.intp, s) for s in shape]
+    kstrides = [context.get_constant(types.intp, s) for s in strides]
+    itemsize = np.dtype(typestr).itemsize
+
+    context.populate_array(
+        ary,
+        data=builder.bitcast(data, ary.data.type),
+        shape=kshape,
+        strides=kstrides,
+        itemsize=context.get_constant(types.intp, itemsize),
+        parent=None,
+        meminfo=None,
+    )
+
+    return ary._getvalue()
diff --git a/python/cuda_cccl/tests/compute/test_select.py b/python/cuda_cccl/tests/compute/test_select.py
index b6b8f311f53..13b3128262e 100644
--- a/python/cuda_cccl/tests/compute/test_select.py
+++ b/python/cuda_cccl/tests/compute/test_select.py
@@ -420,3 +420,126 @@ def condition(pair):
     expected_count = np.sum(h_sums < 70)
 
     assert num_selected == expected_count
+
+
+def test_select_stateful_threshold():
+    """Test stateful select that uses state for threshold"""
+    num_items = 1000
+    h_in = random_array(num_items, np.int32, max_value=100)
+
+    # Create device state containing threshold value
+    threshold_value = 50
+    threshold_state = cp.array([threshold_value], dtype=np.int32)
+
+    # Define condition that references state as closure
+    def threshold_select(x):
+        return x > threshold_state[0]
+
+    d_in = cp.asarray(h_in)
+    d_out = cp.empty_like(d_in)
+    d_num_selected = cp.empty(2, dtype=np.uint64)
+
+    cuda.compute.select(
+        d_in,
+        d_out,
+        d_num_selected,
+        threshold_select,
+        num_items,
+    )
+
+    # Check selected output
+    num_selected = int(d_num_selected[0].get())
+    got = d_out.get()[:num_selected]
+
+    # Verify all output values are > threshold
+    assert np.all(got > threshold_value)
+
+    # Verify we got the expected number of items
+    expected_selected = h_in[h_in > threshold_value]
+    expected_count = len(expected_selected)
+
+    assert num_selected == expected_count
+
+    # Verify exact results
+    assert np.array_equal(got, expected_selected)
+
+
+def test_select_stateful_atomic():
+    """Test stateful select with atomic operations to count rejected items"""
+    from numba import cuda as numba_cuda
+
+    num_items = 1000
+    h_in = random_array(num_items, np.int32, max_value=100)
+
+    # Create device state for counting rejected items
+    reject_counter = cp.zeros(1, dtype=np.int32)
+
+    # Define condition that references state as closure
+    def count_rejects(x):
+        if x > 50:
+            return True
+        else:
+            numba_cuda.atomic.add(reject_counter, 0, 1)
+            return False
+
+    d_in = cp.asarray(h_in)
+    d_out = cp.empty_like(d_in)
+    d_num_selected = cp.empty(2, dtype=np.uint64)
+
+    cuda.compute.select(
+        d_in,
+        d_out,
+        d_num_selected,
+        count_rejects,
+        num_items,
+    )
+
+    # Check selected output
+    num_selected = int(d_num_selected[0].get())
+    got = d_out.get()[:num_selected]
+
+    # Verify all output values are > 50
+    assert np.all(got > 50)
+
+    # Verify we got the expected number of items
+    expected_selected = h_in[h_in > 50]
+    expected_count = len(expected_selected)
+
+    assert num_selected == expected_count
+
+    # Verify exact results
+    assert np.array_equal(got, expected_selected)
+
+    # Verify state contains count of rejected items
+    rejected_count = int(reject_counter[0].get())
+    expected_rejected = len(h_in[h_in <= 50])
+    assert rejected_count == expected_rejected, (
+        f"Expected {expected_rejected} rejections, got {rejected_count}"
+    )
+
+
+def test_select_stateful_counting_rejects():
+    """Test select that counts rejected items using state"""
+    from numba import cuda as numba_cuda
+
+    d_in = cp.arange(100, dtype=np.int32)
+    d_out = cp.empty_like(d_in)
+    d_num_selected = cp.empty(1, dtype=np.uint64)
+
+    reject_count = cp.zeros(1, dtype=np.int32)
+
+    # Define condition that references state as closure
+    def count_rejects(x):
+        if x >= 50:
+            return True
+        else:
+            numba_cuda.atomic.add(reject_count, 0, 1)
+            return False
+
+    cuda.compute.select(d_in, d_out, d_num_selected, count_rejects, len(d_in))
+
+    num_selected = int(d_num_selected.get()[0])
+    num_rejected = int(reject_count.get()[0])
+
+    assert num_selected == 50  # Values 50-99
+    assert num_rejected == 50  # Values 0-49
diff --git a/python/cuda_cccl/tests/compute/test_transform.py b/python/cuda_cccl/tests/compute/test_transform.py
index dc3aceabcb3..a4712e46249 100644
--- a/python/cuda_cccl/tests/compute/test_transform.py
+++ b/python/cuda_cccl/tests/compute/test_transform.py
@@ -380,3 +380,27 @@ def add_vectors(v1: Vec2D, v2: Vec2D) -> Vec2D:
 
     np.testing.assert_equal(result["x"], h_in1["x"] + h_in2["x"])
     np.testing.assert_equal(result["y"], h_in1["y"] + h_in2["y"])
+
+
+def test_unary_transform_stateful_counting():
+    """Test unary_transform with state that counts even numbers."""
+    from numba import cuda as numba_cuda
+
+    d_in = cp.arange(100, dtype=np.int32)
+    d_out = cp.empty_like(d_in)
+
+    even_count = cp.zeros(1, dtype=np.int32)
+
+    # Define op that references state as closure
+    def count_evens(x):
+        if x % 2 == 0:
+            numba_cuda.atomic.add(even_count, 0, 1)
+        return x * 2
+
+    cuda.compute.unary_transform(d_in, d_out, count_evens, len(d_in))
+
+    expected_output = cp.arange(100, dtype=np.int32) * 2
+    np.testing.assert_array_equal(d_out.get(), expected_output.get())
+
+    num_evens = int(even_count.get()[0])
+    assert num_evens == 50  # 0, 2, 4, ..., 98