Improve the ring backend initialization

.circleci/config.yml

@@ -7,6 +7,18 @@ parameters:
   nightly_build:
     type: boolean
     default: false
+  weekly_build:
+    type: boolean
+    default: false
+  test_release:
+    type: boolean
+    default: false
+  linux_release:
+    type: boolean
+    default: false
+  cuda_release:
+    type: boolean
+    default: false
 
 jobs:
   build_documentation:
@@ -61,9 +73,9 @@ jobs:
                  git push -f origin gh-pages
 
   linux_build_and_test:
-    machine:
-      image: ubuntu-2204:current
-      resource_class: large
+    docker:
+      - image: cimg/python:3.9
+
     steps:
       - checkout
       - run:
@@ -75,17 +87,19 @@ jobs:
       - run:
           name: Install dependencies
           command: |
-            export DEBIAN_FRONTEND=noninteractive
-            export NEEDRESTART_MODE=a
-            sudo apt-get update
-            sudo apt-get upgrade -y
             pip install --upgrade cmake
-            sudo apt-get install -y libblas-dev liblapack-dev liblapacke-dev
+            pip install nanobind==2.4.0
+            pip install numpy
+            sudo apt-get update
+            sudo apt-get install libblas-dev liblapack-dev liblapacke-dev
             sudo apt-get install openmpi-bin openmpi-common libopenmpi-dev
       - run:
           name: Install Python package
           command: |
-            pip install -e ".[dev]"
+            CMAKE_ARGS="-DMLX_BUILD_METAL=OFF" \
+              python3 setup.py build_ext --inplace
+            CMAKE_ARGS="-DMLX_BUILD_METAL=OFF" \
+              python3 setup.py develop
       - run:
           name: Generate package stubs
           command: |
@@ -95,14 +109,13 @@ jobs:
       - run:
           name: Run Python tests
           command: |
-            python -m unittest discover python/tests -v
+            python3 -m unittest discover python/tests -v
             mpirun --bind-to none -host localhost:8 -np 8 python python/tests/mpi_test_distributed.py
-            mlx.launch --verbose -n 8 python/tests/ring_test_distributed.py -v 2> >(tee -a stderr.log >&2)
-            if $(grep "\[WARN\]" stderr.log); then echo "Distributed ring test failed"; exit 1; fi
+            mlx.launch --verbose -n 8 python/tests/ring_test_distributed.py
       - run:
           name: Build CPP only
           command: |
-            mkdir -p build && cd build
+            mkdir -p build && cd build 
             cmake .. -DMLX_BUILD_METAL=OFF -DCMAKE_BUILD_TYPE=DEBUG
             make -j `nproc`
       - run:
@@ -157,8 +170,7 @@ jobs:
             LOW_MEMORY=1 DEVICE=cpu python -m xmlrunner discover -v python/tests -o test-results/cpu
             LOW_MEMORY=1 DEVICE=gpu METAL_DEVICE_WRAPPER_TYPE=1 METAL_DEBUG_ERROR_MODE=0 python -m xmlrunner discover -v python/tests -o test-results/gpu
             mpirun --bind-to none -host localhost:8 -np 8 -x DYLD_LIBRARY_PATH=/opt/homebrew/lib/ python python/tests/mpi_test_distributed.py
-            mlx.launch --verbose -n 8 python/tests/ring_test_distributed.py -v 2> >(tee -a stderr.log >&2)
-            if $(grep "\[WARN\]" stderr.log); then echo "Distributed ring test failed"; exit 1; fi
+            mlx.launch --verbose -n 8 python/tests/ring_test_distributed.py
       - run:
           name: Build example extension
           command: |
@@ -201,7 +213,7 @@ jobs:
 
   cuda_build_and_test:
     machine:
-      image: linux-cuda-12:2023.11.1
+      image: linux-cuda-12:default
       resource_class: gpu.nvidia.small.gen2
     steps:
       - checkout
@@ -210,7 +222,7 @@ jobs:
           command: |
             sudo apt-get update
             sudo apt-get install libblas-dev liblapack-dev liblapacke-dev
-            python3 -m venv env
+            python -m venv env
             source env/bin/activate
             CMAKE_ARGS="-DMLX_BUILD_CUDA=ON -DCMAKE_CUDA_COMPILER=`which nvcc`" \
               pip install -e ".[dev]"
@@ -272,18 +284,7 @@ jobs:
           name: Build Python package
           command: |
             source env/bin/activate
-            python setup.py clean --all
-            << parameters.build_env >> MLX_BUILD_STAGE=1 python -m build -w
-      - when:
-          condition:
-            equal: ["3.9", << parameters.python_version >>]
-          steps:
-            - run:
-                name: Build common package
-                command: |
-                  source env/bin/activate
-                  python setup.py clean --all
-                  << parameters.build_env >> MLX_BUILD_STAGE=2 python -m build -w
+            << parameters.build_env >> python -m build -w
       - when:
           condition: << parameters.build_env >>
           steps:
@@ -300,71 +301,59 @@ jobs:
       python_version:
         type: string
         default: "3.9"
-      build_env:
+      extra_env:
         type: string
-        default: ""
-    machine:
-      image: ubuntu-2204:current
-      resource_class: large
+        default: "DEV_RELEASE=1"
+    docker:
+      - image: ubuntu:20.04
     steps:
       - checkout
       - run:
           name: Build wheel
           command: |
             PYTHON=python<< parameters.python_version >>
-            export DEBIAN_FRONTEND=noninteractive
-            export NEEDRESTART_MODE=a
-            sudo apt-get update
-            sudo apt-get upgrade -y
-            TZ=Etc/UTC sudo apt-get -y install tzdata
-            sudo apt-get install -y apt-utils
-            sudo apt-get install -y software-properties-common
-            sudo add-apt-repository -y ppa:deadsnakes/ppa
-            sudo apt-get install -y $PYTHON $PYTHON-dev $PYTHON-full
-            sudo apt-get install -y libblas-dev liblapack-dev liblapacke-dev
-            sudo apt-get install -y build-essential git
+            apt-get update
+            apt-get upgrade -y
+            DEBIAN_FRONTEND=noninteractive TZ=Etc/UTC apt-get -y install tzdata
+            apt-get install -y apt-utils
+            apt-get install -y software-properties-common
+            add-apt-repository -y ppa:deadsnakes/ppa
+            apt-get install -y $PYTHON $PYTHON-dev $PYTHON-full
+            apt-get install -y libblas-dev liblapack-dev liblapacke-dev
+            apt-get install -y build-essential git
             $PYTHON -m venv env
             source env/bin/activate
             pip install --upgrade pip
             pip install --upgrade cmake
+            pip install nanobind==2.4.0
+            pip install --upgrade setuptools
+            pip install numpy
             pip install auditwheel
             pip install patchelf
             pip install build
             pip install twine
-            << parameters.build_env >> pip install ".[dev]" -v
+            << parameters.extra_env >> pip install . -v
             pip install typing_extensions
             python setup.py generate_stubs
-            python setup.py clean --all
-            MLX_BUILD_STAGE=1 << parameters.build_env >> python -m build -w
-            bash python/scripts/repair_linux.sh
-      - when:
-          condition:
-            equal: ["3.9", << parameters.python_version >>]
-          steps:
-            - run:
-                name: Build common package
-                command: |
-                  source env/bin/activate
-                  python setup.py clean --all
-                  << parameters.build_env >> MLX_BUILD_STAGE=2 \
-                    python -m build -w
-                  auditwheel repair dist/mlx_cpu*.whl --plat manylinux_2_35_x86_64
-      - when:
-          condition: << parameters.build_env >>
-          steps:
-            - run:
-                name: Upload packages
-                command: |
-                  source env/bin/activate
-                  twine upload wheelhouse/*.whl
+            << parameters.extra_env >> python -m build --wheel
+            auditwheel show dist/*
+            auditwheel repair dist/* --plat manylinux_2_31_x86_64
+      - run:
+          name: Upload package
+          command: |
+            source env/bin/activate
+            twine upload wheelhouse/*
       - store_artifacts:
           path: wheelhouse/
 
   build_cuda_release:
     parameters:
-      build_env:
+      python_version:
         type: string
-        default: ""
+        default: "3.9"
+      extra_env:
+        type: string
+        default: "DEV_RELEASE=1"
     machine:
       image: linux-cuda-12:default
       resource_class: gpu.nvidia.small.gen2
@@ -375,25 +364,25 @@ jobs:
           command: |
             sudo apt-get update
             sudo apt-get install libblas-dev liblapack-dev liblapacke-dev
-            sudo apt-get install zip
             python -m venv env
             source env/bin/activate
             pip install auditwheel
             pip install patchelf
             pip install build
             pip install twine
-            << parameters.build_env >> MLX_BUILD_STAGE=2 \
+            << parameters.extra_env >> \
               CMAKE_ARGS="-DMLX_BUILD_CUDA=ON -DCMAKE_CUDA_COMPILER=`which nvcc`" \
-              python -m build -w
+              pip install ".[dev]" -v
+            python setup.py generate_stubs
+            << parameters.extra_env >> \
+              CMAKE_ARGS="-DMLX_BUILD_CUDA=ON -DCMAKE_CUDA_COMPILER=`which nvcc`" \
+              python -m build --wheel
             bash python/scripts/repair_cuda.sh
-      - when:
-          condition: << parameters.build_env >>
-          steps:
-            - run:
-                name: Upload package
-                command: |
-                  source env/bin/activate
-                  twine upload wheelhouse/*.whl
+      - run:
+          name: Upload package
+          command: |
+            source env/bin/activate
+            twine upload wheelhouse/*.whl
       - store_artifacts:
           path: wheelhouse/
 
@@ -405,6 +394,8 @@ workflows:
             pattern: "^(?!pull/)[-\\w]+$"
             value: << pipeline.git.branch >>
         - not: << pipeline.parameters.nightly_build >>
+        - not: << pipeline.parameters.weekly_build >>
+        - not: << pipeline.parameters.test_release >>
     jobs:
       - mac_build_and_test:
           matrix:
@@ -418,6 +409,8 @@ workflows:
     when:
       and:
         - not: << pipeline.parameters.nightly_build >>
+        - not: << pipeline.parameters.weekly_build >>
+        - not: << pipeline.parameters.test_release >>
     jobs:
       - build_release:
           filters:
@@ -508,16 +501,7 @@ workflows:
           matrix:
             parameters:
               python_version: ["3.9", "3.10", "3.11", "3.12", "3.13"]
-              build_env: ["PYPI_RELEASE=1"]
-      - build_cuda_release:
-          filters:
-            tags:
-              only: /^v.*/
-            branches:
-              ignore: /.*/
-          matrix:
-            parameters:
-              build_env: ["PYPI_RELEASE=1"]
+              extra_env: ["PYPI_RELEASE=1"]
 
   prb:
     when:
@@ -596,8 +580,99 @@ workflows:
               - macosx_deployment_target: "15.0"
                 xcode_version: "15.0.0"
                 python_version: "3.13"
+  weekly_build:
+    when:
+      and:
+        - equal: [ main, << pipeline.git.branch >> ]
+        - << pipeline.parameters.weekly_build >>
+    jobs:
+      - build_release:
+          matrix:
+            parameters:
+              python_version: ["3.9", "3.10", "3.11", "3.12", "3.13"]
+              macosx_deployment_target: ["13.5", "14.0", "15.0"]
+              build_env: ["DEV_RELEASE=1"]
+              xcode_version: ["16.2.0", "15.0.0"]
+            exclude:
+              - macosx_deployment_target: "13.5"
+                xcode_version: "16.2.0"
+                python_version: "3.9"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "13.5"
+                xcode_version: "16.2.0"
+                python_version: "3.10"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "13.5"
+                xcode_version: "16.2.0"
+                python_version: "3.11"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "13.5"
+                xcode_version: "16.2.0"
+                python_version: "3.12"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "13.5"
+                xcode_version: "16.2.0"
+                python_version: "3.13"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "14.0"
+                xcode_version: "15.0.0"
+                python_version: "3.9"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "14.0"
+                xcode_version: "15.0.0"
+                python_version: "3.10"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "14.0"
+                xcode_version: "15.0.0"
+                python_version: "3.11"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "14.0"
+                xcode_version: "15.0.0"
+                python_version: "3.12"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "14.0"
+                xcode_version: "15.0.0"
+                python_version: "3.13"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "15.0"
+                xcode_version: "15.0.0"
+                python_version: "3.9"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "15.0"
+                xcode_version: "15.0.0"
+                python_version: "3.10"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "15.0"
+                xcode_version: "15.0.0"
+                python_version: "3.11"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "15.0"
+                xcode_version: "15.0.0"
+                python_version: "3.12"
+                build_env: "DEV_RELEASE=1"
+              - macosx_deployment_target: "15.0"
+                xcode_version: "15.0.0"
+                python_version: "3.13"
+                build_env: "DEV_RELEASE=1"
+  linux_test_release:
+    when:
+      and:
+        - equal: [ main, << pipeline.git.branch >> ]
+        - << pipeline.parameters.linux_release >>
+    jobs:
       - build_linux_release:
           matrix:
             parameters:
               python_version: ["3.9", "3.10", "3.11", "3.12", "3.13"]
-      - build_cuda_release
+              extra_env: ["PYPI_RELEASE=1"]
+  cuda_test_release:
+    when:
+      and:
+        - equal: [ main, << pipeline.git.branch >> ]
+        - << pipeline.parameters.cuda_release >>
+    jobs:
+      - build_cuda_release:
+          matrix:
+            parameters:
+              python_version: ["3.9", "3.10", "3.11", "3.12", "3.13"]
+              extra_env: ["PYPI_RELEASE=1"]

ACKNOWLEDGMENTS.md

@@ -19,7 +19,6 @@ MLX was developed with contributions from the following individuals:
 - Gleb Pobudzey: Added the `where` primitive, and groups in 1D and 2D convolutions.
 - Paul Paczuski: Improved stability of BCE loss calculation
 - Max-Heinrich Laves: Added `conv_transpose1d`, `conv_transpose2d`, and `conv_transpose3d` ops.
-- Gökdeniz Gülmez: Added the `Muon (MomentUm Orthogonalized by Newton-schulz)` optimizer.
 
 <a href="https://github.com/ml-explore/mlx/graphs/contributors">
   <img class="dark-light" src="https://contrib.rocks/image?repo=ml-explore/mlx&anon=0&columns=20&max=100&r=true" />

CMakeLists.txt

@@ -64,8 +64,10 @@ if(${CMAKE_SYSTEM_NAME} MATCHES "Darwin")
       message(WARNING "Building for x86_64 arch is not officially supported.")
     endif()
   endif()
+
 else()
   set(MLX_BUILD_METAL OFF)
+  message(WARNING "MLX is prioritised for Apple silicon systems using macOS.")
 endif()
 
 # ----------------------------- Lib -----------------------------

docs/src/dev/extensions.rst

@@ -138,13 +138,13 @@ more concrete:
         * representing the vectorized computation and the axis which
         * corresponds to the output vectorized dimension.
         */
-        std::pair<std::vector<array>, std::vector<int>> vmap(
+        virtual std::pair<std::vector<array>, std::vector<int>> vmap(
             const std::vector<array>& inputs,
             const std::vector<int>& axes) override;
 
-        /** The name of primitive. */
-        const char* name() const override {
-          return "Axpby";
+        /** Print the primitive. */
+        void print(std::ostream& os) override {
+            os << "Axpby";
         }
 
         /** Equivalence check **/

docs/src/install.rst

@@ -23,6 +23,13 @@ To install from PyPI you must meet the following requirements:
     MLX is only available on devices running macOS >= 13.5
     It is highly recommended to use macOS 14 (Sonoma)
 
+
+MLX is also available on conda-forge. To install MLX with conda do:
+
+.. code-block:: shell
+
+   conda install conda-forge::mlx
+
 CUDA
 ^^^^
 
@@ -31,16 +38,8 @@ and SM 7.0 (Volta) and up. To install MLX with CUDA support, run:
 
 .. code-block:: shell
 
-    pip install "mlx[cuda]"
+    pip install mlx-cuda
 
-CPU-only (Linux)
-^^^^^^^^^^^^^^^^
-
-For a CPU-only version of MLX that runs on Linux use:
-
-.. code-block:: shell
-
-    pip install "mlx[cpu]"
 
 Troubleshooting
 ^^^^^^^^^^^^^^^

docs/src/python/optimizers/common_optimizers.rst

@@ -19,4 +19,3 @@ Common Optimizers
    Adamax
    Lion
    MultiOptimizer
-   Muon

examples/extensions/axpby/axpby.h

@@ -74,9 +74,9 @@ class Axpby : public mx::Primitive {
       const std::vector<mx::array>& inputs,
       const std::vector<int>& axes) override;
 
-  /** The name of primitive. */
-  const char* name() const override {
-    return "Axpby";
+  /** Print the primitive. */
+  void print(std::ostream& os) override {
+    os << "Axpby";
   }
 
   /** Equivalence check **/

mlx/backend/common/utils.cpp

@@ -1,20 +1,14 @@
 // Copyright © 2023-2024 Apple Inc.
 
-#include <dlfcn.h>
-
 #include "mlx/backend/common/utils.h"
+#include "mlx/primitives.h"
 
 namespace mlx::core {
 
-std::filesystem::path current_binary_dir() {
-  static std::filesystem::path binary_dir = []() {
-    Dl_info info;
-    if (!dladdr(reinterpret_cast<void*>(&current_binary_dir), &info)) {
-      throw std::runtime_error("Unable to get current binary dir.");
-    }
-    return std::filesystem::path(info.dli_fname).parent_path();
-  }();
-  return binary_dir;
+std::string get_primitive_string(Primitive* primitive) {
+  std::ostringstream op_t;
+  primitive->print(op_t);
+  return op_t.str();
 }
 
 std::tuple<Shape, std::vector<Strides>> collapse_contiguous_dims(

mlx/backend/common/utils.h

@@ -2,7 +2,6 @@
 
 #pragma once
 
-#include <filesystem>
 #include <tuple>
 #include <vector>
 
@@ -10,8 +9,7 @@
 
 namespace mlx::core {
 
-// Return the directory that contains current shared library.
-std::filesystem::path current_binary_dir();
+std::string get_primitive_string(Primitive* primitive);
 
 inline int64_t
 elem_to_loc(int elem, const Shape& shape, const Strides& strides) {

mlx/backend/cpu/cholesky.cpp

@@ -20,7 +20,7 @@ void cholesky_impl(const array& a, array& factor, bool upper, Stream stream) {
 
   // The decomposition is computed in place, so just copy the input to the
   // output.
-  copy_cpu(
+  copy(
       a,
       factor,
       a.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cpu/compiled.cpp

@@ -231,7 +231,7 @@ inline void build_kernel(
       os << "static_cast<" << get_type_string(x.dtype()) << ">(tmp_"
          << namer.get_name(x.inputs()[0]) << ");" << std::endl;
     } else {
-      os << x.primitive().name();
+      x.primitive().print(os);
       os << "()(";
       for (int i = 0; i < x.inputs().size() - 1; i++) {
         os << "tmp_" << namer.get_name(x.inputs()[i]) << ", ";

mlx/backend/cpu/conv.cpp

@@ -883,7 +883,7 @@ void explicit_gemm_conv_1D_cpu(
   // Fill with zeros
   std::vector<array> temps;
   temps.push_back(array(0, conv_dtype));
-  copy_cpu(temps.back(), in_padded, CopyType::Scalar, stream);
+  copy(temps.back(), in_padded, CopyType::Scalar, stream);
 
   // Pick input slice from padded
   size_t data_offset = padding_lo[0] * in_padded.strides()[1];
@@ -895,7 +895,7 @@ void explicit_gemm_conv_1D_cpu(
       in_padded_slice.size(),
       data_offset);
   // Copy input values into the slice
-  copy_cpu_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
+  copy_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
   temps.push_back(in_padded_slice);
 
   // Make strided view
@@ -920,7 +920,7 @@ void explicit_gemm_conv_1D_cpu(
   // Materialize strided view
   Shape strided_reshape = {N * oH, wH * C};
   array in_strided(strided_reshape, in_strided_view.dtype(), nullptr, {});
-  copy_cpu(in_strided_view, in_strided, CopyType::General, stream);
+  copy(in_strided_view, in_strided, CopyType::General, stream);
   temps.push_back(in_strided);
 
   // Check wt dtype and prepare
@@ -938,13 +938,13 @@ void explicit_gemm_conv_1D_cpu(
         wt.size(),
         0);
     gemm_wt = array(wt_transpose.shape(), float32, nullptr, {});
-    copy_cpu(wt_transpose, gemm_wt, CopyType::General, stream);
+    copy(wt_transpose, gemm_wt, CopyType::General, stream);
     temps.push_back(gemm_wt);
   } else if (wt.dtype() != float32 || !wt.flags().row_contiguous) {
     auto ctype =
         wt.flags().row_contiguous ? CopyType::Vector : CopyType::General;
     gemm_wt = array(wt.shape(), float32, nullptr, {});
-    copy_cpu(wt, gemm_wt, ctype, stream);
+    copy(wt, gemm_wt, ctype, stream);
     temps.push_back(gemm_wt);
   }
 
@@ -991,7 +991,7 @@ void explicit_gemm_conv_1D_cpu(
 
   // Copy results if needed
   if (out.dtype() != float32) {
-    copy_cpu_inplace(gemm_out, out, CopyType::Vector, stream);
+    copy_inplace(gemm_out, out, CopyType::Vector, stream);
   }
   encoder.add_temporaries(std::move(temps));
 }
@@ -1029,7 +1029,7 @@ void explicit_gemm_conv_2D_cpu(
   // Fill with zeros
   std::vector<array> temps;
   temps.push_back(array(0, conv_dtype));
-  copy_cpu(temps.back(), in_padded, CopyType::Scalar, stream);
+  copy(temps.back(), in_padded, CopyType::Scalar, stream);
 
   // Pick input slice from padded
   size_t data_offset = padding_lo[0] * in_padded.strides()[1] +
@@ -1044,7 +1044,7 @@ void explicit_gemm_conv_2D_cpu(
   temps.push_back(in_padded_slice);
 
   // Copy input values into the slice
-  copy_cpu_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
+  copy_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
 
   // Make strided view
   Shape strided_shape = {N, oH, oW, wH, wW, C};
@@ -1065,7 +1065,7 @@ void explicit_gemm_conv_2D_cpu(
   // Materialize strided view
   Shape strided_reshape = {N * oH * oW, wH * wW * C};
   array in_strided(strided_reshape, in_strided_view.dtype(), nullptr, {});
-  copy_cpu(in_strided_view, in_strided, CopyType::General, stream);
+  copy(in_strided_view, in_strided, CopyType::General, stream);
   temps.push_back(in_strided);
 
   // Check wt dtype and prepare
@@ -1076,7 +1076,7 @@ void explicit_gemm_conv_2D_cpu(
     auto ctype =
         wt.flags().row_contiguous ? CopyType::Vector : CopyType::General;
     gemm_wt = array(wt.shape(), float32, nullptr, {});
-    copy_cpu(wt, gemm_wt, ctype, stream);
+    copy(wt, gemm_wt, ctype, stream);
     temps.push_back(gemm_wt);
   }
 
@@ -1116,7 +1116,7 @@ void explicit_gemm_conv_2D_cpu(
 
   // Copy results if needed
   if (out.dtype() != float32) {
-    copy_cpu_inplace(gemm_out, out, CopyType::Vector, stream);
+    copy_inplace(gemm_out, out, CopyType::Vector, stream);
   }
   encoder.add_temporaries(std::move(temps));
 }
@@ -1156,7 +1156,7 @@ void explicit_gemm_conv_ND_cpu(
 
   // Fill with zeros
   std::vector<array> temps = {array(0, conv_dtype)};
-  copy_cpu(temps.back(), in_padded, CopyType::Scalar, stream);
+  copy(temps.back(), in_padded, CopyType::Scalar, stream);
 
   // Pick input slice from padded
   size_t data_offset = 0;
@@ -1173,7 +1173,7 @@ void explicit_gemm_conv_ND_cpu(
       data_offset);
 
   // Copy input values into the slice
-  copy_cpu_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
+  copy_inplace(in, in_padded_slice, CopyType::GeneralGeneral, stream);
   temps.push_back(in_padded_slice);
 
   // Make strided view
@@ -1212,7 +1212,7 @@ void explicit_gemm_conv_ND_cpu(
   }
 
   array in_strided(strided_reshape, in_strided_view.dtype(), nullptr, {});
-  copy_cpu(in_strided_view, in_strided, CopyType::General, stream);
+  copy(in_strided_view, in_strided, CopyType::General, stream);
   temps.push_back(in_strided);
 
   // Check wt dtype and prepare
@@ -1223,13 +1223,13 @@ void explicit_gemm_conv_ND_cpu(
     auto ctype =
         wt.flags().row_contiguous ? CopyType::Vector : CopyType::General;
     gemm_wt = array(wt.shape(), float32, nullptr, {});
-    copy_cpu(wt, gemm_wt, ctype, stream);
+    copy(wt, gemm_wt, ctype, stream);
     temps.push_back(gemm_wt);
   }
 
   if (flip) {
     auto gemm_wt_ = array(gemm_wt.shape(), float32, nullptr, {});
-    copy_cpu(gemm_wt, gemm_wt_, CopyType::Vector, stream);
+    copy(gemm_wt, gemm_wt_, CopyType::Vector, stream);
     temps.push_back(gemm_wt_);
 
     // Calculate the total size of the spatial dimensions
@@ -1284,7 +1284,7 @@ void explicit_gemm_conv_ND_cpu(
 
   // Copy results if needed
   if (out.dtype() != float32) {
-    copy_cpu_inplace(gemm_out, out, CopyType::Vector, stream);
+    copy_inplace(gemm_out, out, CopyType::Vector, stream);
   }
   encoder.add_temporaries(std::move(temps));
 }

mlx/backend/cpu/copy.cpp

@@ -295,11 +295,7 @@ inline void copy_inplace_dispatch(
 
 } // namespace
 
-void copy_cpu_inplace(
-    const array& src,
-    array& dst,
-    CopyType ctype,
-    Stream stream) {
+void copy_inplace(const array& src, array& dst, CopyType ctype, Stream stream) {
   auto& encoder = cpu::get_command_encoder(stream);
   encoder.set_input_array(src);
   encoder.set_output_array(dst);
@@ -309,7 +305,7 @@ void copy_cpu_inplace(
        ctype]() mutable { copy_inplace_dispatch(src, dst, ctype); });
 }
 
-void copy_cpu(const array& src, array& dst, CopyType ctype, Stream stream) {
+void copy(const array& src, array& dst, CopyType ctype, Stream stream) {
   bool donated = set_copy_output_data(src, dst, ctype);
   if (donated && src.dtype() == dst.dtype()) {
     // If the output has the same type as the input then there is nothing to
@@ -319,10 +315,10 @@ void copy_cpu(const array& src, array& dst, CopyType ctype, Stream stream) {
   if (ctype == CopyType::GeneralGeneral) {
     ctype = CopyType::General;
   }
-  copy_cpu_inplace(src, dst, ctype, stream);
+  copy_inplace(src, dst, ctype, stream);
 }
 
-void copy_cpu_inplace(
+void copy_inplace(
     const array& src,
     array& dst,
     const Shape& data_shape,

mlx/backend/cpu/copy.h

@@ -10,14 +10,10 @@
 
 namespace mlx::core {
 
-void copy_cpu(const array& src, array& dst, CopyType ctype, Stream stream);
-void copy_cpu_inplace(
-    const array& src,
-    array& dst,
-    CopyType ctype,
-    Stream stream);
+void copy(const array& src, array& dst, CopyType ctype, Stream stream);
+void copy_inplace(const array& src, array& dst, CopyType ctype, Stream stream);
 
-void copy_cpu_inplace(
+void copy_inplace(
     const array& src,
     array& dst,
     const Shape& data_shape,

mlx/backend/cpu/distributed.cpp

@@ -14,7 +14,7 @@ std::pair<array, bool> ensure_row_contiguous(const array& arr, Stream stream) {
     return {arr, false};
   } else {
     array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-    copy_cpu(arr, arr_copy, CopyType::General, stream);
+    copy(arr, arr_copy, CopyType::General, stream);
     return {arr_copy, true};
   }
 };
@@ -35,7 +35,7 @@ void AllReduce::eval_cpu(
       return in;
     } else {
       array arr_copy(in.shape(), in.dtype(), nullptr, {});
-      copy_cpu(in, arr_copy, CopyType::General, s);
+      copy(in, arr_copy, CopyType::General, s);
       out.copy_shared_buffer(arr_copy);
       return arr_copy;
     }

mlx/backend/cpu/eig.cpp

@@ -135,7 +135,7 @@ void Eig::eval_cpu(
       : array(a.shape(), complex64, nullptr, {});
 
   auto a_copy = array(a.shape(), a.dtype(), nullptr, {});
-  copy_cpu(
+  copy(
       a,
       a_copy,
       a.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cpu/eigh.cpp

@@ -196,7 +196,7 @@ void Eigh::eval_cpu(
 
   values.set_data(allocator::malloc(values.nbytes()));
 
-  copy_cpu(
+  copy(
       a,
       vectors,
       a.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cpu/hadamard.cpp

@@ -96,7 +96,7 @@ void Hadamard::eval_cpu(const std::vector<array>& inputs, array& out) {
   if (in.flags().row_contiguous && in.is_donatable()) {
     out.copy_shared_buffer(in);
   } else {
-    copy_cpu(
+    copy(
         in,
         out,
         in.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cpu/indexing.cpp

@@ -517,7 +517,7 @@ void Scatter::eval_cpu(const std::vector<array>& inputs, array& out) {
   // Copy src into out (copy allocates memory for out)
   auto ctype =
       src.flags().row_contiguous ? CopyType::Vector : CopyType::General;
-  copy_cpu(src, out, ctype, stream());
+  copy(src, out, ctype, stream());
 
   auto& encoder = cpu::get_command_encoder(stream());
   std::vector<array> inds;
@@ -686,7 +686,7 @@ void ScatterAxis::eval_cpu(const std::vector<array>& inputs, array& out) {
   // Copy src into out (copy allocates memory for out)
   auto ctype =
       src.flags().row_contiguous ? CopyType::Vector : CopyType::General;
-  copy_cpu(src, out, ctype, stream());
+  copy(src, out, ctype, stream());
 
   auto& encoder = cpu::get_command_encoder(stream());
   encoder.set_input_array(idx);

mlx/backend/cpu/inverse.cpp

@@ -115,7 +115,7 @@ void inverse_impl(
   //   (A⁻¹)ᵀ = (Aᵀ)⁻¹
 
   // The inverse is computed in place, so just copy the input to the output.
-  copy_cpu(
+  copy(
       a,
       inv,
       a.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cpu/logsumexp.cpp

@@ -88,7 +88,7 @@ void LogSumExp::eval_cpu(const std::vector<array>& inputs, array& out) {
       return x;
     } else {
       auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
-      copy_cpu(x, x_copy, CopyType::General, s);
+      copy(x, x_copy, CopyType::General, s);
       encoder.add_temporary(x_copy);
       return x_copy;
     }

mlx/backend/cpu/luf.cpp

@@ -31,7 +31,7 @@ void luf_impl(
   strides[ndim - 1] = M;
   strides[ndim - 2] = 1;
   lu.set_data(allocator::malloc(lu.nbytes()), lu.nbytes(), strides, flags);
-  copy_cpu_inplace(
+  copy_inplace(
       a,
       lu,
       a.shape(),

mlx/backend/cpu/masked_mm.cpp

@@ -124,20 +124,20 @@ void BlockMaskedMM::eval_cpu(const std::vector<array>& inputs, array& out) {
         if (!expand_all && stx == arr.shape(-1) && sty == 1) {
           if (do_copy) {
             array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-            copy_cpu(arr, arr_copy, CopyType::Vector, s);
+            copy(arr, arr_copy, CopyType::Vector, s);
             return std::make_tuple(false, stx, arr_copy, true);
           }
           return std::make_tuple(false, stx, arr, false);
         } else if (!expand_all && stx == 1 && sty == arr.shape(-2)) {
           if (do_copy) {
             array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-            copy_cpu(arr, arr_copy, CopyType::Vector, s);
+            copy(arr, arr_copy, CopyType::Vector, s);
             return std::make_tuple(true, sty, arr_copy, true);
           }
           return std::make_tuple(true, sty, arr, false);
         } else {
           array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-          copy_cpu(arr, arr_copy, CopyType::General, s);
+          copy(arr, arr_copy, CopyType::General, s);
           int64_t stx = arr.shape(-1);
           return std::make_tuple(false, stx, arr_copy, true);
         }
@@ -386,7 +386,7 @@ void GatherMM::eval_cpu(const std::vector<array>& inputs, array& out) {
       return std::make_tuple(true, sty, arr);
     } else {
       temps.push_back(array(arr.shape(), arr.dtype(), nullptr, {}));
-      copy_cpu(arr, temps.back(), CopyType::General, s);
+      copy(arr, temps.back(), CopyType::General, s);
       int64_t stx = arr.shape(-1);
       return std::make_tuple(false, stx, temps.back());
     }
@@ -504,7 +504,7 @@ void SegmentedMM::eval_cpu(const std::vector<array>& inputs, array& out) {
       return std::make_tuple(true, sty, x);
     } else {
       array xc(x.shape(), x.dtype(), nullptr, {});
-      copy_cpu(x, xc, CopyType::General, s);
+      copy(x, xc, CopyType::General, s);
       encoder.add_temporary(xc);
       int64_t stx = x.shape(-1);
       return std::make_tuple(false, stx, xc);

mlx/backend/cpu/matmul.cpp

@@ -81,7 +81,7 @@ void matmul_general(
       return std::make_tuple(true, sty, arr);
     } else {
       temps.push_back(array(arr.shape(), arr.dtype(), nullptr, {}));
-      copy_cpu(arr, temps.back(), CopyType::General, stream);
+      copy(arr, temps.back(), CopyType::General, stream);
       stx = arr.shape(-1);
       return std::make_tuple(false, stx, temps.back());
     }
@@ -142,7 +142,7 @@ void AddMM::eval_cpu(const std::vector<array>& inputs, array& out) {
   CopyType ctype = c.data_size() == 1
       ? CopyType::Scalar
       : (c.flags().row_contiguous ? CopyType::Vector : CopyType::General);
-  copy_cpu(c, out, ctype, stream());
+  copy(c, out, ctype, stream());
   if (inputs[0].shape(-1) == 0) {
     return;
   }

mlx/backend/cpu/primitives.cpp

@@ -22,7 +22,7 @@ void reshape(const array& in, array& out) {
   auto [copy_necessary, out_strides] = prepare_reshape(in, out);
   if (copy_necessary) {
     out.set_data(allocator::malloc(out.nbytes()));
-    copy_cpu_inplace(in, out, CopyType::General, out.primitive().stream());
+    copy_inplace(in, out, CopyType::General, out.primitive().stream());
   } else {
     shared_buffer_reshape(in, out_strides, out);
   }
@@ -175,7 +175,7 @@ void AsType::eval_cpu(const std::vector<array>& inputs, array& out) {
   assert(inputs.size() == 1);
   auto& in = inputs[0];
   CopyType ctype = in.flags().contiguous ? CopyType::Vector : CopyType::General;
-  copy_cpu(in, out, ctype, stream());
+  copy(in, out, ctype, stream());
 }
 
 void Concatenate::eval_cpu(const std::vector<array>& inputs, array& out) {
@@ -198,7 +198,7 @@ void Concatenate::eval_cpu(const std::vector<array>& inputs, array& out) {
     size_t data_offset = strides[axis_] * sizes[i];
     out_slice.copy_shared_buffer(
         out, strides, flags, out_slice.size(), data_offset);
-    copy_cpu_inplace(inputs[i], out_slice, CopyType::GeneralGeneral, stream());
+    copy_inplace(inputs[i], out_slice, CopyType::GeneralGeneral, stream());
   }
 }
 
@@ -211,7 +211,7 @@ void Contiguous::eval_cpu(const std::vector<array>& inputs, array& out) {
        (allow_col_major_ && in.flags().col_contiguous))) {
     out.copy_shared_buffer(in);
   } else {
-    copy_cpu(in, out, CopyType::General, stream());
+    copy(in, out, CopyType::General, stream());
   }
 }
 
@@ -235,7 +235,7 @@ void Full::eval_cpu(const std::vector<array>& inputs, array& out) {
   } else {
     ctype = CopyType::General;
   }
-  copy_cpu(in, out, ctype, stream());
+  copy(in, out, ctype, stream());
 }
 
 void Pad::eval_cpu(const std::vector<array>& inputs, array& out) {
@@ -251,7 +251,7 @@ void Pad::eval_cpu(const std::vector<array>& inputs, array& out) {
   assert(val.dtype() == in.dtype() && in.dtype() == out.dtype());
 
   // Fill output with val
-  copy_cpu(val, out, CopyType::Scalar, stream());
+  copy(val, out, CopyType::Scalar, stream());
 
   // Find offset for start of input values
   size_t data_offset = 0;
@@ -266,7 +266,7 @@ void Pad::eval_cpu(const std::vector<array>& inputs, array& out) {
       out, out.strides(), out.flags(), out_slice.size(), data_offset);
 
   // Copy input values into the slice
-  copy_cpu_inplace(in, out_slice, CopyType::GeneralGeneral, stream());
+  copy_inplace(in, out_slice, CopyType::GeneralGeneral, stream());
 }
 
 void RandomBits::eval_cpu(const std::vector<array>& inputs, array& out) {
@@ -340,7 +340,7 @@ void DynamicSlice::eval_cpu(const std::vector<array>& inputs, array& out) {
   out.set_data(allocator::malloc(out.nbytes()));
   auto [in_offset, donated] =
       compute_dynamic_offset(inputs[1], in.strides(), axes_, stream());
-  copy_cpu_inplace(
+  copy_inplace(
       /* const array& src = */ in,
       /* array& dst = */ out,
       /* const Shape& data_shape = */ out.shape(),
@@ -372,11 +372,11 @@ void DynamicSliceUpdate::eval_cpu(
   auto ctype = in.flags().contiguous && in.size() == in.data_size()
       ? CopyType::Vector
       : CopyType::General;
-  copy_cpu(in, out, in.data_size() == 1 ? CopyType::Scalar : ctype, stream());
+  copy(in, out, in.data_size() == 1 ? CopyType::Scalar : ctype, stream());
 
   auto [out_offset, donated] =
       compute_dynamic_offset(inputs[2], out.strides(), axes_, stream());
-  copy_cpu_inplace(
+  copy_inplace(
       /* const array& src = */ upd,
       /* array& dst = */ out,
       /* const std::vector<int>& data_shape = */ upd.shape(),
@@ -412,14 +412,14 @@ void SliceUpdate::eval_cpu(const std::vector<array>& inputs, array& out) {
   auto ctype = in.flags().contiguous && in.size() == in.data_size()
       ? CopyType::Vector
       : CopyType::General;
-  copy_cpu(in, out, in.data_size() == 1 ? CopyType::Scalar : ctype, stream());
+  copy(in, out, in.data_size() == 1 ? CopyType::Scalar : ctype, stream());
 
   // Calculate out strides, initial offset and if copy needs to be made
   auto [data_offset, out_strides] =
       prepare_slice(out, start_indices_, strides_);
 
   // Do copy
-  copy_cpu_inplace(
+  copy_inplace(
       /* const array& src = */ upd,
       /* array& dst = */ out,
       /* const std::vector<int>& data_shape = */ upd.shape(),
@@ -456,9 +456,9 @@ void View::eval_cpu(const std::vector<array>& inputs, array& out) {
     if (in.dtype() == bool_) {
       auto in_tmp = array(in.shape(), uint8, nullptr, {});
       in_tmp.copy_shared_buffer(in);
-      copy_cpu_inplace(in_tmp, tmp, CopyType::General, stream());
+      copy_inplace(in_tmp, tmp, CopyType::General, stream());
     } else {
-      copy_cpu_inplace(in, tmp, CopyType::General, stream());
+      copy_inplace(in, tmp, CopyType::General, stream());
     }
 
     auto flags = out.flags();

mlx/backend/cpu/qrf.cpp

@@ -26,7 +26,7 @@ void qrf_impl(const array& a, array& q, array& r, Stream stream) {
   strides[in.ndim() - 2] = 1;
   strides[in.ndim() - 1] = M;
   in.set_data(allocator::malloc(in.nbytes()), in.nbytes(), strides, flags);
-  copy_cpu_inplace(a, in, CopyType::GeneralGeneral, stream);
+  copy_inplace(a, in, CopyType::GeneralGeneral, stream);
   auto& encoder = cpu::get_command_encoder(stream);
   q.set_data(allocator::malloc(q.nbytes()));
   r.set_data(allocator::malloc(r.nbytes()));

mlx/backend/cpu/quantized.cpp

@@ -529,7 +529,7 @@ void QuantizedMatmul::eval_cpu(const std::vector<array>& inputs, array& out) {
       return arr;
     } else {
       temps.push_back(array(arr.shape(), arr.dtype(), nullptr, {}));
-      copy_cpu(arr, temps.back(), CopyType::General, s);
+      copy(arr, temps.back(), CopyType::General, s);
       return temps.back();
     }
   };
@@ -579,7 +579,7 @@ void GatherQMM::eval_cpu(const std::vector<array>& inputs, array& out) {
       return arr;
     } else {
       temps.push_back(array(arr.shape(), arr.dtype(), nullptr, {}));
-      copy_cpu(arr, temps.back(), CopyType::General, s);
+      copy(arr, temps.back(), CopyType::General, s);
       return temps.back();
     }
   };
@@ -713,7 +713,7 @@ void fast::AffineQuantize::eval_cpu(
       return std::make_pair(arr, false);
     } else {
       array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-      copy_cpu(arr, arr_copy, CopyType::General, s);
+      copy(arr, arr_copy, CopyType::General, s);
       return std::make_pair(arr_copy, true);
     }
   };

mlx/backend/cpu/scan.cpp

@@ -251,7 +251,7 @@ void Scan::eval_cpu(const std::vector<array>& inputs, array& out) {
   auto in = inputs[0];
   if (!in.flags().row_contiguous) {
     array arr_copy(in.shape(), in.dtype(), nullptr, {});
-    copy_cpu(in, arr_copy, CopyType::General, stream());
+    copy(in, arr_copy, CopyType::General, stream());
     in = arr_copy;
     encoder.add_temporary(arr_copy);
   }

mlx/backend/cpu/softmax.cpp

@@ -132,7 +132,7 @@ void Softmax::eval_cpu(const std::vector<array>& inputs, array& out) {
       return x;
     } else {
       array x_copy(x.shape(), x.dtype(), nullptr, {});
-      copy_cpu(x, x_copy, CopyType::General, s);
+      copy(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }

mlx/backend/cpu/sort.cpp

@@ -334,10 +334,8 @@ void Sort::eval_cpu(const std::vector<array>& inputs, array& out) {
   auto& in = inputs[0];
 
   // Copy input to output
-  CopyType ctype = (in.flags().contiguous && in.strides()[axis_] != 0)
-      ? CopyType::Vector
-      : CopyType::General;
-  copy_cpu(in, out, ctype, stream());
+  CopyType ctype = in.flags().contiguous ? CopyType::Vector : CopyType::General;
+  copy(in, out, ctype, stream());
 
   auto& encoder = cpu::get_command_encoder(stream());
   encoder.set_output_array(out);
@@ -428,10 +426,8 @@ void Partition::eval_cpu(const std::vector<array>& inputs, array& out) {
   auto& in = inputs[0];
 
   // Copy input to output
-  CopyType ctype = (in.flags().contiguous && in.strides()[axis_] != 0)
-      ? CopyType::Vector
-      : CopyType::General;
-  copy_cpu(in, out, ctype, stream());
+  CopyType ctype = in.flags().contiguous ? CopyType::Vector : CopyType::General;
+  copy(in, out, ctype, stream());
 
   auto& encoder = cpu::get_command_encoder(stream());
   encoder.set_output_array(out);

mlx/backend/cpu/svd.cpp

@@ -31,7 +31,7 @@ void svd_impl(
 
   // lapack clobbers the input, so we have to make a copy.
   array in(a.shape(), a.dtype(), nullptr, {});
-  copy_cpu(
+  copy(
       a,
       in,
       a.flags().row_contiguous ? CopyType::Vector : CopyType::General,

mlx/backend/cuda/CMakeLists.txt

@@ -42,7 +42,6 @@ target_sources(
           ${CMAKE_CURRENT_SOURCE_DIR}/ternary.cu
           ${CMAKE_CURRENT_SOURCE_DIR}/unary.cu
           ${CMAKE_CURRENT_SOURCE_DIR}/utils.cpp
-          ${CMAKE_CURRENT_SOURCE_DIR}/quantized.cu
           ${CMAKE_CURRENT_SOURCE_DIR}/worker.cpp)
 
 target_compile_definitions(mlx PRIVATE MLX_USE_CUDA)
@@ -126,7 +125,3 @@ target_link_libraries(mlx PRIVATE CUDA::nvrtc CUDA::cuda_driver)
 # Suppress nvcc warnings on MLX headers.
 target_compile_options(mlx PRIVATE $<$<COMPILE_LANGUAGE:CUDA>:-Xcudafe
                                    --diag_suppress=997>)
-
-# Install CCCL headers for JIT.
-install(DIRECTORY ${cccl_SOURCE_DIR}/include/cuda
-        DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/cccl)

mlx/backend/cuda/binary.cu

@@ -3,6 +3,7 @@
 #include "mlx/backend/common/binary.h"
 #include "mlx/backend/cuda/device.h"
 #include "mlx/backend/cuda/device/binary_ops.cuh"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/kernel_utils.cuh"
 #include "mlx/dtype_utils.h"
 #include "mlx/primitives.h"
@@ -176,7 +177,7 @@ template <typename Op>
 void binary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    std::string_view op,
     const Stream& s) {
   assert(inputs.size() > 1);
   const auto& a = inputs[0];
@@ -290,7 +291,7 @@ template <typename Op>
 void binary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    std::string_view op,
     const Stream& s) {
   auto& a = inputs[0];
   auto& b = inputs[1];
@@ -299,11 +300,11 @@ void binary_op_gpu(
   binary_op_gpu_inplace<Op>(inputs, out, op, s);
 }
 
-#define BINARY_GPU(func)                                              \
-  void func::eval_gpu(const std::vector<array>& inputs, array& out) { \
-    nvtx3::scoped_range r(#func "::eval_gpu");                        \
-    auto& s = out.primitive().stream();                               \
-    binary_op_gpu<cu::func>(inputs, out, name(), s);                  \
+#define BINARY_GPU(func)                                                 \
+  void func::eval_gpu(const std::vector<array>& inputs, array& out) {    \
+    nvtx3::scoped_range r(#func "::eval_gpu");                           \
+    auto& s = out.primitive().stream();                                  \
+    binary_op_gpu<cu::func>(inputs, out, get_primitive_string(this), s); \
   }
 
 BINARY_GPU(Add)
@@ -327,31 +328,33 @@ BINARY_GPU(Subtract)
 void Equal::eval_gpu(const std::vector<array>& inputs, array& out) {
   nvtx3::scoped_range r("Equal::eval_gpu");
   auto& s = out.primitive().stream();
+  auto op = get_primitive_string(this);
   if (equal_nan_) {
-    binary_op_gpu<cu::NaNEqual>(inputs, out, name(), s);
+    binary_op_gpu<cu::NaNEqual>(inputs, out, op, s);
   } else {
-    binary_op_gpu<cu::Equal>(inputs, out, name(), s);
+    binary_op_gpu<cu::Equal>(inputs, out, op, s);
   }
 }
 
 void BitwiseBinary::eval_gpu(const std::vector<array>& inputs, array& out) {
   nvtx3::scoped_range r("BitwiseBinary::eval_gpu");
   auto& s = out.primitive().stream();
+  auto op = get_primitive_string(this);
   switch (op_) {
     case BitwiseBinary::And:
-      binary_op_gpu<cu::BitwiseAnd>(inputs, out, name(), s);
+      binary_op_gpu<cu::BitwiseAnd>(inputs, out, op, s);
       break;
     case BitwiseBinary::Or:
-      binary_op_gpu<cu::BitwiseOr>(inputs, out, name(), s);
+      binary_op_gpu<cu::BitwiseOr>(inputs, out, op, s);
       break;
     case BitwiseBinary::Xor:
-      binary_op_gpu<cu::BitwiseXor>(inputs, out, name(), s);
+      binary_op_gpu<cu::BitwiseXor>(inputs, out, op, s);
       break;
     case BitwiseBinary::LeftShift:
-      binary_op_gpu<cu::LeftShift>(inputs, out, name(), s);
+      binary_op_gpu<cu::LeftShift>(inputs, out, op, s);
       break;
     case BitwiseBinary::RightShift:
-      binary_op_gpu<cu::RightShift>(inputs, out, name(), s);
+      binary_op_gpu<cu::RightShift>(inputs, out, op, s);
       break;
   }
 }

mlx/backend/cuda/binary_two.cu

@@ -3,6 +3,7 @@
 #include "mlx/backend/common/binary.h"
 #include "mlx/backend/cuda/device.h"
 #include "mlx/backend/cuda/device/binary_ops.cuh"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/kernel_utils.cuh"
 #include "mlx/dtype_utils.h"
 #include "mlx/primitives.h"
@@ -183,7 +184,7 @@ template <typename Op>
 void binary_two_op_gpu_inplace(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    std::string_view op,
     const Stream& s) {
   assert(inputs.size() > 1);
   const auto& a = inputs[0];
@@ -313,7 +314,7 @@ template <typename Op>
 void binary_two_op_gpu(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    std::string_view op,
     const Stream& s) {
   auto& a = inputs[0];
   auto& b = inputs[1];
@@ -328,7 +329,7 @@ void DivMod::eval_gpu(
     std::vector<array>& outputs) {
   nvtx3::scoped_range r("DivMod::eval_gpu");
   auto& s = outputs[0].primitive().stream();
-  binary_two_op_gpu<cu::DivMod>(inputs, outputs, name(), s);
+  binary_two_op_gpu<cu::DivMod>(inputs, outputs, get_primitive_string(this), s);
 }
 
 } // namespace mlx::core

mlx/backend/cuda/compiled.cpp

@@ -53,10 +53,9 @@ struct FusedKernelBuilder {
 
     // Build function signature.
     if (contiguous) {
-      os += "template <typename IdxT = uint32_t, int work_per_thread = 1>\n";
+      os += "template <typename IdxT = uint32_t>\n";
     } else {
-      os +=
-          "template <int NDIM, typename IdxT = uint32_t, int work_per_thread = 1>\n";
+      os += "template <int NDIM, typename IdxT = uint32_t>\n";
     }
     os += fmt::format("__global__ void {}(\n", kernel_name + name);
     for (size_t i = 0; i < params.size(); ++i) {
@@ -68,46 +67,12 @@ struct FusedKernelBuilder {
     }
     os += ") {\n";
 
-    // Index. For non contiguous kernels we create a separate index
-    // variable per variable otherwise everyone uses `index`.
+    // Index.
     os +=
-        "  IdxT index = cg::this_grid().thread_rank() * work_per_thread;\n"
+        "  IdxT index = cg::this_grid().thread_rank();\n"
         "  if (index >= size) {\n"
         "    return;\n"
         "  }\n";
-    if (!contiguous) {
-      for (size_t i = 0; i < inputs.size(); ++i) {
-        const auto& x = inputs[i];
-        const std::string& xname = namer.get_name(x);
-        if (is_scalar(x) || is_constant(i)) {
-          continue;
-        }
-        os += "  IdxT " + xname + "_idx = 0;\n";
-      }
-      os += "  {\n";
-      os += "    IdxT loc = index;\n";
-      os +=
-          "    #pragma unroll\n"
-          "    for (int i = NDIM - 1; i >= 0; i--) {\n";
-      for (size_t i = 0; i < inputs.size(); ++i) {
-        const auto& x = inputs[i];
-        const std::string& xname = namer.get_name(x);
-        if (is_scalar(x) || is_constant(i)) {
-          continue;
-        }
-        os += "      " + xname + "_idx += (loc \% shape[i]) * IdxT(" + xname +
-            "_strides[i]);\n";
-      }
-      os +=
-          "      loc /= shape[i];\n"
-          "    }\n"
-          "  }\n";
-    }
-
-    // Work loop
-    os +=
-        "\n"
-        "  for (int i = 0; i < work_per_thread && index < size; i++) {\n";
 
     // Read inputs.
     for (size_t i = 0; i < inputs.size(); ++i) {
@@ -124,9 +89,12 @@ struct FusedKernelBuilder {
       } else if (contiguous) {
         value = fmt::format("{}[index]", xname);
       } else {
-        value = fmt::format("{}[{}_idx]", xname, xname);
+        std::string index = fmt::format(
+            "elem_to_loc_nd<NDIM>(index, shape.data(), {}_strides.data())",
+            xname);
+        value = fmt::format("{}[{}]", xname, index);
       }
-      os += fmt::format("    {} tmp_{} = {};\n", type, xname, value);
+      os += fmt::format("  {} tmp_{} = {};\n", type, xname, value);
     }
 
     // Write tape.
@@ -138,37 +106,23 @@ struct FusedKernelBuilder {
         value = fmt::format(
             "static_cast<{}>(tmp_{})", type, namer.get_name(x.inputs()[0]));
       } else {
-        value = x.primitive().name();
+        std::ostringstream ss;
+        x.primitive().print(ss);
+        value = ss.str();
         value += "{}(";
         for (size_t i = 0; i < x.inputs().size() - 1; ++i) {
           value += fmt::format("tmp_{}, ", namer.get_name(x.inputs()[i]));
         }
         value += fmt::format("tmp_{})", namer.get_name(x.inputs().back()));
       }
-      os += fmt::format("    {} tmp_{} = {};\n", type, xname, value);
+      os += fmt::format("  {} tmp_{} = {};\n", type, xname, value);
     }
 
     // Write output.
     for (const auto& x : outputs) {
-      os += fmt::format("    {0}[index] = tmp_{0};\n", namer.get_name(x));
+      os += fmt::format("  {0}[index] = tmp_{0};\n", namer.get_name(x));
     }
 
-    // End of work loop
-    os +=
-        "\n"
-        "    index++;\n";
-    if (!contiguous) {
-      for (size_t i = 0; i < inputs.size(); ++i) {
-        const auto& x = inputs[i];
-        const std::string& xname = namer.get_name(x);
-        if (is_scalar(x) || is_constant(i)) {
-          continue;
-        }
-        os += "    " + xname + "_idx += " + xname + "_strides[NDIM - 1];\n";
-      }
-    }
-    os += "  }\n";
-
     os += "}\n";
   }
 };
@@ -204,28 +158,15 @@ void Compiled::eval_gpu(
     builder.build("_strided", false);
     builder.os += "\n} // namespace mlx::core::cu\n";
     // Build kernel names.
-    std::vector<std::string> kernel_names;
-    for (auto work_per_thread : std::array<int, 2>{1, 4}) {
+    std::vector<std::string> kernel_names = {
+        fmt::format("mlx::core::cu::{}_contiguous<uint32_t>", lib_name()),
+        fmt::format("mlx::core::cu::{}_contiguous<int64_t>", lib_name()),
+    };
+    for (int i = 1; i <= MAX_NDIM; ++i) {
       kernel_names.push_back(fmt::format(
-          "mlx::core::cu::{}_contiguous<uint32_t, {}>",
-          lib_name(),
-          work_per_thread));
-      kernel_names.push_back(fmt::format(
-          "mlx::core::cu::{}_contiguous<int64_t, {}>",
-          lib_name(),
-          work_per_thread));
-      for (int i = 1; i <= MAX_NDIM; ++i) {
-        kernel_names.push_back(fmt::format(
-            "mlx::core::cu::{}_strided<{}, uint32_t, {}>",
-            lib_name(),
-            i,
-            work_per_thread));
-        kernel_names.push_back(fmt::format(
-            "mlx::core::cu::{}_strided<{}, int64_t, {}>",
-            lib_name(),
-            i,
-            work_per_thread));
-      }
+          "mlx::core::cu::{}_strided<{}, uint32_t>", lib_name(), i));
+      kernel_names.push_back(
+          fmt::format("mlx::core::cu::{}_strided<{}, int64_t>", lib_name(), i));
     }
     return std::make_pair(std::move(builder.os), std::move(kernel_names));
   });
@@ -268,21 +209,13 @@ void Compiled::eval_gpu(
     args.append<uint32_t>(outputs[0].data_size());
   }
 
-  // Choose work per thread
-  int work_per_thread = 4;
-  if (!contiguous && shape.back() % work_per_thread != 0) {
-    work_per_thread = 1;
-  }
-
   // Launch kernel.
   const char* index_type = large ? "int64_t" : "uint32_t";
   std::string kernel_name = fmt::format("mlx::core::cu::{}", lib_name());
   if (contiguous) {
-    kernel_name +=
-        fmt::format("_contiguous<{}, {}>", index_type, work_per_thread);
+    kernel_name += fmt::format("_contiguous<{}>", index_type);
   } else {
-    kernel_name += fmt::format(
-        "_strided<{}, {}, {}>", shape.size(), index_type, work_per_thread);
+    kernel_name += fmt::format("_strided<{}, {}>", shape.size(), index_type);
   }
   auto& encoder = cu::get_command_encoder(s);
   for (const auto& in : inputs) {
@@ -293,8 +226,7 @@ void Compiled::eval_gpu(
   }
 
   auto kernel = mod.get_kernel(kernel_name);
-  auto [num_blocks, block_dims] =
-      get_launch_args(kernel, outputs[0], large, work_per_thread);
+  auto [num_blocks, block_dims] = get_launch_args(kernel, outputs[0], large);
   encoder.add_kernel_node(kernel, num_blocks, block_dims, args.args());
 }
 

mlx/backend/cuda/device.cpp

@@ -66,6 +66,7 @@ CommandEncoder& Device::get_command_encoder(Stream s) {
 }
 
 CommandEncoder::CaptureContext::CaptureContext(CommandEncoder& enc) : enc(enc) {
+  CHECK_CUDA_ERROR(cudaGraphCreate(&graph, 0));
   CHECK_CUDA_ERROR(
       cudaStreamBeginCapture(enc.stream(), cudaStreamCaptureModeGlobal));
 }

mlx/backend/cuda/device/atomic_ops.cuh

@@ -2,7 +2,7 @@
 
 #pragma once
 
-#include "mlx/backend/cuda/device/complex.cuh"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/device/fp16_math.cuh"
 
 #include <cuda/atomic>
@@ -48,7 +48,7 @@ inline __device__ void atomic_add(__half* out, __half val) {
   atomicAdd(out, val);
 }
 
-inline __device__ void atomic_add(complex64_t* out, complex64_t val) {
+inline __device__ void atomic_add(cuComplex* out, cuComplex val) {
 #if __CUDA_ARCH__ < 900
   atomic_add_general(out, val);
 #else
@@ -58,7 +58,12 @@ inline __device__ void atomic_add(complex64_t* out, complex64_t val) {
 
 inline __device__ void atomic_add(__nv_bfloat16* out, __nv_bfloat16 val) {
 #if __CUDA_ARCH__ < 800
+#if CCCL_VERSION >= 2008000
   atomic_add_general(out, val);
+#else
+  bool cccl_version_too_old_for_bfloat16_atomic_add = false;
+  assert(cccl_version_too_old_for_bfloat16_atomic_add);
+#endif
 #else
   atomicAdd(out, val);
 #endif

mlx/backend/cuda/device/binary_ops.cuh

@@ -44,7 +44,7 @@ struct Remainder {
       } else {
         return x % y;
       }
-    } else if constexpr (is_complex_v<T>) {
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
       return x % y;
     } else {
       T r = fmod(x, y);
@@ -66,12 +66,14 @@ struct Equal {
 struct NaNEqual {
   template <typename T>
   __device__ bool operator()(T x, T y) {
-    if constexpr (is_complex_v<T>) {
+    if constexpr (std::is_same_v<T, cuComplex>) {
       return x == y ||
-          (isnan(x.real()) && isnan(y.real()) && isnan(x.imag()) &&
-           isnan(y.imag())) ||
-          (x.real() == y.real() && isnan(x.imag()) && isnan(y.imag())) ||
-          (isnan(x.real()) && isnan(y.real()) && x.imag() == y.imag());
+          (isnan(cuCrealf(x)) && isnan(cuCrealf(y)) && isnan(cuCimagf(x)) &&
+           isnan(cuCimagf(y))) ||
+          (cuCrealf(x) == cuCrealf(y) && isnan(cuCimagf(x)) &&
+           isnan(cuCimagf(y))) ||
+          (isnan(cuCrealf(x)) && isnan(cuCrealf(y)) &&
+           cuCimagf(x) == cuCimagf(y));
     } else {
       return x == y || (isnan(x) && isnan(y));
     }
@@ -109,17 +111,17 @@ struct LessEqual {
 struct LogAddExp {
   template <typename T>
   __device__ T operator()(T x, T y) {
-    if constexpr (is_complex_v<T>) {
-      if (isnan(x.real()) || isnan(x.imag()) || isnan(y.real()) ||
-          isnan(y.imag())) {
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      if (isnan(cuCrealf(x)) || isnan(cuCimagf(x)) || isnan(cuCrealf(y)) ||
+          isnan(cuCimagf(y))) {
         return {
             cuda::std::numeric_limits<float>::quiet_NaN(),
             cuda::std::numeric_limits<float>::quiet_NaN()};
       }
-      auto max = x.real() > y.real() ? x : y;
-      auto min = x.real() < y.real() ? x : y;
-      auto min_real = min.real();
-      auto max_real = max.real();
+      auto max = cuCrealf(x) > cuCrealf(y) ? x : y;
+      auto min = cuCrealf(x) < cuCrealf(y) ? x : y;
+      auto min_real = cuCrealf(min);
+      auto max_real = cuCrealf(max);
       if (!isfinite(min_real) && (min_real == max_real)) {
         if (min_real < 0) {
           return min;
@@ -148,8 +150,8 @@ struct Maximum {
   __device__ T operator()(T x, T y) {
     if constexpr (cuda::std::is_integral_v<T>) {
       return max(x, y);
-    } else if constexpr (is_complex_v<T>) {
-      if (isnan(x.real()) || isnan(x.imag())) {
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      if (isnan(cuCrealf(x)) || isnan(cuCimagf(x))) {
         return x;
       }
       return x > y ? x : y;
@@ -167,8 +169,8 @@ struct Minimum {
   __device__ T operator()(T x, T y) {
     if constexpr (cuda::std::is_integral_v<T>) {
       return min(x, y);
-    } else if constexpr (is_complex_v<T>) {
-      if (isnan(x.real()) || isnan(x.imag())) {
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      if (isnan(cuCrealf(x)) || isnan(cuCimagf(x))) {
         return x;
       }
       return x < y ? x : y;
@@ -191,8 +193,8 @@ struct Multiply {
 struct NotEqual {
   template <typename T>
   __device__ bool operator()(T x, T y) {
-    if constexpr (is_complex_v<T>) {
-      return x.real() != y.real() || x.imag() != y.imag();
+    if constexpr (std::is_same_v<T, cuComplex>) {
+      return cuCrealf(x) != cuCrealf(y) || cuCimagf(x) != cuCimagf(y);
     } else {
       return x != y;
     }
@@ -212,8 +214,19 @@ struct Power {
         base *= base;
       }
       return res;
-    } else if constexpr (is_complex_v<T>) {
-      return pow(base, exp);
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      if (base.y == 0 && base.x == 0) {
+        if (isnan(exp.x) || isnan(exp.y)) {
+          auto nan = cuda::std::numeric_limits<float>::quiet_NaN();
+          return make_cuFloatComplex(nan, nan);
+        }
+        return make_cuFloatComplex(0.0, 0.0);
+      }
+      auto x_theta = atan2f(base.y, base.x);
+      auto x_ln_r = 0.5 * logf(base.x * base.x + base.y * base.y);
+      auto mag = expf(exp.x * x_ln_r - exp.y * x_theta);
+      auto phase = exp.y * x_ln_r + exp.x * x_theta;
+      return make_cuFloatComplex(mag * cosf(phase), mag * sinf(phase));
     } else {
       return powf(base, exp);
     }

mlx/backend/cuda/device/cast_op.cuh

@@ -2,8 +2,7 @@
 
 #pragma once
 
-#include "mlx/backend/cuda/device/complex.cuh"
-
+#include <cuComplex.h>
 #include <cuda_bf16.h>
 #include <cuda_fp16.h>
 #include <thrust/iterator/transform_iterator.h>
@@ -21,43 +20,50 @@ struct CastOp {
 };
 
 // Castings between complex and boolean.
-template <typename T>
-struct CastOp<complex_t<T>, bool> {
+// TODO: Should make a custom complex type.
+template <>
+struct CastOp<cuComplex, bool> {
   static constexpr bool is_castable = true;
 
-  __device__ bool operator()(complex_t<T> x) {
-    return x.real() != 0 && x.imag() != 0;
+  __device__ bool operator()(cuComplex x) {
+    return x.x != 0 && x.y != 0;
   }
 };
 
-template <typename T>
-struct CastOp<bool, complex_t<T>> {
+template <>
+struct CastOp<bool, cuComplex> {
   static constexpr bool is_castable = true;
 
-  __device__ complex_t<T> operator()(bool x) {
-    return x ? complex_t<T>{1, 1} : complex_t<T>{0, 0};
+  __device__ cuComplex operator()(bool x) {
+    return x ? make_cuFloatComplex(1, 1) : make_cuFloatComplex(0, 0);
   }
 };
 
 // Converting a complex number to real number discards the imaginary part.
-template <typename T, typename DstT>
-struct CastOp<complex_t<T>, DstT, cuda::std::enable_if_t<!is_complex_v<DstT>>> {
-  static constexpr bool is_castable = cuda::std::is_convertible_v<T, DstT>;
+template <typename DstT>
+struct CastOp<
+    cuComplex,
+    DstT,
+    cuda::std::enable_if_t<!cuda::std::is_same_v<cuComplex, DstT>>> {
+  static constexpr bool is_castable = cuda::std::is_convertible_v<float, DstT>;
 
-  __device__ DstT operator()(complex_t<T> x) {
-    static_assert(!is_complex_v<DstT>);
-    return static_cast<DstT>(x.real());
+  __device__ DstT operator()(cuComplex x) {
+    static_assert(!cuda::std::is_same_v<cuComplex, DstT>);
+    return static_cast<DstT>(cuCrealf(x));
   }
 };
 
 // Allow converting a real number to complex number.
-template <typename SrcT, typename T>
-struct CastOp<SrcT, complex_t<T>, cuda::std::enable_if_t<!is_complex_v<SrcT>>> {
-  static constexpr bool is_castable = cuda::std::is_convertible_v<SrcT, T>;
+template <typename SrcT>
+struct CastOp<
+    SrcT,
+    cuComplex,
+    cuda::std::enable_if_t<!cuda::std::is_same_v<SrcT, cuComplex>>> {
+  static constexpr bool is_castable = cuda::std::is_convertible_v<SrcT, float>;
 
-  __device__ complex_t<T> operator()(SrcT x) {
-    static_assert(!is_complex_v<SrcT>);
-    return complex_t<T>{static_cast<T>(x), 0};
+  __device__ cuComplex operator()(SrcT x) {
+    static_assert(!cuda::std::is_same_v<SrcT, cuComplex>);
+    return cuComplex{static_cast<float>(x), 0};
   }
 };
 
@@ -82,7 +88,8 @@ struct CastOp<
     SrcT,
     DstT,
     cuda::std::enable_if_t<
-        !cuda::std::is_convertible_v<SrcT, DstT> && !is_complex_v<SrcT> &&
+        !cuda::std::is_convertible_v<SrcT, DstT> &&
+        !cuda::std::is_same_v<SrcT, cuComplex> &&
         (cuda::std::is_same_v<DstT, __half> ||
          cuda::std::is_same_v<DstT, __nv_bfloat16>)>> {
   static constexpr bool is_castable = true;
@@ -97,7 +104,8 @@ struct CastOp<
     SrcT,
     DstT,
     cuda::std::enable_if_t<
-        !cuda::std::is_convertible_v<SrcT, DstT> && !is_complex_v<SrcT> &&
+        !cuda::std::is_convertible_v<SrcT, DstT> &&
+        !cuda::std::is_same_v<DstT, cuComplex> &&
         !cuda::std::is_same_v<DstT, __half> &&
         !cuda::std::is_same_v<DstT, __nv_bfloat16> &&
         (cuda::std::is_same_v<SrcT, __half> ||

mlx/backend/cuda/device/cexpf.cuh

@@ -0,0 +1,138 @@
+// Copyright © 2025 Apple Inc.
+// Copyright © 2008-2013 NVIDIA Corporation
+// Copyright © 2013 Filipe RNC Maia
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+//
+// Forked from
+// https://github.com/NVIDIA/cccl/blob/main/thrust/thrust/detail/complex/cexpf.h
+
+// TODO: We should use thrust::exp but the thrust header in old CUDA versions
+// can not be used in JIT.
+
+#pragma once
+
+#include <cuComplex.h>
+#include <cuda/std/cstdint>
+
+namespace mlx::core::cu::detail {
+
+using ieee_float_shape_type = union {
+  float value;
+  uint32_t word;
+};
+
+inline __device__ void get_float_word(uint32_t& i, float d) {
+  ieee_float_shape_type gf_u;
+  gf_u.value = (d);
+  (i) = gf_u.word;
+}
+
+inline __device__ void get_float_word(int32_t& i, float d) {
+  ieee_float_shape_type gf_u;
+  gf_u.value = (d);
+  (i) = gf_u.word;
+}
+
+inline __device__ void set_float_word(float& d, uint32_t i) {
+  ieee_float_shape_type sf_u;
+  sf_u.word = (i);
+  (d) = sf_u.value;
+}
+
+inline __device__ float frexp_expf(float x, int* expt) {
+  const uint32_t k = 235;
+  const float kln2 = 162.88958740F;
+
+  float exp_x;
+  uint32_t hx;
+
+  exp_x = expf(x - kln2);
+  get_float_word(hx, exp_x);
+  *expt = (hx >> 23) - (0x7f + 127) + k;
+  set_float_word(exp_x, (hx & 0x7fffff) | ((0x7f + 127) << 23));
+  return exp_x;
+}
+
+inline __device__ cuComplex ldexp_cexpf(cuComplex z, int expt) {
+  float x, y, exp_x, scale1, scale2;
+  int ex_expt, half_expt;
+
+  x = cuCrealf(z);
+  y = cuCimagf(z);
+  exp_x = frexp_expf(x, &ex_expt);
+  expt += ex_expt;
+
+  half_expt = expt / 2;
+  set_float_word(scale1, (0x7f + half_expt) << 23);
+  half_expt = expt - half_expt;
+  set_float_word(scale2, (0x7f + half_expt) << 23);
+
+  return cuComplex{
+      cosf(y) * exp_x * scale1 * scale2, sinf(y) * exp_x * scale1 * scale2};
+}
+
+inline __device__ cuComplex cexpf(const cuComplex& z) {
+  float x, y, exp_x;
+  uint32_t hx, hy;
+
+  const uint32_t exp_ovfl = 0x42b17218, cexp_ovfl = 0x43400074;
+
+  x = cuCrealf(z);
+  y = cuCimagf(z);
+
+  get_float_word(hy, y);
+  hy &= 0x7fffffff;
+
+  /* cexp(x + I 0) = exp(x) + I 0 */
+  if (hy == 0) {
+    return cuComplex{expf(x), y};
+  }
+  get_float_word(hx, x);
+  /* cexp(0 + I y) = cos(y) + I sin(y) */
+  if ((hx & 0x7fffffff) == 0) {
+    return cuComplex{cosf(y), sinf(y)};
+  }
+  if (hy >= 0x7f800000) {
+    if ((hx & 0x7fffffff) != 0x7f800000) {
+      /* cexp(finite|NaN +- I Inf|NaN) = NaN + I NaN */
+      return cuComplex{y - y, y - y};
+    } else if (hx & 0x80000000) {
+      /* cexp(-Inf +- I Inf|NaN) = 0 + I 0 */
+      return cuComplex{0.0, 0.0};
+    } else {
+      /* cexp(+Inf +- I Inf|NaN) = Inf + I NaN */
+      return cuComplex{x, y - y};
+    }
+  }
+
+  if (hx >= exp_ovfl && hx <= cexp_ovfl) {
+    /*
+     * x is between 88.7 and 192, so we must scale to avoid
+     * overflow in expf(x).
+     */
+    return ldexp_cexpf(z, 0);
+  } else {
+    /*
+     * Cases covered here:
+     *  -  x < exp_ovfl and exp(x) won't overflow (common case)
+     *  -  x > cexp_ovfl, so exp(x) * s overflows for all s > 0
+     *  -  x = +-Inf (generated by exp())
+     *  -  x = NaN (spurious inexact exception from y)
+     */
+    exp_x = expf(x);
+    return cuComplex{exp_x * cosf(y), exp_x * sinf(y)};
+  }
+}
+
+} // namespace mlx::core::cu::detail

mlx/backend/cuda/device/complex.cuh

@@ -1,60 +0,0 @@
-// Copyright © 2025 Apple Inc.
-
-#pragma once
-
-// Make multiplication and division faster.
-#define LIBCUDACXX_ENABLE_SIMPLIFIED_COMPLEX_OPERATIONS
-
-#include <cuda/std/complex>
-#include <cuda/std/type_traits>
-
-namespace mlx::core::cu {
-
-// TODO: Consider using a faster implementation as cuda::std::complex has to
-// conform to C++ standard.
-template <typename T>
-using complex_t = cuda::std::complex<T>;
-
-using complex64_t = complex_t<float>;
-using complex128_t = complex_t<double>;
-
-template <typename T>
-struct is_complex : cuda::std::false_type {};
-
-template <typename T>
-struct is_complex<cuda::std::complex<T>> : cuda::std::true_type {};
-
-template <typename T>
-inline constexpr bool is_complex_v = is_complex<T>::value;
-
-// cuda::std::complex is missing some operators.
-template <typename T>
-inline __host__ __device__ complex_t<T> operator%(
-    complex_t<T> a,
-    complex_t<T> b) {
-  T r = a.real() - floor(a.real() / b.real()) * b.real();
-  T i = a.imag() - floor(a.imag() / b.imag()) * b.imag();
-  return complex_t<T>{r, i};
-}
-
-template <typename T>
-inline __host__ __device__ bool operator>(complex_t<T> a, complex_t<T> b) {
-  return (a.real() > b.real()) || (a.real() == b.real() && a.imag() > b.imag());
-}
-
-template <typename T>
-inline __host__ __device__ bool operator<(complex_t<T> a, complex_t<T> b) {
-  return operator>(b, a);
-}
-
-template <typename T>
-inline __host__ __device__ bool operator<=(complex_t<T> a, complex_t<T> b) {
-  return !(a > b);
-}
-
-template <typename T>
-inline __host__ __device__ bool operator>=(complex_t<T> a, complex_t<T> b) {
-  return !(a < b);
-}
-
-} // namespace mlx::core::cu

mlx/backend/cuda/device/cucomplex_math.cuh

@@ -0,0 +1,240 @@
+// Copyright © 2025 Apple Inc.
+// Copyright © 2017-2024 The Simons Foundation, Inc.
+//
+// FINUFFT is licensed under the Apache License, Version 2.0 (the
+// "License"); you may not use this file except in compliance with the
+// License.  You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+//
+// Forked from
+// https://github.com/flatironinstitute/finufft/blob/main/include/cufinufft/contrib/helper_math.h
+
+#pragma once
+
+#include <cuComplex.h>
+
+// This header provides some helper functions for cuComplex types.
+// It mainly wraps existing CUDA implementations to provide operator overloads
+// e.g. cuAdd, cuSub, cuMul, cuDiv, cuCreal, cuCimag, cuCabs, cuCarg, cuConj are
+// all provided by CUDA
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator+(const cuDoubleComplex& a, const cuDoubleComplex& b) {
+  return cuCadd(a, b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator-(const cuDoubleComplex& a, const cuDoubleComplex& b) {
+  return cuCsub(a, b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator*(const cuDoubleComplex& a, const cuDoubleComplex& b) {
+  return cuCmul(a, b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator/(const cuDoubleComplex& a, const cuDoubleComplex& b) {
+  return cuCdiv(a, b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator%(const cuDoubleComplex& a, const cuDoubleComplex& b) {
+  double r = cuCreal(a) - (floorf(cuCreal(a) / cuCreal(b)) * cuCreal(b));
+  double i = cuCimag(a) - (floorf(cuCimag(a) / cuCimag(b)) * cuCimag(b));
+  return make_cuDoubleComplex(r, i);
+}
+
+__forceinline__ __host__ __device__ bool operator==(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  return cuCreal(a) == cuCreal(b) && cuCimag(a) == cuCimag(b);
+}
+
+__forceinline__ __host__ __device__ bool operator!=(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  return !(a == b);
+}
+
+__forceinline__ __host__ __device__ bool operator>(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  double mag_a = sqrt(cuCreal(a) * cuCreal(a) + cuCimag(a) * cuCimag(a));
+  double mag_b = sqrt(cuCreal(b) * cuCreal(b) + cuCimag(b) * cuCimag(b));
+  return mag_a > mag_b;
+}
+
+__forceinline__ __host__ __device__ bool operator>=(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  return a > b || a == b;
+}
+
+__forceinline__ __host__ __device__ bool operator<(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  return b > a;
+}
+
+__forceinline__ __host__ __device__ bool operator<=(
+    const cuDoubleComplex& a,
+    const cuDoubleComplex& b) {
+  return b > a || a == b;
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator+(const cuDoubleComplex& a, double b) {
+  return make_cuDoubleComplex(cuCreal(a) + b, cuCimag(a));
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator+(double a, const cuDoubleComplex& b) {
+  return make_cuDoubleComplex(a + cuCreal(b), cuCimag(b));
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator-(const cuDoubleComplex& a, double b) {
+  return make_cuDoubleComplex(cuCreal(a) - b, cuCimag(a));
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator-(double a, const cuDoubleComplex& b) {
+  return make_cuDoubleComplex(a - cuCreal(b), -cuCimag(b));
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator*(const cuDoubleComplex& a, double b) {
+  return make_cuDoubleComplex(cuCreal(a) * b, cuCimag(a) * b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator*(double a, const cuDoubleComplex& b) {
+  return make_cuDoubleComplex(a * cuCreal(b), a * cuCimag(b));
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator/(const cuDoubleComplex& a, double b) {
+  return make_cuDoubleComplex(cuCreal(a) / b, cuCimag(a) / b);
+}
+
+__forceinline__ __host__ __device__ cuDoubleComplex
+operator/(double a, const cuDoubleComplex& b) {
+  double denom = cuCreal(b) * cuCreal(b) + cuCimag(b) * cuCimag(b);
+  return make_cuDoubleComplex(
+      (a * cuCreal(b)) / denom, (-a * cuCimag(b)) / denom);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator+(const cuFloatComplex& a, const cuFloatComplex& b) {
+  return cuCaddf(a, b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator-(const cuFloatComplex& a, const cuFloatComplex& b) {
+  return cuCsubf(a, b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator*(const cuFloatComplex& a, const cuFloatComplex& b) {
+  return cuCmulf(a, b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator/(const cuFloatComplex& a, const cuFloatComplex& b) {
+  return cuCdivf(a, b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator%(const cuFloatComplex& a, const cuFloatComplex& b) {
+  float r = cuCrealf(a) - (floorf(cuCrealf(a) / cuCrealf(b)) * cuCrealf(b));
+  float i = cuCimagf(a) - (floorf(cuCimagf(a) / cuCimagf(b)) * cuCimagf(b));
+  return make_cuFloatComplex(r, i);
+}
+
+__forceinline__ __host__ __device__ bool operator==(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  return cuCrealf(a) == cuCrealf(b) && cuCimagf(a) == cuCimagf(b);
+}
+
+__forceinline__ __host__ __device__ bool operator!=(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  return !(a == b);
+}
+
+__forceinline__ __host__ __device__ bool operator>(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  float mag_a = sqrt(cuCrealf(a) * cuCrealf(a) + cuCimagf(a) * cuCimagf(a));
+  float mag_b = sqrt(cuCrealf(b) * cuCrealf(b) + cuCimagf(b) * cuCimagf(b));
+  return mag_a > mag_b;
+}
+
+__forceinline__ __host__ __device__ bool operator>=(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  return a > b || a == b;
+}
+
+__forceinline__ __host__ __device__ bool operator<(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  return b > a;
+}
+
+__forceinline__ __host__ __device__ bool operator<=(
+    const cuFloatComplex& a,
+    const cuFloatComplex& b) {
+  return b > a || a == b;
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator+(const cuFloatComplex& a, float b) {
+  return make_cuFloatComplex(cuCrealf(a) + b, cuCimagf(a));
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator+(float a, const cuFloatComplex& b) {
+  return make_cuFloatComplex(a + cuCrealf(b), cuCimagf(b));
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator-(const cuFloatComplex& a, float b) {
+  return make_cuFloatComplex(cuCrealf(a) - b, cuCimagf(a));
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator-(float a, const cuFloatComplex& b) {
+  return make_cuFloatComplex(a - cuCrealf(b), -cuCimagf(b));
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator*(const cuFloatComplex& a, float b) {
+  return make_cuFloatComplex(cuCrealf(a) * b, cuCimagf(a) * b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator*(float a, const cuFloatComplex& b) {
+  return make_cuFloatComplex(a * cuCrealf(b), a * cuCimagf(b));
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator/(const cuFloatComplex& a, float b) {
+  return make_cuFloatComplex(cuCrealf(a) / b, cuCimagf(a) / b);
+}
+
+__forceinline__ __host__ __device__ cuFloatComplex
+operator/(float a, const cuFloatComplex& b) {
+  float denom = cuCrealf(b) * cuCrealf(b) + cuCimagf(b) * cuCimagf(b);
+  return make_cuFloatComplex(
+      (a * cuCrealf(b)) / denom, (-a * cuCimagf(b)) / denom);
+}

mlx/backend/cuda/device/unary_ops.cuh

@@ -2,6 +2,8 @@
 
 #pragma once
 
+#include "mlx/backend/cuda/device/cexpf.cuh"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/device/fp16_math.cuh"
 #include "mlx/backend/cuda/device/utils.cuh"
 
@@ -14,6 +16,8 @@ struct Abs {
   __device__ T operator()(T x) {
     if constexpr (cuda::std::is_unsigned_v<T>) {
       return x;
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {sqrt(cuCrealf(x) * cuCrealf(x) + cuCimagf(x) * cuCimagf(x)), 0};
     } else {
       return abs(x);
     }
@@ -25,6 +29,8 @@ struct ArcCos {
   __device__ T operator()(T x) {
     return acos(x);
   }
+
+  __device__ cuComplex operator()(cuComplex x);
 };
 
 struct ArcCosh {
@@ -39,6 +45,8 @@ struct ArcSin {
   __device__ T operator()(T x) {
     return asin(x);
   }
+
+  __device__ cuComplex operator()(cuComplex x);
 };
 
 struct ArcSinh {
@@ -53,6 +61,8 @@ struct ArcTan {
   __device__ T operator()(T x) {
     return atan(x);
   }
+
+  __device__ cuComplex operator()(cuComplex x);
 };
 
 struct ArcTanh {
@@ -74,8 +84,6 @@ struct Ceil {
   __device__ T operator()(T x) {
     if constexpr (cuda::std::is_integral_v<T>) {
       return x;
-    } else if constexpr (is_complex_v<T>) {
-      return T{ceil(x.real()), ceil(x.imag())};
     } else {
       return ceil(x);
     }
@@ -83,23 +91,34 @@ struct Ceil {
 };
 
 struct Conjugate {
-  template <typename T>
-  __device__ complex_t<T> operator()(complex_t<T> x) {
-    return conj(x);
+  __device__ cuComplex operator()(cuComplex x) {
+    return {cuCrealf(x), -cuCimagf(x)};
   }
 };
 
 struct Cos {
   template <typename T>
   __device__ T operator()(T x) {
-    return cos(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {
+          cos(cuCrealf(x)) * cosh(cuCimagf(x)),
+          -sin(cuCrealf(x)) * sinh(cuCimagf(x))};
+    } else {
+      return cos(x);
+    }
   }
 };
 
 struct Cosh {
   template <typename T>
   __device__ T operator()(T x) {
-    return cosh(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {
+          cosh(cuCrealf(x)) * cos(cuCimagf(x)),
+          sinh(cuCrealf(x)) * sin(cuCimagf(x))};
+    } else {
+      return cosh(x);
+    }
   }
 };
 
@@ -132,7 +151,11 @@ struct ErfInv {
 struct Exp {
   template <typename T>
   __device__ T operator()(T x) {
-    return exp(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return detail::cexpf(x);
+    } else {
+      return exp(x);
+    }
   }
 };
 
@@ -154,8 +177,6 @@ struct Floor {
   __device__ T operator()(T x) {
     if constexpr (cuda::std::is_integral_v<T>) {
       return x;
-    } else if constexpr (is_complex_v<T>) {
-      return T{floor(x.real()), floor(x.imag())};
     } else {
       return floor(x);
     }
@@ -163,25 +184,30 @@ struct Floor {
 };
 
 struct Imag {
-  template <typename T>
-  __device__ auto operator()(complex_t<T> x) {
-    return x.imag();
+  __device__ float operator()(cuComplex x) {
+    return cuCimagf(x);
   }
 };
 
 struct Log {
   template <typename T>
   __device__ T operator()(T x) {
-    return log(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      auto r = log(cuCrealf(Abs{}(x)));
+      auto i = atan2f(cuCimagf(x), cuCrealf(x));
+      return {r, i};
+    } else {
+      return log(x);
+    }
   }
 };
 
 struct Log2 {
   template <typename T>
   __device__ T operator()(T x) {
-    if constexpr (is_complex_v<T>) {
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
       auto y = Log{}(x);
-      return {y.real() / CUDART_LN2_F, y.imag() / CUDART_LN2_F};
+      return {cuCrealf(y) / CUDART_LN2_F, cuCimagf(y) / CUDART_LN2_F};
     } else {
       return log2(x);
     }
@@ -191,17 +217,23 @@ struct Log2 {
 struct Log10 {
   template <typename T>
   __device__ T operator()(T x) {
-    return log10(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      auto y = Log{}(x);
+      return {cuCrealf(y) / CUDART_LNT_F, cuCimagf(y) / CUDART_LNT_F};
+      return y;
+    } else {
+      return log10(x);
+    }
   }
 };
 
 struct Log1p {
   template <typename T>
   __device__ T operator()(T z) {
-    if constexpr (is_complex_v<T>) {
-      float x = z.real();
-      float y = z.imag();
-      float zabs = Abs{}(z).real();
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      float x = cuCrealf(z);
+      float y = cuCimagf(z);
+      float zabs = cuCrealf(Abs{}(z));
       float theta = atan2f(y, x + 1);
       if (zabs < 0.5f) {
         float r = x * (2 + x) + y * y;
@@ -228,8 +260,8 @@ struct LogicalNot {
 struct Negative {
   template <typename T>
   __device__ T operator()(T x) {
-    if constexpr (is_complex_v<T>) {
-      return T{0, 0} - x;
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return 0 - x;
     } else {
       return -x;
     }
@@ -237,17 +269,16 @@ struct Negative {
 };
 
 struct Real {
-  template <typename T>
-  __device__ auto operator()(complex_t<T> x) {
-    return x.real();
+  __device__ float operator()(cuComplex x) {
+    return cuCrealf(x);
   }
 };
 
 struct Round {
   template <typename T>
   __device__ T operator()(T x) {
-    if constexpr (is_complex_v<T>) {
-      return {rint(x.real()), rint(x.imag())};
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {rint(cuCrealf(x)), rint(cuCimagf(x))};
     } else {
       return rint(x);
     }
@@ -267,8 +298,8 @@ struct Sign {
   __device__ T operator()(T x) {
     if constexpr (cuda::std::is_unsigned_v<T>) {
       return x != 0;
-    } else if constexpr (is_complex_v<T>) {
-      if (x.real() == 0 && x.imag() == 0) {
+    } else if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      if (cuCrealf(x) == 0 && cuCimagf(x) == 0) {
         return x;
       } else {
         return x / Abs()(x);
@@ -284,14 +315,26 @@ struct Sign {
 struct Sin {
   template <typename T>
   __device__ T operator()(T x) {
-    return sin(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {
+          sin(cuCrealf(x)) * cosh(cuCimagf(x)),
+          cos(cuCrealf(x)) * sinh(cuCimagf(x))};
+    } else {
+      return sin(x);
+    }
   }
 };
 
 struct Sinh {
   template <typename T>
   __device__ T operator()(T x) {
-    return sinh(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      return {
+          sinh(cuCrealf(x)) * cos(cuCimagf(x)),
+          cosh(cuCrealf(x)) * sin(cuCimagf(x))};
+    } else {
+      return sinh(x);
+    }
   }
 };
 
@@ -307,31 +350,77 @@ struct Sqrt {
   __device__ T operator()(T x) {
     return sqrt(x);
   }
+
+  __device__ cuComplex operator()(cuComplex x) {
+    auto xr = cuCrealf(x);
+    auto xi = cuCimagf(x);
+    if (xr == 0.0f && xi == 0.0f) {
+      return {0.0f, 0.0f};
+    }
+    auto r = cuCrealf(Abs{}(x));
+    auto a = sqrt((r + xr) / 2.0f);
+    auto b_abs = sqrt((r - xr) / 2.0f);
+    auto b = copysign(b_abs, xi);
+    return {a, b};
+  }
 };
 
 struct Rsqrt {
   template <typename T>
   __device__ T operator()(T x) {
-    if constexpr (is_complex_v<T>) {
-      return 1.0f / Sqrt{}(x);
-    } else {
-      return rsqrt(x);
-    }
+    return rsqrt(x);
+  }
+  __device__ cuComplex operator()(cuComplex x) {
+    return 1.0f / Sqrt{}(x);
   }
 };
 
 struct Tan {
   template <typename T>
   __device__ T operator()(T x) {
-    return tan(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      float tan_a = tan(cuCrealf(x));
+      float tanh_b = tanh(cuCimagf(x));
+      float t1 = tan_a * tanh_b;
+      float denom = 1. + t1 * t1;
+      return {(tan_a - tanh_b * t1) / denom, (tanh_b + tan_a * t1) / denom};
+    } else {
+      return tan(x);
+    }
   }
 };
 
 struct Tanh {
   template <typename T>
   __device__ T operator()(T x) {
-    return tanh(x);
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
+      float tanh_a = tanh(cuCrealf(x));
+      float tan_b = tan(cuCimagf(x));
+      float t1 = tanh_a * tan_b;
+      float denom = 1. + t1 * t1;
+      return {(tanh_a + tan_b * t1) / denom, (tan_b - tanh_a * t1) / denom};
+    } else {
+      return tanh(x);
+    }
   }
 };
 
+inline __device__ cuComplex ArcCos::operator()(cuComplex x) {
+  auto i = cuComplex{0.0, 1.0};
+  auto y = Log{}(x + i * Sqrt{}(1.0 - x * x));
+  return {cuCimagf(y), -cuCrealf(y)};
+};
+
+inline __device__ cuComplex ArcSin::operator()(cuComplex x) {
+  auto i = cuComplex{0.0f, 1.0f};
+  auto y = Log{}(i * x + Sqrt{}(1.0f - x * x));
+  return {cuCimagf(y), -cuCrealf(y)};
+};
+
+inline __device__ cuComplex ArcTan::operator()(cuComplex x) {
+  auto i = cuComplex{0.0f, 1.0f};
+  auto ix = i * x;
+  return (1.0f / cuComplex{0.0f, 2.0f}) * Log{}((1.0f + ix) / (1.0f - ix));
+};
+
 } // namespace mlx::core::cu

mlx/backend/cuda/device/utils.cuh

@@ -8,9 +8,9 @@
 
 #pragma once
 
-#include "mlx/backend/cuda/device/complex.cuh"
 #include "mlx/backend/cuda/device/config.h"
 
+#include <cuComplex.h>
 #include <cuda_bf16.h>
 #include <cuda_fp16.h>
 #include <cuda/std/array>
@@ -127,13 +127,13 @@ struct Limits<bool> {
   }
 };
 
-template <typename T>
-struct Limits<complex_t<T>> {
-  static constexpr __host__ __device__ complex_t<T> max() {
-    return {Limits<T>::max(), Limits<T>::max()};
+template <>
+struct Limits<cuComplex> {
+  static constexpr __host__ __device__ cuComplex max() {
+    return {Limits<float>::max(), Limits<float>::max()};
   }
-  static constexpr __host__ __device__ complex_t<T> min() {
-    return {Limits<T>::min(), Limits<T>::min()};
+  static constexpr __host__ __device__ cuComplex min() {
+    return {Limits<float>::min(), Limits<float>::min()};
   }
 };
 

mlx/backend/cuda/event.cu

@@ -90,6 +90,8 @@ bool CudaEvent::completed() const {
 // SharedEvent implementations
 ///////////////////////////////////////////////////////////////////////////////
 
+namespace {
+
 __host__ __device__ void event_wait(SharedEvent::Atomic* ac, uint64_t value) {
   uint64_t current;
   while ((current = ac->load()) < value) {
@@ -110,6 +112,8 @@ __global__ void event_signal_kernel(SharedEvent::Atomic* ac, uint64_t value) {
   event_signal(ac, value);
 }
 
+} // namespace
+
 SharedEvent::SharedEvent() {
   // Allocate cuda::atomic on managed memory.
   Atomic* ac;

mlx/backend/cuda/jit_module.cpp

@@ -13,7 +13,6 @@
 
 #include <fmt/format.h>
 #include <nvrtc.h>
-#include <unistd.h>
 
 namespace mlx::core::cu {
 
@@ -51,32 +50,6 @@ const std::string& cuda_home() {
   return home;
 }
 
-// Return the location of CCCL headers shipped with the distribution.
-const std::string& cccl_dir() {
-  static std::string dir = []() {
-    std::filesystem::path path;
-#if defined(MLX_CCCL_DIR)
-    // First search the install dir if defined.
-    path = MLX_CCCL_DIR;
-    if (std::filesystem::exists(path)) {
-      return path.string();
-    }
-#endif
-    // Then search dynamically from the dir of libmlx.so file.
-    path = current_binary_dir().parent_path() / "include" / "cccl";
-    if (std::filesystem::exists(path)) {
-      return path.string();
-    }
-    // Finally check the environment variable.
-    path = std::getenv("MLX_CCCL_DIR");
-    if (!path.empty() && std::filesystem::exists(path)) {
-      return path.string();
-    }
-    return std::string();
-  }();
-  return dir;
-}
-
 // Get the cache directory for storing compiled results.
 const std::filesystem::path& ptx_cache_dir() {
   static std::filesystem::path cache = []() -> std::filesystem::path {
@@ -137,8 +110,7 @@ void write_cached_ptx(
     const std::filesystem::path& cache_dir,
     const std::string& module_name,
     const std::vector<char>& ptx,
-    const std::vector<std::pair<std::string, std::string>>& ptx_kernels,
-    const std::string& source_code) {
+    const std::vector<std::pair<std::string, std::string>>& ptx_kernels) {
   if (cache_dir.empty()) {
     return;
   }
@@ -151,9 +123,6 @@ void write_cached_ptx(
   for (const auto& [name, mangled] : ptx_kernels) {
     txt_file << name << "\t" << mangled << std::endl;
   }
-
-  std::ofstream source_file(cache_dir / (module_name + ".cu"));
-  source_file << source_code;
 }
 
 // Return if |device|'s version is not newer than |major|.|minor| version.
@@ -192,8 +161,9 @@ constexpr const char* g_include_names[] = {
     INCLUDE_PREFIX "atomic_ops.cuh",
     INCLUDE_PREFIX "binary_ops.cuh",
     INCLUDE_PREFIX "cast_op.cuh",
+    INCLUDE_PREFIX "cexpf.cuh",
     INCLUDE_PREFIX "config.h",
-    INCLUDE_PREFIX "complex.cuh",
+    INCLUDE_PREFIX "cucomplex_math.cuh",
     INCLUDE_PREFIX "fp16_math.cuh",
     INCLUDE_PREFIX "indexing.cuh",
     INCLUDE_PREFIX "scatter_ops.cuh",
@@ -208,8 +178,9 @@ constexpr const char* g_headers[] = {
     jit_source_atomic_ops,
     jit_source_binary_ops,
     jit_source_cast_op,
+    jit_source_cexpf,
     jit_source_config,
-    jit_source_complex,
+    jit_source_cucomplex_math,
     jit_source_fp16_math,
     jit_source_indexing,
     jit_source_scatter_ops,
@@ -246,24 +217,16 @@ JitModule::JitModule(
     }
 
     // Compile program.
-    std::vector<const char*> args;
     bool use_sass = compiler_supports_device_sass(device);
     std::string compute = fmt::format(
         "--gpu-architecture={}_{}{}",
         use_sass ? "sm" : "compute",
         device.compute_capability_major(),
         device.compute_capability_minor());
-    args.push_back(compute.c_str());
-    std::string cccl_include = cccl_dir();
-    if (!cccl_include.empty()) {
-      cccl_include = fmt::format("--include-path={}", cccl_include);
-      args.push_back(cccl_include.c_str());
-    }
-    std::string cuda_include =
-        fmt::format("--include-path={}/include", cuda_home());
-    args.push_back(cuda_include.c_str());
+    std::string include = fmt::format("--include-path={}/include", cuda_home());
+    const char* args[] = {compute.c_str(), include.c_str()};
     nvrtcResult compile_result =
-        nvrtcCompileProgram(prog, args.size(), args.data());
+        nvrtcCompileProgram(prog, std::size(args), args);
     if (compile_result != NVRTC_SUCCESS) {
       size_t log_size;
       CHECK_NVRTC_ERROR(nvrtcGetProgramLogSize(prog, &log_size));
@@ -293,8 +256,7 @@ JitModule::JitModule(
     } else {
       CHECK_NVRTC_ERROR(nvrtcGetPTX(prog, ptx.data()));
     }
-    write_cached_ptx(
-        ptx_cache_dir(), module_name, ptx, ptx_kernels, source_code);
+    write_cached_ptx(ptx_cache_dir(), module_name, ptx, ptx_kernels);
   }
 
   // Load module.

mlx/backend/cuda/kernel_utils.cuh

@@ -11,6 +11,7 @@
 #include "mlx/array.h"
 #include "mlx/backend/cuda/device/utils.cuh"
 
+#include <cuComplex.h>
 #include <cuda.h>
 #include <cuda_bf16.h>
 #include <cuda_fp16.h>
@@ -78,7 +79,7 @@ struct CTypeToCudaType<bfloat16_t> {
 
 template <>
 struct CTypeToCudaType<complex64_t> {
-  using type = cu::complex64_t;
+  using type = cuComplex;
 };
 
 template <typename T>
@@ -90,14 +91,10 @@ inline constexpr bool is_floating_v =
     cuda::std::is_same_v<T, float> || cuda::std::is_same_v<T, double> ||
     cuda::std::is_same_v<T, float16_t> || cuda::std::is_same_v<T, bfloat16_t>;
 
-// Type traits for detecting complex numbers.
-template <typename T>
-inline constexpr bool is_complex_v = cuda::std::is_same_v<T, complex64_t> ||
-    cuda::std::is_same_v<T, complex128_t>;
-
 // Type traits for detecting complex or real floating point numbers.
 template <typename T>
-inline constexpr bool is_inexact_v = is_floating_v<T> || is_complex_v<T>;
+inline constexpr bool is_inexact_v =
+    is_floating_v<T> || cuda::std::is_same_v<T, complex64_t>;
 
 // Utility to copy data from vector to array in host.
 template <int NDIM = MAX_NDIM, typename T = int32_t>

mlx/backend/cuda/layer_norm.cu

@@ -237,7 +237,8 @@ void LayerNorm::eval_gpu(
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }
@@ -294,7 +295,9 @@ void LayerNormVJP::eval_gpu(
       return x;
     }
     copied = true;
-    return contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
+    return x_copy;
   };
   bool donate_x = inputs[0].is_donatable();
   bool donate_g = inputs[3].is_donatable();

mlx/backend/cuda/logsumexp.cu

@@ -108,7 +108,8 @@ void LogSumExp::eval_gpu(const std::vector<array>& inputs, array& out) {
     if (x.flags().contiguous && x.strides()[x.ndim() - 1] == 1) {
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       encoder.add_temporary(x_copy);
       return x_copy;
     }

mlx/backend/cuda/matmul.cpp

@@ -27,35 +27,6 @@ void check_cublas_error(const char* name, cublasStatus_t err) {
   }
 }
 
-struct CublasPreference {
-  CublasPreference(Device& device) {
-    // The recommended cublas workspace size is 4 MiB for pre-Hopper and 32 MiB
-    // for Hopper+:
-    // https://docs.nvidia.com/cuda/cublas/#cublassetworkspace
-    uint64_t MiB = 1024 * 1024;
-    uint64_t workspace_size =
-        device.compute_capability_major() >= 9 ? 32 * MiB : 4 * MiB;
-
-    CHECK_CUBLAS_ERROR(cublasLtMatmulPreferenceCreate(&pref_));
-    CHECK_CUBLAS_ERROR(cublasLtMatmulPreferenceSetAttribute(
-        pref_,
-        CUBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES,
-        &workspace_size,
-        sizeof(uint64_t)));
-  }
-
-  ~CublasPreference() {
-    CHECK_CUBLAS_ERROR(cublasLtMatmulPreferenceDestroy(pref_));
-  }
-
-  cublasLtMatmulPreference_t pref_{nullptr};
-};
-
-cublasLtMatmulPreference_t cublas_preference(Device& device) {
-  static CublasPreference pref(device);
-  return pref.pref_;
-}
-
 class MatMul {
  public:
   MatMul(
@@ -72,7 +43,7 @@ class MatMul {
       int32_t batch_count,
       int64_t a_batch_stride,
       int64_t b_batch_stride)
-      : handle_(device.lt_handle()), pref_(cublas_preference(device)) {
+      : handle_(device.lt_handle()) {
     heuristic_.state = CUBLAS_STATUS_NOT_INITIALIZED;
 
     auto scale_type = dtype_to_cuda_type(dtype);
@@ -106,6 +77,20 @@ class MatMul {
         type, b_rows, b_cols, b_transposed, ldb, batch_count, b_batch_stride);
     out_desc_ = create_matrix_layout(
         type, a_rows, b_cols, false, b_cols, batch_count, a_rows * b_cols);
+
+    // The recommended cublas workspace size is 4 MiB for pre-Hopper and 32 MiB
+    // for Hopper+:
+    // https://docs.nvidia.com/cuda/cublas/#cublassetworkspace
+    uint64_t MiB = 1024 * 1024;
+    uint64_t workspace_size =
+        device.compute_capability_major() >= 9 ? 32 * MiB : 4 * MiB;
+
+    CHECK_CUBLAS_ERROR(cublasLtMatmulPreferenceCreate(&pref_));
+    CHECK_CUBLAS_ERROR(cublasLtMatmulPreferenceSetAttribute(
+        pref_,
+        CUBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES,
+        &workspace_size,
+        sizeof(uint64_t)));
   }
 
   MatMul(
@@ -119,6 +104,7 @@ class MatMul {
       uint64_t b_rows,
       uint64_t b_cols,
       int64_t ldb,
+      bool c_transposed,
       int64_t ldc,
       int32_t batch_count,
       int64_t a_batch_stride,
@@ -140,15 +126,15 @@ class MatMul {
             b_batch_stride) {
     auto type = dtype_to_cuda_type(dtype);
     c_desc_ = create_matrix_layout(
-        type, a_rows, b_cols, false, ldc, batch_count, c_batch_stride);
+        type, a_rows, b_cols, c_transposed, ldc, batch_count, c_batch_stride);
   }
 
   ~MatMul() {
-    CHECK_CUBLAS_ERROR(cublasLtMatrixLayoutDestroy(a_desc_));
-    CHECK_CUBLAS_ERROR(cublasLtMatrixLayoutDestroy(b_desc_));
-    CHECK_CUBLAS_ERROR(cublasLtMatrixLayoutDestroy(c_desc_));
-    CHECK_CUBLAS_ERROR(cublasLtMatrixLayoutDestroy(out_desc_));
-    CHECK_CUBLAS_ERROR(cublasLtMatmulDescDestroy(matmul_desc_));
+    cublasLtMatrixLayoutDestroy(a_desc_);
+    cublasLtMatrixLayoutDestroy(b_desc_);
+    cublasLtMatrixLayoutDestroy(c_desc_);
+    cublasLtMatrixLayoutDestroy(out_desc_);
+    cublasLtMatmulDescDestroy(matmul_desc_);
   }
 
   void run(
@@ -273,9 +259,9 @@ class MatMul {
     return desc;
   }
 
-  cublasLtMatmulPreference_t pref_{nullptr};
   cublasLtHandle_t handle_{nullptr};
   cublasLtMatmulDesc_t matmul_desc_{nullptr};
+  cublasLtMatmulPreference_t pref_{nullptr};
   cublasLtMatrixLayout_t a_desc_{nullptr};
   cublasLtMatrixLayout_t b_desc_{nullptr};
   cublasLtMatrixLayout_t c_desc_{nullptr};
@@ -296,7 +282,8 @@ check_transpose(cu::CommandEncoder& enc, const Stream& s, const array& arr) {
   } else if (stx == 1 && sty == arr.shape(-2)) {
     return std::make_tuple(true, sty, arr);
   } else {
-    array arr_copy = contiguous_copy_gpu(arr, s);
+    array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
+    copy_gpu(arr, arr_copy, CopyType::General, s);
     enc.add_temporary(arr_copy);
     return std::make_tuple(false, arr.shape(-1), arr_copy);
   }
@@ -402,7 +389,9 @@ void AddMM::eval_gpu(const std::vector<array>& inputs, array& out) {
   assert(inputs.size() == 3);
   auto& a_pre = inputs[0];
   auto& b_pre = inputs[1];
-  auto c = inputs[2];
+  auto& c_pre = inputs[2];
+
+  out.set_data(allocator::malloc(out.nbytes()));
 
   /////////////////////////////////////////////////////////////////////////////
   // Init checks and prep
@@ -415,24 +404,7 @@ void AddMM::eval_gpu(const std::vector<array>& inputs, array& out) {
   // the arrays
   auto [a_transposed, lda, a] = check_transpose(encoder, s, a_pre);
   auto [b_transposed, ldb, b] = check_transpose(encoder, s, b_pre);
-
-  int64_t ldc;
-  {
-    auto stx = c.strides()[c.ndim() - 2];
-    auto sty = c.strides()[c.ndim() - 1];
-    if (sty == 1 && stx == c.shape(-1)) {
-      ldc = stx;
-      out.set_data(allocator::malloc(out.nbytes()));
-    } else if (sty == 1 && stx == 0) {
-      ldc = 0;
-      out.set_data(allocator::malloc(out.nbytes()));
-    } else {
-      // Copy C into out and set C to out
-      ldc = c.shape(-1);
-      copy_gpu(c, out, CopyType::General, s);
-      c = out;
-    }
-  }
+  auto [c_transposed, ldc, c] = check_transpose(encoder, s, c_pre);
 
   /////////////////////////////////////////////////////////////////////////////
   // Check and collapse batch dimensions
@@ -470,6 +442,7 @@ void AddMM::eval_gpu(const std::vector<array>& inputs, array& out) {
       K,
       N,
       ldb,
+      c_transposed,
       ldc,
       batch_shape.back(),
       a_batch_strides.back(),

mlx/backend/cuda/primitives.cu

@@ -91,6 +91,7 @@ NO_GPU_MULTI(Eigh)
 
 namespace fast {
 NO_GPU(ScaledDotProductAttention)
+NO_GPU_MULTI(AffineQuantize)
 NO_GPU_MULTI(CustomKernel)
 } // namespace fast
 

mlx/backend/cuda/quantized.cu

@@ -1,386 +0,0 @@
-// Copyright © 2025 Apple Inc.
-
-#include "mlx/backend/cuda/device.h"
-#include "mlx/backend/cuda/kernel_utils.cuh"
-#include "mlx/backend/gpu/copy.h"
-#include "mlx/dtype_utils.h"
-#include "mlx/fast_primitives.h"
-
-#include <cooperative_groups.h>
-#include <cooperative_groups/reduce.h>
-#include <nvtx3/nvtx3.hpp>
-
-namespace mlx::core {
-namespace cu {
-
-namespace cg = cooperative_groups;
-
-template <int bits, int wsize = 8>
-inline constexpr __device__ short get_pack_factor() {
-  return (bits == 3 || bits == 5) ? 8 : (bits == 6 ? 4 : wsize / bits);
-}
-
-template <int bits, int wsize = 8>
-inline constexpr __device__ short get_bytes_per_pack() {
-  constexpr int power_of_2_bits = (bits & (bits - 1)) == 0;
-  return power_of_2_bits ? (wsize / 8) : (bits == 5 ? 5 : 3);
-}
-
-template <typename T, int group_size, int bits>
-__global__ void
-affine_quantize(const T* w, uint8_t* out, T* scales, T* biases, size_t size) {
-  auto block_size = cg::this_thread_block().dim_threads();
-  auto block_idx = cg::this_thread_block().group_index();
-  auto idx_in_block = cg::this_thread_block().thread_index();
-
-  auto tidx = block_idx.x * block_size.x + idx_in_block.x;
-  auto tidy = block_idx.y * block_size.y + idx_in_block.y;
-
-  auto grid_dim_x =
-      cg::this_grid().dim_blocks().x * cg::this_grid().block_index().x;
-  constexpr float eps = 1e-7;
-  constexpr int simd_size = WARP_SIZE;
-  constexpr float n_bins = (1 << bits) - 1;
-  constexpr int pack_factor = get_pack_factor<bits, 8>();
-  constexpr int bytes_per_pack = get_bytes_per_pack<bits>();
-  constexpr int values_per_reduce = group_size / simd_size;
-  constexpr int writes_per_reduce = pack_factor / values_per_reduce;
-  constexpr int writes_per_pack =
-      writes_per_reduce > 1 ? 1 : values_per_reduce / pack_factor;
-  constexpr int power_of_2_bits = (bits & (bits - 1)) == 0;
-
-  size_t offset = tidx + grid_dim_x * size_t(tidy);
-  size_t in_index = offset * values_per_reduce;
-  if (in_index >= size) {
-    return;
-  }
-  size_t out_index = power_of_2_bits
-      ? offset * writes_per_pack
-      : offset * bytes_per_pack / writes_per_reduce;
-
-  float w_thread[values_per_reduce];
-  float w_min = Limits<float>::max();
-  float w_max = 0;
-
-#pragma clang loop unroll(full)
-  for (int i = 0; i < values_per_reduce; i++) {
-    float val = w[in_index + i];
-    w_thread[i] = val;
-    w_min = min(w_min, val);
-    w_max = max(w_max, val);
-  }
-
-  cg::greater<float> max_op;
-  cg::less<float> min_op;
-  auto warp = cg::tiled_partition<WARP_SIZE>(cg::this_thread_block());
-
-  w_min = cg::reduce(warp, w_min, min_op);
-  w_max = cg::reduce(warp, w_max, max_op);
-
-  float scale = max((w_max - w_min) / n_bins, eps);
-  bool side = abs(w_min) > abs(w_max);
-  scale = side ? scale : -scale;
-  float edge = side ? w_min : w_max;
-  float q0 = round(edge / scale);
-  bool at_zero = q0 == 0.0f;
-  scale = at_zero ? scale : edge / q0;
-  float bias = at_zero ? 0 : edge;
-
-  // Write out the scales and biases
-  size_t gindex = in_index / group_size;
-  if (in_index % group_size == 0) {
-    scales[gindex] = static_cast<T>(scale);
-    biases[gindex] = static_cast<T>(bias);
-  }
-
-  using OutType = std::conditional_t<bits == 5, uint64_t, uint32_t>;
-  OutType output = 0;
-
-#pragma clang loop unroll(full)
-  for (int i = 0; i < values_per_reduce; i++) {
-    uint8_t val = min(round((w_thread[i] - bias) / scale), n_bins);
-    if (bits == 8) {
-      output = val;
-    } else {
-      output |= val << (bits * (i % pack_factor));
-    }
-
-    if (pack_factor < values_per_reduce && i % pack_factor == pack_factor - 1) {
-      out[out_index + i / pack_factor] = output;
-      output = 0;
-    } else {
-#pragma clang loop unroll(full)
-      for (int j = 1; j < writes_per_reduce; j++) {
-        uint8_t sval = warp.shfl_down(val, j);
-        output |= static_cast<OutType>(sval)
-            << (bits * (j * values_per_reduce + i));
-      }
-    }
-  }
-  if constexpr (bits == 3 || bits == 6) {
-    if (in_index % pack_factor == 0 && out_index % bytes_per_pack == 0) {
-      out[out_index] = output & 0xff;
-      out[out_index + 1] = (output & 0xff00) >> 8;
-      out[out_index + 2] = (output & 0xff0000) >> 16;
-    }
-  } else if constexpr (bits == 5) {
-    if (in_index % pack_factor == 0 && out_index % bytes_per_pack == 0) {
-      out[out_index] = output & 0xff;
-      out[out_index + 1] = (output & 0xff00) >> 8;
-      out[out_index + 2] = (output & 0xff0000) >> 16;
-      out[out_index + 3] = (output & 0xff000000) >> 24;
-      out[out_index + 4] = (output & 0xff00000000) >> 32;
-    }
-  } else {
-    if constexpr (writes_per_reduce > 0) {
-      if (out_index % writes_per_reduce == 0) {
-        out[out_index / writes_per_reduce] = output;
-      }
-    }
-  }
-}
-
-template <typename T, int group_size, int bits>
-__global__ void affine_dequantize(
-    const uint8_t* w,
-    const T* scales,
-    const T* biases,
-    T* out,
-    size_t size) {
-  auto block_size = cg::this_thread_block().dim_threads();
-  auto block_idx = cg::this_thread_block().group_index();
-  auto idx_in_block = cg::this_thread_block().thread_index();
-
-  auto tidx = block_idx.x * block_size.x + idx_in_block.x;
-  auto tidy = block_idx.y * block_size.y + idx_in_block.y;
-
-  auto grid_dim_x =
-      cg::this_grid().dim_blocks().x * cg::this_grid().block_index().x;
-
-  constexpr int pack_factor = get_pack_factor<bits, 8>();
-  constexpr int bytes_per_pack = get_bytes_per_pack<bits>();
-
-  size_t offset = tidx + grid_dim_x * size_t(tidy);
-  size_t oindex = offset * pack_factor;
-
-  if (oindex >= size) {
-    return;
-  }
-
-  size_t gindex = oindex / group_size;
-  T scale = scales[gindex];
-  T bias = biases[gindex];
-  out += oindex;
-
-  if constexpr (bits == 3) {
-    w += offset * bytes_per_pack;
-    out[0] = static_cast<T>(w[0] & 0x7) * scale + bias;
-    out[1] = static_cast<T>((w[0] & 0x38) >> 3) * scale + bias;
-    out[2] = (static_cast<T>((w[0] & 0xc0) >> 6) +
-              static_cast<T>((w[1] & 0x1) << 2)) *
-            scale +
-        bias;
-    out[3] = static_cast<T>((w[1] & 0xe) >> 1) * scale + bias;
-    out[4] = static_cast<T>((w[1] & 0x70) >> 4) * scale + bias;
-    out[5] = (static_cast<T>((w[1] & 0x80) >> 7) +
-              static_cast<T>((w[2] & 0x3) << 1)) *
-            scale +
-        bias;
-    out[6] = static_cast<T>((w[2] & 0x1c) >> 2) * scale + bias;
-    out[7] = static_cast<T>((w[2] & 0xe0) >> 5) * scale + bias;
-  } else if constexpr (bits == 5) {
-    w += offset * bytes_per_pack;
-    out[0] = static_cast<T>(w[0] & 0x1f) * scale + bias;
-    out[1] = (static_cast<T>((w[0] & 0xe0) >> 5) +
-              static_cast<T>((w[1] & 0x3) << 3)) *
-            scale +
-        bias;
-    out[2] = static_cast<T>((w[1] & 0x7c) >> 2) * scale + bias;
-    out[3] = (static_cast<T>((w[1] & 0x80) >> 7) +
-              static_cast<T>((w[2] & 0xf) << 1)) *
-            scale +
-        bias;
-    out[4] = (static_cast<T>((w[2] & 0xf0) >> 4) +
-              static_cast<T>((w[3] & 0x1) << 4)) *
-            scale +
-        bias;
-    out[5] = static_cast<T>((w[3] & 0x3e) >> 1) * scale + bias;
-    out[6] = (static_cast<T>((w[3] & 0xc0) >> 6) +
-              static_cast<T>((w[4] & 0x7) << 2)) *
-            scale +
-        bias;
-    out[7] = static_cast<T>((w[4] & 0xf8) >> 3) * scale + bias;
-  } else if constexpr (bits == 6) {
-    w += offset * bytes_per_pack;
-    out[0] = static_cast<T>(w[0] & 0x3f) * scale + bias;
-    out[1] = (static_cast<T>((w[0] >> 6) & 0x03) +
-              static_cast<T>((w[1] & 0x0f) << 2)) *
-            scale +
-        bias;
-    out[2] = (static_cast<T>((w[1] >> 4) & 0x0f) +
-              static_cast<T>((w[2] & 0x03) << 4)) *
-            scale +
-        bias;
-    out[3] = static_cast<T>((w[2] >> 2) & 0x3f) * scale + bias;
-  } else {
-    uint val = w[offset];
-#pragma clang loop unroll(full)
-    for (int i = 0; i < pack_factor; i++) {
-      uint8_t d;
-      if (bits == 2) {
-        d = (val >> (bits * i)) & 0x03;
-      } else if (bits == 4) {
-        d = (val >> (bits * i)) & 0x0f;
-      } else if (bits == 8) {
-        d = val;
-      }
-      out[i] = scale * static_cast<T>(d) + bias;
-    }
-  }
-}
-
-} // namespace cu
-namespace {
-
-inline array ensure_row_contiguous(
-    const array& x,
-    cu::CommandEncoder& enc,
-    const Stream& s) {
-  if (!x.flags().row_contiguous) {
-    array x_copy = contiguous_copy_gpu(x, s);
-    enc.add_temporary(x_copy);
-    return x_copy;
-  } else {
-    return x;
-  }
-}
-
-} // namespace
-
-template <typename F>
-void dispatch_groups(int group_size, F&& f) {
-  switch (group_size) {
-    case 32:
-      f(std::integral_constant<int, 32>{});
-      break;
-    case 64:
-      f(std::integral_constant<int, 64>{});
-      break;
-    case 128:
-      f(std::integral_constant<int, 128>{});
-      break;
-  }
-}
-
-template <typename F>
-void dispatch_bits(int bits, F&& f) {
-  switch (bits) {
-    case 2:
-      f(std::integral_constant<int, 2>{});
-      break;
-    case 3:
-      f(std::integral_constant<int, 3>{});
-      break;
-    case 4:
-      f(std::integral_constant<int, 4>{});
-      break;
-    case 5:
-      f(std::integral_constant<int, 5>{});
-      break;
-    case 6:
-      f(std::integral_constant<int, 6>{});
-      break;
-    case 8:
-      f(std::integral_constant<int, 8>{});
-      break;
-  }
-}
-
-void fast::AffineQuantize::eval_gpu(
-    const std::vector<array>& inputs,
-    std::vector<array>& outputs) {
-  auto& w_pre = inputs[0];
-  auto& out = outputs[0];
-  out.set_data(allocator::malloc(out.nbytes()));
-
-  auto& s = stream();
-  auto& d = cu::device(s.device);
-  auto& enc = d.get_command_encoder(s);
-
-  auto w = ensure_row_contiguous(w_pre, enc, s);
-  enc.set_input_array(w);
-  if (dequantize_) {
-    auto scales = ensure_row_contiguous(inputs[1], enc, s);
-    auto biases = ensure_row_contiguous(inputs[2], enc, s);
-    enc.set_input_array(scales);
-    enc.set_input_array(biases);
-    enc.set_output_array(out);
-  } else {
-    auto& scales = outputs[1];
-    auto& biases = outputs[2];
-    scales.set_data(allocator::malloc(scales.nbytes()));
-    biases.set_data(allocator::malloc(biases.nbytes()));
-    enc.set_output_array(out);
-    enc.set_output_array(scales);
-    enc.set_output_array(biases);
-  }
-
-  auto dtype = dequantize_ ? outputs[0].dtype() : inputs[0].dtype();
-
-  // Treat uint32 as uint8 in kernel
-  int uint8_per_uint32 = 4;
-  int packs_per_int = (bits_ == 3 || bits_ == 5) ? 8
-      : bits_ == 6                               ? 4
-                                                 : 8 / bits_;
-  int per_thread = dequantize_ ? packs_per_int : group_size_ / WARP_SIZE;
-  size_t size =
-      dequantize_ ? out.size() / packs_per_int : w.size() / per_thread;
-
-  bool large = size > UINT_MAX;
-  auto grid_shape = w.shape();
-
-  if (dequantize_) {
-    grid_shape.back() *= uint8_per_uint32;
-  } else {
-    grid_shape.back() /= per_thread;
-  }
-
-  dispatch_float_types(dtype, "affine_quantize", [&](auto type_tag) {
-    dispatch_groups(group_size_, [&](auto group_size) {
-      dispatch_bits(bits_, [&](auto bits) {
-        using DataType = cuda_type_t<MLX_GET_TYPE(type_tag)>;
-        if (dequantize_) {
-          auto kernel =
-              cu::affine_dequantize<DataType, group_size.value, bits.value>;
-          auto [num_blocks, block_dims] =
-              get_launch_args(kernel, size, grid_shape, w.strides(), large);
-          enc.add_kernel_node(
-              kernel,
-              num_blocks,
-              block_dims,
-              w.data<uint8_t>(),
-              inputs[1].data<DataType>(),
-              inputs[2].data<DataType>(),
-              out.data<DataType>(),
-              out.size());
-        } else {
-          auto kernel =
-              cu::affine_quantize<DataType, group_size.value, bits.value>;
-          auto [num_blocks, block_dims] =
-              get_launch_args(kernel, size, grid_shape, w.strides(), large);
-          enc.add_kernel_node(
-              kernel,
-              num_blocks,
-              block_dims,
-              w.data<DataType>(),
-              out.data<uint8_t>(),
-              outputs[1].data<DataType>(),
-              outputs[2].data<DataType>(),
-              w.size());
-        }
-      });
-    });
-  });
-}
-
-} // namespace mlx::core

mlx/backend/cuda/reduce.cu

@@ -47,7 +47,8 @@ void Reduce::eval_gpu(const std::vector<array>& inputs, array& out) {
     }
   }
   if (plan.type == GeneralReduce || broadcasted || !in.flags().contiguous) {
-    array in_copy = contiguous_copy_gpu(in, s);
+    array in_copy(in.shape(), in.dtype(), nullptr, {});
+    copy_gpu(in, in_copy, CopyType::General, s);
     encoder.add_temporary(in_copy);
     in = in_copy;
     plan = get_reduction_plan(in, axes_);

mlx/backend/cuda/reduce/reduce.cuh

@@ -3,6 +3,7 @@
 #include <type_traits>
 
 #include "mlx/backend/common/reduce.h"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/kernel_utils.cuh"
 #include "mlx/backend/cuda/reduce/reduce_ops.cuh"
 #include "mlx/dtype_utils.h"

mlx/backend/cuda/reduce/reduce_ops.cuh

@@ -69,18 +69,6 @@ struct Prod {
 struct Min {
   template <typename T>
   __device__ __forceinline__ T operator()(T a, T b) {
-    if constexpr (is_complex_v<T>) {
-      if (isnan(a.real()) || isnan(a.imag())) {
-        return a;
-      }
-      if (isnan(b.real()) || isnan(b.imag())) {
-        return b;
-      }
-    } else if constexpr (!cuda::std::is_integral_v<T>) {
-      if (isnan(a) || isnan(b)) {
-        return cuda::std::numeric_limits<float>::quiet_NaN();
-      }
-    }
     return a < b ? a : b;
   }
 
@@ -93,18 +81,6 @@ struct Min {
 struct Max {
   template <typename T>
   __device__ __forceinline__ T operator()(T a, T b) {
-    if constexpr (is_complex_v<T>) {
-      if (isnan(a.real()) || isnan(a.imag())) {
-        return a;
-      }
-      if (isnan(b.real()) || isnan(b.imag())) {
-        return b;
-      }
-    } else if constexpr (!cuda::std::is_integral_v<T>) {
-      if (isnan(a) || isnan(b)) {
-        return cuda::std::numeric_limits<float>::quiet_NaN();
-      }
-    }
     return a > b ? a : b;
   }
 
@@ -175,7 +151,7 @@ struct ReduceInit<Or, T> {
 template <typename T>
 struct ReduceInit<Sum, T> {
   static constexpr __host__ __device__ auto value() {
-    if constexpr (is_complex_v<T>) {
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
       return T{0, 0};
     } else {
       return cast_to<typename ReduceResult<Sum, T>::type>(0);
@@ -186,7 +162,7 @@ struct ReduceInit<Sum, T> {
 template <typename T>
 struct ReduceInit<Prod, T> {
   static constexpr __host__ __device__ auto value() {
-    if constexpr (is_complex_v<T>) {
+    if constexpr (cuda::std::is_same_v<T, cuComplex>) {
       return T{1, 0};
     } else {
       return cast_to<typename ReduceResult<Prod, T>::type>(1);

mlx/backend/cuda/rms_norm.cu

@@ -206,7 +206,8 @@ void RMSNorm::eval_gpu(
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }
@@ -258,7 +259,9 @@ void RMSNormVJP::eval_gpu(
       return x;
     }
     copied = true;
-    return contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
+    return x_copy;
   };
   bool donate_x = inputs[0].is_donatable();
   bool donate_g = inputs[2].is_donatable();

mlx/backend/cuda/scan.cu

@@ -379,7 +379,9 @@ void Scan::eval_gpu(const std::vector<array>& inputs, array& out) {
           in.flags());
     }
   } else {
-    in = contiguous_copy_gpu(in, s);
+    array arr_copy(in.shape(), in.dtype(), nullptr, {});
+    copy_gpu(in, arr_copy, CopyType::General, s);
+    in = std::move(arr_copy);
     out.copy_shared_buffer(in);
   }
 

mlx/backend/cuda/softmax.cu

@@ -125,7 +125,8 @@ void Softmax::eval_gpu(const std::vector<array>& inputs, array& out) {
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }

mlx/backend/cuda/sort.cu

@@ -72,7 +72,8 @@ void gpu_sort(const Stream& s, array in, array& out_, int axis, bool argsort) {
   bool is_segmented_sort = in.flags().contiguous && in.strides()[axis] == 1;
   if (!is_segmented_sort) {
     array trans = swapaxes_in_eval(in, axis, last_dim);
-    in = contiguous_copy_gpu(trans, s);
+    in = array(trans.shape(), trans.dtype(), nullptr, {});
+    copy_gpu(trans, in, CopyType::General, s);
     encoder.add_temporary(in);
     out = array(allocator::malloc(out.nbytes()), in.shape(), out.dtype());
     encoder.add_temporary(out);

mlx/backend/cuda/unary.cu

@@ -2,6 +2,7 @@
 
 #include "mlx/backend/common/unary.h"
 #include "mlx/backend/cuda/device.h"
+#include "mlx/backend/cuda/device/cucomplex_math.cuh"
 #include "mlx/backend/cuda/device/unary_ops.cuh"
 #include "mlx/backend/cuda/iterators/general_iterator.cuh"
 #include "mlx/backend/cuda/kernel_utils.cuh"
@@ -70,10 +71,10 @@ constexpr bool supports_unary_op() {
         !std::is_same_v<In, bool>;
   }
   if (std::is_same_v<Op, Ceil> || std::is_same_v<Op, Floor>) {
-    return std::is_same_v<In, Out> && !mlx::core::is_complex_v<In>;
+    return std::is_same_v<In, Out> && !std::is_same_v<In, complex64_t>;
   }
   if (std::is_same_v<Op, Conjugate>) {
-    return std::is_same_v<In, Out> && mlx::core::is_complex_v<In>;
+    return std::is_same_v<In, Out> && std::is_same_v<In, complex64_t>;
   }
   if (std::is_same_v<Op, ArcCos> || std::is_same_v<Op, ArcSin> ||
       std::is_same_v<Op, ArcTan> || std::is_same_v<Op, Cos> ||
@@ -87,7 +88,7 @@ constexpr bool supports_unary_op() {
     return std::is_same_v<In, Out> && is_inexact_v<In>;
   }
   if (std::is_same_v<Op, Imag> || std::is_same_v<Op, Real>) {
-    return mlx::core::is_complex_v<In> && std::is_same_v<Out, float>;
+    return std::is_same_v<In, complex64_t> && std::is_same_v<Out, float>;
   }
   if (std::is_same_v<Op, LogicalNot>) {
     return std::is_same_v<In, Out> && std::is_same_v<In, bool>;
@@ -101,7 +102,7 @@ template <typename Op>
 void unary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   auto& in = inputs[0];
   if (in.size() == 0) {
@@ -177,17 +178,17 @@ template <typename Op>
 void unary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   set_unary_output_data(inputs[0], out);
   unary_op_gpu_inplace<Op>(inputs, out, op, s);
 }
 
-#define UNARY_GPU(func)                                               \
-  void func::eval_gpu(const std::vector<array>& inputs, array& out) { \
-    nvtx3::scoped_range r(#func "::eval_gpu");                        \
-    auto& s = out.primitive().stream();                               \
-    unary_op_gpu<cu::func>(inputs, out, name(), s);                   \
+#define UNARY_GPU(func)                                                 \
+  void func::eval_gpu(const std::vector<array>& inputs, array& out) {   \
+    nvtx3::scoped_range r(#func "::eval_gpu");                          \
+    auto& s = out.primitive().stream();                                 \
+    unary_op_gpu<cu::func>(inputs, out, get_primitive_string(this), s); \
   }
 
 UNARY_GPU(Abs)
@@ -223,15 +224,16 @@ UNARY_GPU(Tanh)
 void Log::eval_gpu(const std::vector<array>& inputs, array& out) {
   nvtx3::scoped_range r("Log::eval_gpu");
   auto& s = out.primitive().stream();
+  auto op = get_primitive_string(this);
   switch (base_) {
     case Base::e:
-      unary_op_gpu<cu::Log>(inputs, out, name(), s);
+      unary_op_gpu<cu::Log>(inputs, out, op, s);
       break;
     case Base::two:
-      unary_op_gpu<cu::Log2>(inputs, out, name(), s);
+      unary_op_gpu<cu::Log2>(inputs, out, op, s);
       break;
     case Base::ten:
-      unary_op_gpu<cu::Log10>(inputs, out, name(), s);
+      unary_op_gpu<cu::Log10>(inputs, out, op, s);
       break;
   }
 }
@@ -242,7 +244,7 @@ void Round::eval_gpu(const std::vector<array>& inputs, array& out) {
   const auto& in = inputs[0];
   auto& s = out.primitive().stream();
   if (issubdtype(in.dtype(), inexact)) {
-    unary_op_gpu<cu::Round>(inputs, out, name(), s);
+    unary_op_gpu<cu::Round>(inputs, out, get_primitive_string(this), s);
   } else {
     // No-op integer types
     out.copy_shared_buffer(in);

mlx/backend/cuda/utils.cpp

@@ -61,7 +61,7 @@ const char* dtype_to_cuda_type(const Dtype& dtype) {
     case float64:
       return "double";
     case complex64:
-      return "complex64_t";
+      return "cuComplex";
     default:
       return "unknown";
   }

mlx/backend/gpu/copy.cpp

@@ -46,10 +46,4 @@ void copy_gpu_inplace(
       in, out, in.shape(), i_strides, out.strides(), i_offset, 0, ctype, s);
 }
 
-array contiguous_copy_gpu(const array& arr, const Stream& s) {
-  array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
-  copy_gpu(arr, arr_copy, CopyType::General, s);
-  return arr_copy;
-}
-
 } // namespace mlx::core

mlx/backend/gpu/copy.h

@@ -43,7 +43,4 @@ void copy_gpu_inplace(
 // Fill the output with the scalar val
 void fill_gpu(const array& val, array& out, const Stream& s);
 
-// Return a contiguous array with same shape that copies the data of |arr|.
-array contiguous_copy_gpu(const array& arr, const Stream& s);
-
 } // namespace mlx::core

mlx/backend/metal/binary.cpp

@@ -7,20 +7,20 @@
 
 #define BINARY_GPU(func)                                              \
   void func::eval_gpu(const std::vector<array>& inputs, array& out) { \
-    binary_op_gpu(inputs, out, name());                               \
+    binary_op_gpu(inputs, out, get_primitive_string(this));           \
   }
 
 #define BINARY_GPU_MULTI(func)                                         \
   void func::eval_gpu(                                                 \
       const std::vector<array>& inputs, std::vector<array>& outputs) { \
-    binary_op_gpu(inputs, outputs, name());                            \
+    binary_op_gpu(inputs, outputs, get_primitive_string(this));        \
   }
 
 namespace mlx::core {
 
 std::string get_kernel_name(
     BinaryOpType bopt,
-    const char* op,
+    const std::string& op,
     const array& a,
     bool large,
     int ndim,
@@ -65,7 +65,7 @@ std::string get_kernel_name(
 void binary_op_gpu_inplace(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   auto& a = inputs[0];
   auto& b = inputs[1];
@@ -165,7 +165,7 @@ void binary_op_gpu_inplace(
 void binary_op_gpu(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   assert(inputs.size() == 2);
   auto& a = inputs[0];
@@ -179,7 +179,7 @@ void binary_op_gpu(
 void binary_op_gpu(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op) {
+    const std::string& op) {
   auto& s = outputs[0].primitive().stream();
   binary_op_gpu(inputs, outputs, op, s);
 }
@@ -187,7 +187,7 @@ void binary_op_gpu(
 void binary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   std::vector<array> outputs = {out};
   binary_op_gpu_inplace(inputs, outputs, op, s);
@@ -196,7 +196,7 @@ void binary_op_gpu_inplace(
 void binary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s) {
   assert(inputs.size() == 2);
   auto& a = inputs[0];
@@ -209,7 +209,7 @@ void binary_op_gpu(
 void binary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op) {
+    const std::string& op) {
   auto& s = out.primitive().stream();
   binary_op_gpu(inputs, out, op, s);
 }
@@ -237,19 +237,19 @@ BINARY_GPU(Subtract)
 void BitwiseBinary::eval_gpu(const std::vector<array>& inputs, array& out) {
   switch (op_) {
     case BitwiseBinary::And:
-      binary_op_gpu(inputs, out, name());
+      binary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case BitwiseBinary::Or:
-      binary_op_gpu(inputs, out, name());
+      binary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case BitwiseBinary::Xor:
-      binary_op_gpu(inputs, out, name());
+      binary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case BitwiseBinary::LeftShift:
-      binary_op_gpu(inputs, out, name());
+      binary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case BitwiseBinary::RightShift:
-      binary_op_gpu(inputs, out, name());
+      binary_op_gpu(inputs, out, get_primitive_string(this));
       break;
   }
 }

mlx/backend/metal/binary.h

@@ -9,25 +9,25 @@ namespace mlx::core {
 void binary_op_gpu(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    const std::string& op,
     const Stream& s);
 
 void binary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s);
 
 void binary_op_gpu_inplace(
     const std::vector<array>& inputs,
     std::vector<array>& outputs,
-    const char* op,
+    const std::string& op,
     const Stream& s);
 
 void binary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string& op,
     const Stream& s);
 
 } // namespace mlx::core

mlx/backend/metal/compiled.cpp

@@ -212,7 +212,9 @@ inline void build_kernel(
           get_type_string(x.dtype()),
           namer.get_name(x.inputs()[0]));
     } else {
-      os += x.primitive().name();
+      std::ostringstream ss;
+      x.primitive().print(ss);
+      os += ss.str();
       os += "()(";
       for (int i = 0; i < x.inputs().size() - 1; i++) {
         os += fmt::format("tmp_{0}, ", namer.get_name(x.inputs()[i]));

mlx/backend/metal/conv.cpp

@@ -149,7 +149,8 @@ void explicit_gemm_conv_group_ND_gpu(
       wt, {wt.strides(0), 1, C_per_group}, wt.flags(), wt.size());
 
   // Materialize
-  array wt_transpose = contiguous_copy_gpu(wt_view, s);
+  auto wt_transpose = array(wt_view.shape(), wt_view.dtype(), nullptr, {});
+  copy_gpu(wt_view, wt_transpose, CopyType::General, s);
 
   // Perform gemm
   std::vector<array> copies = {in_unfolded, wt_transpose};
@@ -960,12 +961,16 @@ void Convolution::eval_gpu(const std::vector<array>& inputs, array& out) {
   auto in = inputs[0];
   auto wt = inputs[1];
   if (!in.flags().row_contiguous) {
-    in = contiguous_copy_gpu(in, s);
-    copies.push_back(in);
+    array arr_copy(in.shape(), in.dtype(), nullptr, {});
+    copy_gpu(in, arr_copy, CopyType::General, s);
+    copies.push_back(arr_copy);
+    in = arr_copy;
   }
   if (!wt.flags().row_contiguous) {
-    wt = contiguous_copy_gpu(wt, s);
-    copies.push_back(wt);
+    array arr_copy(wt.shape(), wt.dtype(), nullptr, {});
+    copy_gpu(wt, arr_copy, CopyType::General, s);
+    copies.push_back(arr_copy);
+    wt = arr_copy;
   }
 
   // 3D conv

mlx/backend/metal/copy.cpp

@@ -86,7 +86,7 @@ void copy_gpu_inplace(
     }
   } else {
     work_per_thread = get_work_per_thread(out.dtype(), out.data_size());
-    if (!large && work_per_thread > 1) {
+    if (work_per_thread > 1) {
       kernel_name += "n";
     }
   }

mlx/backend/metal/device.cpp

@@ -1,18 +1,20 @@
 // Copyright © 2023-2024 Apple Inc.
 
 #include <cstdlib>
+#include <filesystem>
 #include <sstream>
 
 #define NS_PRIVATE_IMPLEMENTATION
 #define CA_PRIVATE_IMPLEMENTATION
 #define MTL_PRIVATE_IMPLEMENTATION
 
-#include "mlx/backend/common/utils.h"
 #include "mlx/backend/metal/device.h"
 #include "mlx/backend/metal/metal.h"
 #include "mlx/backend/metal/utils.h"
 #include "mlx/utils.h"
 
+namespace fs = std::filesystem;
+
 namespace mlx::core::metal {
 
 namespace {
@@ -78,7 +80,12 @@ MTL::Library* try_load_bundle(
 std::pair<MTL::Library*, NS::Error*> load_colocated_library(
     MTL::Device* device,
     const std::string& relative_path) {
-  auto path = current_binary_dir() / relative_path;
+  std::string binary_dir = get_binary_directory();
+  if (binary_dir.size() == 0) {
+    return {nullptr, nullptr};
+  }
+
+  auto path = fs::path(binary_dir) / relative_path;
   if (!path.has_extension()) {
     path.replace_extension(".metallib");
   }
@@ -190,7 +197,7 @@ MTL::Library* load_library(
 
   std::ostringstream msg;
   msg << "Failed to load the metallib " << lib_name << ".metallib. "
-      << "We attempted to load it from <" << current_binary_dir() << "/"
+      << "We attempted to load it from <" << get_binary_directory() << "/"
       << lib_name << ".metallib" << ">";
 #ifdef SWIFTPM_BUNDLE
   msg << " and from the Swift PM bundle.";

mlx/backend/metal/device.h

@@ -3,6 +3,8 @@
 #pragma once
 
 #include <Metal/Metal.hpp>
+#include <dlfcn.h>
+#include <filesystem>
 #include <functional>
 #include <mutex>
 #include <shared_mutex>
@@ -13,8 +15,22 @@
 #include "mlx/array.h"
 #include "mlx/device.h"
 
+namespace fs = std::filesystem;
+
 namespace mlx::core::metal {
 
+// Note, this function must be left inline in a header so that it is not
+// dynamically linked.
+inline std::string get_binary_directory() {
+  Dl_info info;
+  std::string directory;
+  int success = dladdr((void*)get_binary_directory, &info);
+  if (success) {
+    directory = fs::path(info.dli_fname).remove_filename().c_str();
+  }
+  return directory;
+}
+
 using MTLFCList =
     std::vector<std::tuple<const void*, MTL::DataType, NS::UInteger>>;
 

mlx/backend/metal/jit_kernels.cpp

@@ -8,6 +8,12 @@ using namespace fmt::literals;
 
 namespace mlx::core {
 
+std::string op_name(const array& arr) {
+  std::ostringstream op_t;
+  arr.primitive().print(op_t);
+  return op_t.str();
+}
+
 MTL::ComputePipelineState* get_arange_kernel(
     metal::Device& d,
     const std::string& kernel_name,
@@ -27,7 +33,7 @@ MTL::ComputePipelineState* get_unary_kernel(
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op) {
+    const std::string op) {
   std::string lib_name = kernel_name.substr(kernel_name.find("_") + 1);
   auto lib = d.get_library(lib_name, [&]() {
     auto in_t = get_type_string(in_type);
@@ -52,10 +58,10 @@ MTL::ComputePipelineState* get_unary_kernel(
 }
 
 void append_binary_kernels(
-    const std::string& lib_name,
+    const std::string lib_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op,
+    const std::string op,
     std::string& kernel_source) {
   const std::array<std::pair<std::string, std::string>, 7> kernel_types = {{
       {"ss", "binary_ss"},
@@ -106,7 +112,7 @@ MTL::ComputePipelineState* get_binary_kernel(
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op) {
+    const std::string op) {
   std::string lib_name = kernel_name.substr(kernel_name.find("_") + 1);
   auto lib = d.get_library(lib_name, [&]() {
     std::string kernel_source;
@@ -123,7 +129,7 @@ MTL::ComputePipelineState* get_binary_two_kernel(
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op) {
+    const std::string op) {
   std::string lib_name = kernel_name.substr(kernel_name.find("_") + 1);
   auto lib = d.get_library(lib_name, [&]() {
     std::string kernel_source = metal::utils();
@@ -138,7 +144,7 @@ MTL::ComputePipelineState* get_ternary_kernel(
     metal::Device& d,
     const std::string& kernel_name,
     Dtype type,
-    const char* op) {
+    const std::string op) {
   std::string lib_name = kernel_name.substr(kernel_name.find("_") + 1);
   auto lib = d.get_library(lib_name, [&]() {
     auto t_str = get_type_string(type);

mlx/backend/metal/kernels.h

@@ -19,27 +19,27 @@ MTL::ComputePipelineState* get_unary_kernel(
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op);
+    const std::string op);
 
 MTL::ComputePipelineState* get_binary_kernel(
     metal::Device& d,
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op);
+    const std::string op);
 
 MTL::ComputePipelineState* get_binary_two_kernel(
     metal::Device& d,
     const std::string& kernel_name,
     Dtype in_type,
     Dtype out_type,
-    const char* op);
+    const std::string op);
 
 MTL::ComputePipelineState* get_ternary_kernel(
     metal::Device& d,
     const std::string& kernel_name,
     Dtype type,
-    const char* op);
+    const std::string op);
 
 MTL::ComputePipelineState* get_copy_kernel(
     metal::Device& d,
@@ -257,10 +257,8 @@ MTL::ComputePipelineState* get_gather_qmm_kernel(
 
 // Create a GPU kernel template definition for JIT compilation
 template <typename... Args>
-std::string get_template_definition(
-    std::string_view name,
-    std::string_view func,
-    Args... args) {
+std::string
+get_template_definition(std::string name, std::string func, Args... args) {
   std::ostringstream s;
   s << func << "<";
   bool first = true;

mlx/backend/metal/logsumexp.cpp

@@ -25,7 +25,8 @@ void LogSumExp::eval_gpu(const std::vector<array>& inputs, array& out) {
     if (x.flags().contiguous && x.strides()[x.ndim() - 1] == 1) {
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       d.add_temporary(x_copy, s.index);
       return x_copy;
     }

mlx/backend/metal/matmul.cpp

@@ -33,7 +33,8 @@ std::tuple<bool, int64_t, array> check_transpose(
   } else if (stx == 1 && (!is_vector || sty == arr.shape(-2))) {
     return std::make_tuple(true, sty, arr);
   } else {
-    array arr_copy = contiguous_copy_gpu(arr, s);
+    array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
+    copy_gpu(arr, arr_copy, CopyType::General, s);
     copies.push_back(arr_copy);
     return std::make_tuple(false, arr.shape(-1), arr_copy);
   }
@@ -42,7 +43,8 @@ std::tuple<bool, int64_t, array> check_transpose(
 inline array
 ensure_row_contiguous(const array& x, metal::Device& d, const Stream& s) {
   if (!x.flags().row_contiguous) {
-    array x_copy = contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     d.add_temporary(x_copy, s.index);
     return x_copy;
   } else {
@@ -73,7 +75,8 @@ ensure_batch_contiguous(const array& x, metal::Device& d, const Stream& s) {
     }
   }
 
-  array x_copy = contiguous_copy_gpu(x, s);
+  array x_copy(x.shape(), x.dtype(), nullptr, {});
+  copy_gpu(x, x_copy, CopyType::General, s);
   d.add_temporary(x_copy, s.index);
   return std::make_tuple(false, x_copy.strides()[x_copy.ndim() - 2], x_copy);
 }
@@ -1891,7 +1894,8 @@ void segmented_mm(
       return std::make_tuple(false, x);
     }
 
-    array x_copy = contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     d.add_temporary(x_copy, s.index);
     return std::make_tuple(true, x_copy);
   };

mlx/backend/metal/nojit_kernels.cpp

@@ -18,7 +18,7 @@ MTL::ComputePipelineState* get_unary_kernel(
     const std::string& kernel_name,
     Dtype,
     Dtype,
-    const char*) {
+    const std::string) {
   return d.get_kernel(kernel_name);
 }
 
@@ -27,7 +27,7 @@ MTL::ComputePipelineState* get_binary_kernel(
     const std::string& kernel_name,
     Dtype,
     Dtype,
-    const char*) {
+    const std::string) {
   return d.get_kernel(kernel_name);
 }
 
@@ -36,7 +36,7 @@ MTL::ComputePipelineState* get_binary_two_kernel(
     const std::string& kernel_name,
     Dtype,
     Dtype,
-    const char*) {
+    const std::string) {
   return d.get_kernel(kernel_name);
 }
 
@@ -44,7 +44,7 @@ MTL::ComputePipelineState* get_ternary_kernel(
     metal::Device& d,
     const std::string& kernel_name,
     Dtype,
-    const char*) {
+    const std::string) {
   return d.get_kernel(kernel_name);
 }
 

mlx/backend/metal/normalization.cpp

@@ -40,7 +40,8 @@ void RMSNorm::eval_gpu(
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }
@@ -106,7 +107,9 @@ void RMSNormVJP::eval_gpu(
     if (x.flags().row_contiguous) {
       return {x, false};
     }
-    array x_copy = contiguous_copy_gpu(x, s);
+
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     return {x_copy, true};
   };
   bool donate_x = inputs[0].is_donatable();
@@ -238,7 +241,8 @@ void LayerNorm::eval_gpu(
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }
@@ -315,7 +319,8 @@ void LayerNormVJP::eval_gpu(
     if (x.flags().row_contiguous) {
       return {x, false};
     }
-    array x_copy = contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     return {x_copy, true};
   };
   bool donate_x = inputs[0].is_donatable();

mlx/backend/metal/quantized.cpp

@@ -20,7 +20,8 @@ namespace {
 inline array
 ensure_row_contiguous(const array& x, metal::Device& d, const Stream& s) {
   if (!x.flags().row_contiguous) {
-    array x_copy = contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     d.add_temporary(x_copy, s.index);
     return x_copy;
   } else {
@@ -37,7 +38,8 @@ inline array ensure_row_contiguous_matrix(
   if (stride_0 == x.shape(-1) && stride_1 == 1) {
     return x;
   } else {
-    array x_copy = contiguous_copy_gpu(x, s);
+    array x_copy(x.shape(), x.dtype(), nullptr, {});
+    copy_gpu(x, x_copy, CopyType::General, s);
     d.add_temporary(x_copy, s.index);
     return x_copy;
   }

mlx/backend/metal/reduce.cpp

@@ -989,7 +989,8 @@ void Reduce::eval_gpu(const std::vector<array>& inputs, array& out) {
     //       input for the axes with stride smaller than the minimum reduction
     //       stride.
     if (plan.type == GeneralReduce) {
-      array in_copy = contiguous_copy_gpu(in, s);
+      array in_copy(in.shape(), in.dtype(), nullptr, {});
+      copy_gpu(in, in_copy, CopyType::General, s);
       d.add_temporary(in_copy, s.index);
       in = in_copy;
       plan = get_reduction_plan(in, axes_);

mlx/backend/metal/scaled_dot_product_attention.cpp

@@ -398,7 +398,8 @@ void ScaledDotProductAttention::eval_gpu(
   auto copy_unless = [&copies, &s](
                          auto predicate, const array& arr) -> const array& {
     if (!predicate(arr)) {
-      array arr_copy = contiguous_copy_gpu(arr, s);
+      array arr_copy(arr.shape(), arr.dtype(), nullptr, {});
+      copy_gpu(arr, arr_copy, CopyType::General, s);
       copies.push_back(std::move(arr_copy));
       return copies.back();
     } else {

mlx/backend/metal/scan.cpp

@@ -30,7 +30,9 @@ void Scan::eval_gpu(const std::vector<array>& inputs, array& out) {
           in.flags());
     }
   } else {
-    in = contiguous_copy_gpu(in, s);
+    array arr_copy(in.shape(), in.dtype(), nullptr, {});
+    copy_gpu(in, arr_copy, CopyType::General, s);
+    in = std::move(arr_copy);
     out.copy_shared_buffer(in);
   }
 

mlx/backend/metal/softmax.cpp

@@ -35,7 +35,8 @@ void Softmax::eval_gpu(const std::vector<array>& inputs, array& out) {
       }
       return x;
     } else {
-      array x_copy = contiguous_copy_gpu(x, s);
+      auto x_copy = array(x.shape(), x.dtype(), nullptr, {});
+      copy_gpu(x, x_copy, CopyType::General, s);
       out.copy_shared_buffer(x_copy);
       return x_copy;
     }

mlx/backend/metal/ternary.cpp

@@ -11,7 +11,7 @@ namespace mlx::core {
 void ternary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s) {
   assert(inputs.size() == 3);
   auto& a = inputs[0];
@@ -128,7 +128,7 @@ void ternary_op_gpu_inplace(
 void ternary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s) {
   auto& a = inputs[0];
   auto& b = inputs[1];
@@ -141,13 +141,13 @@ void ternary_op_gpu(
 void ternary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op) {
+    const std::string op) {
   auto& s = out.primitive().stream();
   ternary_op_gpu(inputs, out, op, s);
 }
 
 void Select::eval_gpu(const std::vector<array>& inputs, array& out) {
-  ternary_op_gpu(inputs, out, name());
+  ternary_op_gpu(inputs, out, get_primitive_string(this));
 }
 
 } // namespace mlx::core

mlx/backend/metal/ternary.h

@@ -9,13 +9,13 @@ namespace mlx::core {
 void ternary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s);
 
 void ternary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s);
 
 } // namespace mlx::core

mlx/backend/metal/unary.cpp

@@ -8,7 +8,7 @@
 
 #define UNARY_GPU(func)                                               \
   void func::eval_gpu(const std::vector<array>& inputs, array& out) { \
-    unary_op_gpu(inputs, out, name());                                \
+    unary_op_gpu(inputs, out, get_primitive_string(this));            \
   }
 
 namespace mlx::core {
@@ -16,7 +16,7 @@ namespace mlx::core {
 void unary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s) {
   auto& in = inputs[0];
   bool contig = in.flags().contiguous;
@@ -98,7 +98,7 @@ void unary_op_gpu_inplace(
 void unary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s) {
   set_unary_output_data(inputs[0], out);
   unary_op_gpu_inplace(inputs, out, op, s);
@@ -107,7 +107,7 @@ void unary_op_gpu(
 void unary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op) {
+    const std::string op) {
   auto& s = out.primitive().stream();
   unary_op_gpu(inputs, out, op, s);
 }
@@ -146,13 +146,13 @@ UNARY_GPU(Tanh)
 void Log::eval_gpu(const std::vector<array>& inputs, array& out) {
   switch (base_) {
     case Base::e:
-      unary_op_gpu(inputs, out, name());
+      unary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case Base::two:
-      unary_op_gpu(inputs, out, name());
+      unary_op_gpu(inputs, out, get_primitive_string(this));
       break;
     case Base::ten:
-      unary_op_gpu(inputs, out, name());
+      unary_op_gpu(inputs, out, get_primitive_string(this));
       break;
   }
 }
@@ -161,7 +161,7 @@ void Round::eval_gpu(const std::vector<array>& inputs, array& out) {
   assert(inputs.size() == 1);
   const auto& in = inputs[0];
   if (issubdtype(in.dtype(), inexact)) {
-    unary_op_gpu(inputs, out, name());
+    unary_op_gpu(inputs, out, get_primitive_string(this));
   } else {
     // No-op integer types
     out.copy_shared_buffer(in);

mlx/backend/metal/unary.h

@@ -9,13 +9,13 @@ namespace mlx::core {
 void unary_op_gpu(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s);
 
 void unary_op_gpu_inplace(
     const std::vector<array>& inputs,
     array& out,
-    const char* op,
+    const std::string op,
     const Stream& s);
 
 } // namespace mlx::core

mlx/backend/metal/utils.h

@@ -40,7 +40,7 @@ inline void debug_set_primitive_buffer_label(
   if (auto cbuf_label = command_buffer->label(); cbuf_label) {
     label << cbuf_label->utf8String();
   }
-  label << primitive.name();
+  primitive.print(label);
   command_buffer->setLabel(make_string(label));
 #endif
 }

mlx/compile.cpp

@@ -107,7 +107,7 @@ Compiled::Compiled(
     // name and type of output
     os << namer.get_name(a) << kindof(a.dtype()) << a.itemsize();
     // computation performed
-    os << a.primitive().name();
+    a.primitive().print(os);
     // name of inputs to the function
     for (auto& inp : a.inputs()) {
       os << namer.get_name(inp);
@@ -170,16 +170,11 @@ bool Compiled::is_equivalent(const Primitive& other) const {
       });
 }
 
-const char* Compiled::name() const {
-  if (name_.empty()) {
-    std::ostringstream os;
-    os << "Compiled";
-    for (auto& a : tape_) {
-      os << a.primitive().name();
-    }
-    name_ = os.str();
+void Compiled::print(std::ostream& os) {
+  os << "Compiled";
+  for (auto& a : tape_) {
+    a.primitive().print(os);
   }
-  return name_.c_str();
 }
 
 std::vector<Shape> Compiled::output_shapes(const std::vector<array>& inputs) {

mlx/distributed/primitives.h

@@ -45,22 +45,27 @@ class AllReduce : public DistPrimitive {
       const std::vector<int>& argnums,
       const std::vector<array>& outputs) override;
 
-  const char* name() const override {
+  void print(std::ostream& os) override {
     switch (reduce_type_) {
       case And:
-        return "And AllReduce";
+        os << "And";
       case Or:
-        return "Or AllReduce";
+        os << "And";
+        break;
       case Sum:
-        return "Sum AllReduce";
+        os << "Sum";
+        break;
       case Prod:
-        return "Prod AllReduce";
+        os << "Prod";
+        break;
       case Min:
-        return "Min AllReduce";
+        os << "Min";
+        break;
       case Max:
-        return "Max AllReduce";
+        os << "Max";
+        break;
     }
-    return "<unknwon AllReduce>";
+    os << " AllReduce";
   }
 
  private:
@@ -89,7 +94,7 @@ class AllGather : public DistPrimitive {
       const std::vector<int>& argnums,
       const std::vector<array>& outputs) override;
 
-  DEFINE_NAME(AllGather);
+  DEFINE_PRINT(AllGather);
 };
 
 class Send : public DistPrimitive {
@@ -105,7 +110,7 @@ class Send : public DistPrimitive {
       const std::vector<array>& inputs,
       const std::vector<int>& axes) override;
 
-  DEFINE_NAME(Send);
+  DEFINE_PRINT(Send);
 
  private:
   int dst_;
@@ -121,7 +126,7 @@ class Recv : public DistPrimitive {
   void eval_gpu(const std::vector<array>& inputs, std::vector<array>& outputs)
       override;
 
-  DEFINE_NAME(Recv);
+  DEFINE_PRINT(Recv);
 
  private:
   int src_;

mlx/distributed/ring/ring.cpp

@@ -22,78 +22,20 @@
 #include "mlx/backend/cpu/encoder.h"
 #include "mlx/distributed/distributed.h"
 #include "mlx/distributed/distributed_impl.h"
+#include "mlx/dtype_utils.h"
 #include "mlx/threadpool.h"
 
 #ifndef SOL_TCP
 #define SOL_TCP IPPROTO_TCP
 #endif
 
-#define SWITCH_TYPE(x, ...)  \
-  switch ((x).dtype()) {     \
-    case bool_: {            \
-      using T = bool;        \
-      __VA_ARGS__;           \
-    } break;                 \
-    case int8: {             \
-      using T = int8_t;      \
-      __VA_ARGS__;           \
-    } break;                 \
-    case int16: {            \
-      using T = int16_t;     \
-      __VA_ARGS__;           \
-    } break;                 \
-    case int32: {            \
-      using T = int32_t;     \
-      __VA_ARGS__;           \
-    } break;                 \
-    case int64: {            \
-      using T = int64_t;     \
-      __VA_ARGS__;           \
-    } break;                 \
-    case uint8: {            \
-      using T = uint8_t;     \
-      __VA_ARGS__;           \
-    } break;                 \
-    case uint16: {           \
-      using T = uint16_t;    \
-      __VA_ARGS__;           \
-    } break;                 \
-    case uint32: {           \
-      using T = uint32_t;    \
-      __VA_ARGS__;           \
-    } break;                 \
-    case uint64: {           \
-      using T = uint64_t;    \
-      __VA_ARGS__;           \
-    } break;                 \
-    case bfloat16: {         \
-      using T = bfloat16_t;  \
-      __VA_ARGS__;           \
-    } break;                 \
-    case float16: {          \
-      using T = float16_t;   \
-      __VA_ARGS__;           \
-    } break;                 \
-    case float32: {          \
-      using T = float;       \
-      __VA_ARGS__;           \
-    } break;                 \
-    case float64: {          \
-      using T = double;      \
-      __VA_ARGS__;           \
-    } break;                 \
-    case complex64: {        \
-      using T = complex64_t; \
-      __VA_ARGS__;           \
-    } break;                 \
-  }
-
 namespace mlx::core::distributed::ring {
 
 constexpr const size_t ALL_SUM_SIZE = 8 * 1024 * 1024;
 constexpr const size_t ALL_SUM_BUFFERS = 2;
 constexpr const int CONN_ATTEMPTS = 5;
 constexpr const int CONN_WAIT = 1000;
+constexpr const int INIT_TIMEOUT = 20000;
 
 using GroupImpl = mlx::core::distributed::detail::GroupImpl;
 using json = nlohmann::json;
@@ -503,6 +445,7 @@ std::vector<int> make_connections(
 
   return sockets;
 }
+
 template <typename T>
 struct SumOp {
   void operator()(const T* input, T* output, size_t N) {
@@ -550,19 +493,27 @@ class RingGroup : public GroupImpl {
     size_ = nodes.size();
     int connect_to = (rank_ + 1) % size_;
 
-    // We define the connection order by having the rank_ == size_ - 1 connect
-    // first and accept after.
-    if (rank_ < connect_to) {
-      log_info(verbose_, "Rank", rank_, "accepting");
-      sockets_left_ = std::move(accept_connections(nodes[rank_]));
-      log_info(verbose_, "Rank", rank_, "connecting to", connect_to);
-      sockets_right_ = std::move(make_connections(nodes[connect_to], verbose));
-    } else {
-      log_info(verbose_, "Rank", rank_, "connecting to", connect_to);
-      sockets_right_ = std::move(make_connections(nodes[connect_to], verbose));
-      log_info(verbose_, "Rank", rank_, "accepting");
-      sockets_left_ = std::move(accept_connections(nodes[rank_]));
+    // Initialize the ring by making all the connections
+    log_info(verbose_, "Rank", rank_, "accepting");
+    log_info(verbose_, "Rank", rank_, "connecting to", connect_to);
+    auto sl = std::async(std::launch::async, accept_connections, nodes[rank_]);
+    auto sr = std::async(
+        std::launch::async, make_connections, nodes[connect_to], verbose);
+    std::future_status status_sl, status_sr;
+    for (int i = 0; i < 10; i++) {
+      status_sl = sl.wait_for(std::chrono::milliseconds(INIT_TIMEOUT / 10));
+      status_sr = sl.wait_for(std::chrono::milliseconds(INIT_TIMEOUT / 10));
+      if (status_sl == std::future_status::ready &&
+          status_sr == std::future_status::ready) {
+        break;
+      }
     }
+    if (status_sl != std::future_status::ready ||
+        status_sr != std::future_status::ready) {
+      throw std::runtime_error("[ring] Ring initialization timed out");
+    }
+    sockets_left_ = std::move(sl.get());
+    sockets_right_ = std::move(sr.get());
 
     // Failure if we couldn't make right or left sockets
     if (sockets_right_.empty()) {
@@ -628,18 +579,24 @@ class RingGroup : public GroupImpl {
   }
 
   void all_sum(const array& input, array& output, Stream stream) override {
-    SWITCH_TYPE(
-        output, all_reduce<T, SumOp<T>>(input, output, stream, SumOp<T>()));
+    dispatch_all_types(output.dtype(), [&](auto type_tag) {
+      using T = MLX_GET_TYPE(type_tag);
+      all_reduce<T, SumOp<T>>(input, output, stream, SumOp<T>());
+    });
   }
 
   void all_max(const array& input, array& output, Stream stream) override {
-    SWITCH_TYPE(
-        output, all_reduce<T, MaxOp<T>>(input, output, stream, MaxOp<T>()));
+    dispatch_all_types(output.dtype(), [&](auto type_tag) {
+      using T = MLX_GET_TYPE(type_tag);
+      all_reduce<T, MaxOp<T>>(input, output, stream, MaxOp<T>());
+    });
   }
 
   void all_min(const array& input, array& output, Stream stream) override {
-    SWITCH_TYPE(
-        output, all_reduce<T, MinOp<T>>(input, output, stream, MinOp<T>()));
+    dispatch_all_types(output.dtype(), [&](auto type_tag) {
+      using T = MLX_GET_TYPE(type_tag);
+      all_reduce<T, MinOp<T>>(input, output, stream, MinOp<T>());
+    });
   }
 
   std::shared_ptr<GroupImpl> split(int color, int key = -1) override {

mlx/export.cpp

@@ -354,7 +354,9 @@ struct PrimitiveFactory {
 
   void save(Writer& os, const std::shared_ptr<Primitive>& p) {
     serialize(os, p->stream());
-    std::string name = p->name();
+    std::ostringstream pout;
+    p->print(pout);
+    auto name = pout.str();
     name = name.substr(0, name.find(' '));
     if (auto it = name_remap.find(name); it != name_remap.end()) {
       name = it->second;

mlx/fast_primitives.h

@@ -58,7 +58,7 @@ class RMSNorm : public Custom {
       const std::vector<int>& argnums,
       const std::vector<array>& outputs) override;
 
-  DEFINE_NAME(RMSNorm)
+  DEFINE_PRINT(RMSNorm)
   bool is_equivalent(const Primitive& other) const override;
   DEFINE_INPUT_OUTPUT_SHAPE()
 
@@ -85,7 +85,7 @@ class RMSNormVJP : public Custom {
   void eval_gpu(const std::vector<array>& inputs, std::vector<array>& outputs)
       override;
 
-  DEFINE_NAME(RMSNormVJP)
+  DEFINE_PRINT(RMSNormVJP)
   bool is_equivalent(const Primitive& other) const override;
   auto state() const {
     return std::make_pair(nullptr, eps_);
@@ -118,7 +118,7 @@ class LayerNorm : public Custom {
       const std::vector<int>& argnums,
       const std::vector<array>& outputs) override;
 
-  DEFINE_NAME(LayerNorm)
+  DEFINE_PRINT(LayerNorm)
   bool is_equivalent(const Primitive& other) const override;
   DEFINE_INPUT_OUTPUT_SHAPE()
   auto state() const {
@@ -144,7 +144,7 @@ class LayerNormVJP : public Custom {
   void eval_gpu(const std::vector<array>& inputs, std::vector<array>& outputs)
       override;
 
-  DEFINE_NAME(LayerNormVJP)
+  DEFINE_PRINT(LayerNormVJP)
   bool is_equivalent(const Primitive& other) const override;
   auto state() const {
     return std::make_pair(nullptr, eps_);
@@ -186,7 +186,7 @@ class RoPE : public Custom {
       const std::vector<int>& argnums,
       const std::vector<array>& outputs) override;
 
-  DEFINE_NAME(RoPE)
+  DEFINE_PRINT(RoPE)
   bool is_equivalent(const Primitive& other) const override;
   DEFINE_INPUT_OUTPUT_SHAPE()
   auto state() const {
@@ -233,7 +233,7 @@ class ScaledDotProductAttention : public Custom {
   void eval_gpu(const std::vector<array>& inputs, array& out);
   bool is_equivalent(const Primitive& other) const override;
 
-  DEFINE_NAME(ScaledDotProductAttention);
+  DEFINE_PRINT(ScaledDotProductAttention);
   DEFINE_INPUT_OUTPUT_SHAPE()
   auto state() const {
     return std::make_tuple(nullptr, scale_, do_causal_);
@@ -263,7 +263,7 @@ class AffineQuantize : public Custom {
   void eval_gpu(const std::vector<array>& inputs, std::vector<array>& outputs)
       override;
 
-  DEFINE_NAME(AffineQuantize);
+  DEFINE_PRINT(AffineQuantize);
 
   bool is_equivalent(const Primitive& other) const override;
   std::vector<Shape> output_shapes(const std::vector<array>& inputs) override;
@@ -311,7 +311,7 @@ class CustomKernel : public Primitive {
   void eval_gpu(const std::vector<array>& inputs, std::vector<array>& outputs)
       override;
 
-  DEFINE_NAME(CustomKernel);
+  DEFINE_PRINT(CustomKernel);
 
  private:
   std::string source_;

mlx/graph_utils.cpp

@@ -93,7 +93,7 @@ void print_graph(
   os << "\n";
 
   for (auto& arr : tape) {
-    os << arr.primitive().name();
+    arr.primitive().print(os);
     os << " ";
     print_arrs(arr.inputs());
     os << " -> ";
@@ -143,7 +143,7 @@ void export_to_dot(
           os << "{ ";
           os << x.primitive_id();
           os << " [label =\"";
-          os << x.primitive().name();
+          x.primitive().print(os);
           os << "\", shape=rectangle]";
           os << "; }" << std::endl;
           // Arrows to primitive's inputs

mlx/linalg.cpp

@@ -500,7 +500,7 @@ array cross(
 void validate_eig(
     const array& a,
     const StreamOrDevice& stream,
-    const std::string& fname) {
+    const std::string fname) {
   check_cpu_stream(stream, fname);
   check_float_or_complex(a.dtype(), fname);
 

mlx/primitives.cpp

@@ -181,7 +181,7 @@ std::vector<array> Primitive::jvp(
     const std::vector<int>&) {
   std::ostringstream msg;
   msg << "[Primitive::jvp] Not implemented for ";
-  msg << name();
+  print(msg);
   msg << ".";
   throw std::invalid_argument(msg.str());
 }
@@ -193,7 +193,7 @@ std::vector<array> Primitive::vjp(
     const std::vector<array>&) {
   std::ostringstream msg;
   msg << "[Primitive::vjp] Not implemented for ";
-  msg << name();
+  print(msg);
   msg << ".";
   throw std::invalid_argument(msg.str());
 }
@@ -203,7 +203,7 @@ std::pair<std::vector<array>, std::vector<int>> Primitive::vmap(
     const std::vector<int>&) {
   std::ostringstream msg;
   msg << "[Primitive::vmap] Not implemented for ";
-  msg << name();
+  print(msg);
   msg << ".";
   throw std::invalid_argument(msg.str());
 }
@@ -211,7 +211,7 @@ std::pair<std::vector<array>, std::vector<int>> Primitive::vmap(
 std::vector<Shape> Primitive::output_shapes(const std::vector<array>&) {
   std::ostringstream msg;
   msg << "[Primitive::output_shapes] ";
-  msg << name();
+  this->print(msg);
   msg << " cannot infer output shapes.";
   throw std::invalid_argument(msg.str());
 }
@@ -743,6 +743,26 @@ bool BitwiseBinary::is_equivalent(const Primitive& other) const {
   return op_ == a_other.op_;
 }
 
+void BitwiseBinary::print(std::ostream& os) {
+  switch (op_) {
+    case BitwiseBinary::And:
+      os << "BitwiseAnd";
+      break;
+    case BitwiseBinary::Or:
+      os << "BitwiseOr";
+      break;
+    case BitwiseBinary::Xor:
+      os << "BitwiseXor";
+      break;
+    case BitwiseBinary::LeftShift:
+      os << "LeftShift";
+      break;
+    case BitwiseBinary::RightShift:
+      os << "RightShift";
+      break;
+  }
+}
+
 std::pair<std::vector<array>, std::vector<int>> BitwiseBinary::vmap(
     const std::vector<array>& inputs,
     const std::vector<int>& axes) {
@@ -2439,7 +2459,7 @@ std::vector<array> Imag::vjp(
   assert(primals.size() == 1);
   assert(argnums.size() == 1);
   return {multiply(
-      array(complex64_t{0.0f, 1.0f}, primals[0].dtype()),
+      array(complex64_t{0.0f, -1.0f}, primals[0].dtype()),
       cotangents[0],
       stream())};
 }
@@ -2768,19 +2788,15 @@ std::vector<array> Matmul::vjp(
   std::vector<int> reorder(cotan.ndim());
   std::iota(reorder.begin(), reorder.end(), 0);
   std::iter_swap(reorder.end() - 1, reorder.end() - 2);
-  auto& s = stream();
-
-  auto complex_transpose = [&](const array& x) {
-    return transpose(conjugate(x, s), reorder, s);
-  };
-
   for (auto arg : argnums) {
     if (arg == 0) {
       // M X N * (K X N).T -> M X K
-      vjps.push_back(matmul(cotan, complex_transpose(primals[1]), s));
+      vjps.push_back(
+          matmul(cotan, transpose(primals[1], reorder, stream()), stream()));
     } else {
       // (M X K).T * M X N -> K X N
-      vjps.push_back(matmul(complex_transpose(primals[0]), cotan, s));
+      vjps.push_back(
+          matmul(transpose(primals[0], reorder, stream()), cotan, stream()));
     }
   }
   return vjps;
@@ -5355,13 +5371,8 @@ std::pair<std::vector<array>, std::vector<int>> View::vmap(
   return {{view(inputs[0], dtype_, stream())}, axes};
 }
 
-const char* View::name() const {
-  if (name_.empty()) {
-    std::ostringstream os;
-    os << "View " << dtype_;
-    name_ = os.str();
-  }
-  return name_.c_str();
+void View::print(std::ostream& os) {
+  os << "View " << dtype_;
 }
 
 bool View::is_equivalent(const Primitive& other) const {

mlx/random.cpp

@@ -92,6 +92,29 @@ T below_one() {
   return f;
 }
 
+// Get the next representable value above -1.0 for half precision
+// floating point types (fp16, bf16)
+template <typename T>
+T above_minus_one() {
+  T f = T(-1.0);
+  uint16_t* m = (uint16_t*)&f;
+  *m -= 1;
+  return f;
+}
+
+// Get the next representable value above -1.0 for half precision
+// use std::nextafter as default case.
+array above_minus_one_with_default(Dtype dtype) {
+  switch (dtype) {
+    case float16:
+      return array(above_minus_one<float16_t>(), dtype);
+    case bfloat16:
+      return array(above_minus_one<bfloat16_t>(), dtype);
+    default:
+      return array(std::nextafter(-1.0f, 0.0f), dtype);
+  }
+}
+
 array uniform(
     const array& low,
     const array& high,
@@ -116,27 +139,31 @@ array uniform(
         << " from broadcasted shape " << out_shape << ".";
     throw std::invalid_argument(msg.str());
   }
-
-  // Get random values between [0, nextafter(1.0, 0.0)] since samples must
+  // Get random values between [0, nextafter(maxval, 0.0f)] since samples must
   // be in [low, high)
-  auto get_upper = [&dtype]() {
+  auto get_limits = [&dtype]() {
     switch (dtype) {
       case float32:
-        return array(std::nextafter(1.0f, 0.0f), float32);
+        return std::make_pair(
+            array(std::nextafter(1.0f, 0.0f), float32),
+            array(std::numeric_limits<uint32_t>::max(), float32));
       case float16:
-        return array(below_one<float16_t>(), float32);
+        return std::make_pair(
+            array(below_one<float16_t>(), float16),
+            array(std::numeric_limits<uint16_t>::max(), float32));
       case bfloat16:
-        return array(below_one<bfloat16_t>(), float32);
+        return std::make_pair(
+            array(below_one<bfloat16_t>(), bfloat16),
+            array(std::numeric_limits<uint16_t>::max(), float32));
       default:
         throw std::runtime_error("[uniform] Unsupported type.");
     }
   };
 
-  auto upper = get_upper();
-  auto maxval = array(std::numeric_limits<uint32_t>::max(), float32);
-  auto out = bits(shape, size_of(float32), key, stream);
-  out = divide(out, maxval, stream);
-  out = astype(minimum(out, upper, stream), dtype, stream);
+  auto [upper, maxval] = get_limits();
+  auto out = bits(shape, size_of(dtype), key, stream);
+  out = astype(divide(out, maxval, stream), dtype, stream);
+  out = minimum(out, upper, stream);
   return add(multiply(range, out, stream), lo, stream);
 }
 
@@ -156,7 +183,7 @@ inline array complex_normal(
     const std::optional<array>& key,
     StreamOrDevice s) {
   auto stream = to_stream(s);
-  auto low = array(std::nextafter(-1.0f, 0.0f), float32);
+  auto low = above_minus_one_with_default(float32);
   auto high = array(1.0f, float32);
   shape.push_back(2);
   auto samples =
@@ -180,23 +207,18 @@ array normal(
     StreamOrDevice s /* = {} */) {
   if (dtype == complex64) {
     return complex_normal(shape, loc, scale, key, s);
-  } else if (!issubdtype(dtype, floating)) {
-    throw std::invalid_argument(
-        "[normal] Can only generate uniform numbers with "
-        "floating point type.");
   }
 
   auto stream = to_stream(s);
-  auto low = array(std::nextafter(-1.0f, 0.0f), float32);
-  auto high = array(1.0f, float32);
-  auto samples = uniform(low, high, shape, float32, key, stream);
+  auto low = above_minus_one_with_default(dtype);
+  auto high = array(1.0f, dtype);
+  auto samples = uniform(low, high, shape, dtype, key, stream);
   auto applied_scale = array(std::sqrt(2.0), dtype);
   if (scale.has_value()) {
     applied_scale =
         multiply(applied_scale, astype(*scale, dtype, stream), stream);
   }
-  samples = astype(erfinv(samples, stream), dtype, stream);
-  samples = multiply(applied_scale, samples, stream);
+  samples = multiply(applied_scale, erfinv(samples, stream), stream);
   if (loc.has_value()) {
     samples = add(astype(*loc, dtype, stream), samples, stream);
   }
@@ -447,23 +469,16 @@ array laplace(
     const float scale /* = 1.0 */,
     const std::optional<array>& key /*= nullopt */,
     StreamOrDevice s /* = {} */) {
-  if (!issubdtype(dtype, floating)) {
-    throw std::invalid_argument(
-        "[laplace] Can only generate uniform numbers with real"
-        "floating point type.");
-  }
-
   auto stream = to_stream(s);
-  auto low = array(std::nextafter(-1.0f, 0.0f), float32);
-  auto high = array(1.0f, float32);
-  auto samples = uniform(low, high, shape, float32, key, stream);
+  auto low = above_minus_one_with_default(dtype);
+  auto high = array(1.0f, dtype);
+  auto samples = uniform(low, high, shape, dtype, key, stream);
   // Use inverse CDF to generate Laplacian noise
   samples = multiply(
       sign(samples, stream),
       log1p(
           multiply(array(-1.0f, dtype), abs(samples, stream), stream), stream),
       stream);
-  samples = astype(samples, dtype, stream);
 
   if (scale != 1.0) {
     samples = multiply(array(scale, dtype), samples, stream);

mlx/transforms.cpp

@@ -33,7 +33,7 @@ class Synchronizer : public Primitive {
   void eval_cpu(const std::vector<array>&, std::vector<array>&) override {}
   void eval_gpu(const std::vector<array>&, std::vector<array>&) override {}
 
-  DEFINE_NAME(Synchronize);
+  DEFINE_PRINT(Synchronize);
 };
 
 // Initialize the static tracing members from transforms_impl.h

mlx/version.h

@@ -4,7 +4,7 @@
 
 #define MLX_VERSION_MAJOR 0
 #define MLX_VERSION_MINOR 26
-#define MLX_VERSION_PATCH 5
+#define MLX_VERSION_PATCH 3
 #define MLX_VERSION_NUMERIC \
   (100000 * MLX_VERSION_MAJOR + 1000 * MLX_VERSION_MINOR + MLX_VERSION_PATCH)
 

pyproject.toml

@@ -1,6 +1,6 @@
 [build-system]
 requires = [
-  "setuptools>=80",
+  "setuptools>=42",
   "nanobind==2.4.0",
   "cmake>=3.25",
 ]

python/mlx/optimizers/optimizers.py

@@ -848,106 +848,6 @@ class Adafactor(Optimizer):
         return parameter - update
 
 
-class Muon(Optimizer):
-    r"""The Muon optimizer.
-
-    Our Muon (MomentUm Orthogonalized by Newton-schulz) optimizer follows the
-    original implementation: `Muon: An optimizer for hidden layers in neural
-    networks <https://kellerjordan.github.io/posts/muon/>`_
-
-    Note:
-        - Muon may be sub-optimal for the embedding layer, the final fully
-          connected layer, or any 0D/1D parameters. Those should be optimized
-          by a different method (e.g., :class:`AdamW`).
-        - For 4D convolutional filters, it works by flattening their last
-          dimensions.
-
-    Args:
-        learning_rate (float or callable): The learning rate.
-        momentum (float, optional): The momentum strength. Default: ``0.95``
-        weight_decay (float, optional): The weight decay (L2 penalty).
-            Default: ``0.01``
-        nesterov (bool, optional): Enables Nesterov momentum. Recommended for
-            better performance.  Default: ``True``
-        ns_steps (int, optional): Number of Newton-Schulz iteration steps for
-            orthogonalization.  Default: ``5``
-    """
-
-    def __init__(
-        self,
-        learning_rate: Union[float, Callable[[mx.array], mx.array]],
-        momentum: float = 0.95,
-        weight_decay: float = 0.01,
-        nesterov: bool = True,
-        ns_steps: int = 5,
-    ):
-        super().__init__()
-
-        self._maybe_schedule("learning_rate", learning_rate)
-        self.momentum = momentum
-        self.weight_decay = weight_decay
-        self.nesterov = nesterov
-        self.ns_steps = ns_steps
-
-    def init_single(self, parameter: mx.array, state: dict):
-        """Initialize optimizer state"""
-        state["v"] = mx.zeros_like(parameter)
-
-    def _zeropower_via_newtonschulz5(self, X, steps: int):
-        assert (
-            X.ndim == 2
-        ), f"Expected a 2D array for Newton-Schulz iteration, got shape {X.shape} instead."
-        a, b, c = (3.4445, -4.7750, 2.0315)
-        transpose_needed = X.shape[-2] > X.shape[-1]
-
-        if transpose_needed:
-            X = X.T
-
-        X = X / (mx.linalg.norm(X, keepdims=True) + 1e-7)
-
-        for _ in range(steps):
-            A = X @ X.T
-            B = mx.addmm(b * A, A, A, beta=1.0, alpha=c)
-            X = mx.addmm(a * X, B, X, beta=1.0, alpha=1.0)
-
-        if transpose_needed:
-            X = X.T
-        return X
-
-    def apply_single(self, gradient: mx.array, parameter: mx.array, state: dict):
-        """Performs the Muon parameter update"""
-
-        if self.weight_decay != 0:
-            gradient = gradient + self.weight_decay * parameter
-
-        v = self.momentum * state["v"]
-        v = v + (1 - self.momentum) * gradient
-        state["v"] = v
-
-        if self.nesterov:
-            update = gradient * (1 - self.momentum) + v * self.momentum
-        else:
-            update = v
-
-        lr = self.learning_rate.astype(gradient.dtype)
-
-        if update.ndim >= 2:
-            original_shape = update.shape
-            reshape_needed = update.ndim > 2
-
-            if reshape_needed:
-                update = mx.reshape(update, (update.shape[0], -1))
-
-            update = self._zeropower_via_newtonschulz5(update, steps=self.ns_steps)
-
-            if reshape_needed:
-                update = mx.reshape(update, original_shape)
-
-            lr *= max(1, update.shape[-2] / update.shape[-1]) ** 0.5
-
-        return parameter - lr * update
-
-
 def clip_grad_norm(grads, max_norm):
     """Clips the global norm of the gradients.