Compare commits

..

13 Commits

Author SHA1 Message Date
Angelos Katharopoulos 2ecf184f91 Add a test 2026-05-07 23:54:35 -07:00
Angelos Katharopoulos 1ea24e11f0 Fix qvm split k with batch dim 2026-05-07 23:54:35 -07:00
Cheng c9aa560577 Make device_count() return 0 when there is no GPU (#3486) 2026-05-08 08:33:55 +09:00
Angelos Katharopoulos ff57d875ea Fix indexing bug in slice update with op (#3483) 2026-05-06 17:04:50 -07:00
Cheng 80bcd1c658 [CUDA] Fix half type matmul in cutlass kernels (#3469) 2026-05-06 08:35:53 +09:00
serenposh 1fdd4e23c2 Clearer error when shape dimension overflows int32 (#3425)
Co-authored-by: Kanishk <kanishk.chores@gmail.com>
Co-authored-by: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
2026-05-05 09:53:36 +09:00
Pedro Cuenca b43965925f Define ST_F8_E8M0 (#3448) 2026-05-05 09:22:23 +09:00
Abhilash Shankarampeta 0938db7e54 Add determinant and sign-log-determinant functions to mlx.core.linalg (#3416)
Co-authored-by: Lucas Fernandes Martins <Lucas-Fernandes-Martins@users.noreply.github.com>
2026-05-05 09:06:23 +09:00
Irakli Salia e8ebdebeeb Add barrier to JACCL (#3459) 2026-04-28 09:39:56 -07:00
Cheng d7d0992d75 Reuse nightly build's ccache for release (#3458) 2026-04-28 10:54:41 +09:00
Kimon N. bdb6ff8881 Keep gguflib input-validation asserts active in release builds (#3436) 2026-04-27 08:46:57 +09:00
Long Yixing 894c948773 [CUDA] Fix qmm_naive K-tail dispatch for FP quantized kernels (#3445) 2026-04-27 08:40:14 +09:00
Angelos Katharopoulos 211e57be53 Bump minor (#3438) 2026-04-22 11:09:30 -07:00
45 changed files with 909 additions and 514 deletions
@@ -20,7 +20,7 @@ runs:
run: |
pip install auditwheel "build<=1.4.2" patchelf setuptools
python setup.py clean --all
MLX_DISABLE_SM90A_KERNELS=1 MLX_BUILD_STAGE=2 python -m build -w
MLX_BUILD_STAGE=2 python -m build -w
auditwheel repair dist/mlx_cuda*.whl \
--plat manylinux_2_35_${{ inputs.arch }} \
+4 -1
View File
@@ -14,6 +14,9 @@ inputs:
description: 'Whether to enable ccache'
required: false
default: 'true'
ccache-key:
required: false
default: 'ccache'
runs:
using: "composite"
@@ -33,7 +36,7 @@ runs:
if: ${{ inputs.use-ccache == 'true' }}
uses: hendrikmuhs/ccache-action@v1.2
with:
key: ccache-${{ runner.os }}-${{ runner.arch }}-${{ inputs.toolkit }}
key: ${{ inputs.ccache-key }}-${{ runner.os }}-${{ runner.arch }}-${{ inputs.toolkit }}
max-size: 1GB
# ccache-action bug: running "apt-get update" fails on large arm runner.
update-package-index: false
+9 -5
View File
@@ -85,20 +85,24 @@ jobs:
build_cuda_release:
if: github.repository == 'ml-explore/mlx'
runs-on: ubuntu-22-large
strategy:
matrix:
arch: ['x86_64', 'aarch64']
toolkit: ['cuda-12.9', 'cuda-13.0']
runs-on: ${{ matrix.arch == 'x86_64' && 'ubuntu-22-large' || 'ubuntu-22-large-arm' }}
steps:
- uses: actions/checkout@v6
- uses: ./.github/actions/setup-linux
with:
toolkit: 'cuda-12.9'
toolkit: ${{ matrix.toolkit }}
ccache-key: 'ccache-release'
- name: Build Python package
uses: ./.github/actions/build-cuda-release
with:
toolkit: 'cuda-12.9'
arch: 'x86_64'
arch: ${{ matrix.arch }}
- name: Upload artifacts
uses: actions/upload-artifact@v7
with:
name: mlx-cuda
name: mlx-${{ matrix.toolkit }}-${{ matrix.arch }}
path: wheelhouse/mlx_cuda_*.whl
retention-days: 7
+1 -1
View File
@@ -146,7 +146,7 @@ jobs:
- uses: ./.github/actions/setup-linux
with:
toolkit: ${{ matrix.toolkit }}
use-ccache: false
ccache-key: 'ccache-release'
- name: Build Python package
uses: ./.github/actions/build-cuda-release
with:
+2
View File
@@ -14,6 +14,7 @@ Linear Algebra
cholesky
cholesky_inv
cross
det
qr
svd
eigvals
@@ -23,5 +24,6 @@ Linear Algebra
lu
lu_factor
pinv
slogdet
solve
solve_triangular
+1 -10
View File
@@ -98,15 +98,6 @@ void Recv::eval_cpu(
void ReduceScatter::eval_cpu(
const std::vector<array>& inputs,
std::vector<array>& outputs) {
assert(inputs.size() == 1);
assert(outputs.size() == 1);
auto [in, copied] = ensure_row_contiguous(inputs[0], stream());
outputs[0].set_data(allocator::malloc(outputs[0].nbytes()));
distributed::detail::sum_scatter(group(), in, outputs[0], stream());
if (copied) {
auto& enc = cpu::get_command_encoder(stream());
enc.add_temporary(in);
}
throw std::runtime_error("[ReduceScatter] Not implemented yet.");
}
} // namespace mlx::core::distributed
+2 -3
View File
@@ -67,11 +67,10 @@ void luf_impl(
/* ipiv */ reinterpret_cast<int*>(pivots_ptr),
/* info */ &info);
if (info != 0) {
if (info < 0) {
std::stringstream ss;
ss << "[LUF::eval_cpu] sgetrf_ failed with code " << info
<< ((info > 0) ? " because matrix is singular"
: " because argument had an illegal value");
<< " because argument had an illegal value";
throw std::runtime_error(ss.str());
}
+2 -3
View File
@@ -168,9 +168,8 @@ set_target_properties(mlx PROPERTIES CUDA_ARCHITECTURES
"${MLX_CUDA_ARCHITECTURES}")
# Skip Hopper-only kernels when not building for sm90a.
if(NOT DEFINED ENV{MLX_DISABLE_SM90A_KERNELS}
AND (("90a" IN_LIST MLX_CUDA_ARCHITECTURES) OR ("90a-real" IN_LIST
MLX_CUDA_ARCHITECTURES)))
if(("90a" IN_LIST MLX_CUDA_ARCHITECTURES) OR ("90a-real" IN_LIST
MLX_CUDA_ARCHITECTURES))
target_compile_definitions(mlx PRIVATE MLX_CUDA_SM90A_ENABLED)
endif()
+10 -5
View File
@@ -43,6 +43,12 @@ class GatherGemm {
using ElementD = typename CollectiveEpilogue::ElementD;
using StrideD = typename CollectiveEpilogue::StrideD;
static_assert(
cute::is_same_v<
ElementAccumulator,
typename CollectiveEpilogue::ElementAccumulator>,
"Mainloop and epilogue do not agree on accumulator value type.");
static constexpr int SharedStorageSize = static_cast<int>(cute::max(
sizeof(typename CollectiveMainloop::SharedStorage),
sizeof(typename CollectiveEpilogue::SharedStorage)));
@@ -98,7 +104,9 @@ class GatherGemm {
CUTLASS_DEVICE void operator()(const Params& params, char* smem_buf) {
int thread_idx = int(threadIdx.x);
auto [m_coord, n_coord, l_coord] = uint3(blockIdx);
int m_coord = int(blockIdx.x);
int n_coord = int(blockIdx.y);
int l_coord = int(blockIdx.z);
auto shape_MNKL = append<4>(params.problem_shape, Int<1>{});
auto cta_tile = TileShape{};
@@ -220,7 +228,7 @@ void gather_mm(
using TileShape = Shape<_128, _128, _8>;
using DispatchPolicy = cutlass::gemm::MainloopSm70TwoStage;
using TiledMma = TiledMMA<
MMA_Atom<UniversalFMA<Accumulator, Element, Element, Element>>,
MMA_Atom<UniversalFMA<Accumulator, Element, Element, Accumulator>>,
Layout<Shape<_16, _16, _1>>>;
using CopyTraitsA = SimtCopyTraits<Element, k_major_a.value>;
@@ -296,9 +304,6 @@ void cutlass_gather_mm(
int n = out.shape(-1);
int k = a.shape(-1);
int l = out.size() / (m * n);
if (m < 16 || n < 16) {
throw std::invalid_argument("[gather_mm] M/N is too small.");
}
encoder.set_input_array(a);
encoder.set_input_array(b);
@@ -245,7 +245,7 @@ void grouped_gemm_v2(
LayoutB,
cutlass::ComplexTransform::kNone,
GemmConfiguration::kAlignmentAB,
typename GemmConfiguration::Element,
typename GemmConfiguration::Accumulator,
cutlass::layout::RowMajor,
typename GemmConfiguration::Accumulator,
typename GemmConfiguration::OpClass,
+7 -6
View File
@@ -52,7 +52,7 @@ bool supports_qmm_sm90(
if (!biases) {
return false;
}
if (!x.flags().row_contiguous || !is_last_2_dims_row_contiguous(w) ||
if (!is_last_2_dims_row_contiguous(w) ||
!is_last_2_dims_row_contiguous(scales) ||
!is_last_2_dims_row_contiguous(*biases)) {
return false;
@@ -139,7 +139,7 @@ bool supports_qmm_sm80(
if ((n % 128 != 0) || (k % std::max(64, group_size) != 0)) {
return false;
}
if (!x.flags().row_contiguous || !is_last_2_dims_row_contiguous(w) ||
if (!is_last_2_dims_row_contiguous(w) ||
!is_last_2_dims_row_contiguous(scales)) {
return false;
}
@@ -224,7 +224,7 @@ bool supports_qmm_naive(
if (transpose && (k % std::max(64, group_size) != 0)) {
return false;
}
if (!x.flags().row_contiguous || !is_last_2_dims_row_contiguous(w) ||
if (!is_last_2_dims_row_contiguous(w) ||
!is_last_2_dims_row_contiguous(scales)) {
return false;
}
@@ -265,7 +265,7 @@ void qmm_naive(
if constexpr (k_major.value) {
if (has_k_residue) {
throw std::invalid_argument(
"[quantized_matmul] K must be multiples of group_size.");
"[quantized_matmul] K must be multiples of max(64, group_size).");
}
f.template operator()<false>();
} else {
@@ -276,7 +276,8 @@ void qmm_naive(
};
int m = out.ndim() > 1 ? out.shape(-2) : 1;
int k = x.shape(-1);
bool has_k_residue = k % group_size != 0;
int tile_k = std::max(64, group_size);
bool has_k_residue = k % tile_k != 0;
bool sm80 = encoder.device().compute_capability_major() >= 8;
dispatch_bool(transpose, [&](auto k_major) {
dispatch_k(k_major, has_k_residue, [&]<bool HasKResidue>() {
@@ -342,7 +343,7 @@ bool supports_qmv(
if (k % 8 != 0) {
return false;
}
if (!x.flags().row_contiguous || !is_last_2_dims_row_contiguous(w) ||
if (!is_last_2_dims_row_contiguous(w) ||
!is_last_2_dims_row_contiguous(scales)) {
return false;
}
+3 -1
View File
@@ -60,7 +60,9 @@ __global__ void qmm_naive_kernel(
CUTE_STATIC_ASSERT_V(congruent(select<0,1,3>(shape_MNKL), dC));
int thread_idx = int(threadIdx.x);
auto [m_coord, n_coord, l_coord] = static_cast<uint3>(blockIdx);
int m_coord = int(blockIdx.x);
int n_coord = int(blockIdx.y);
int l_coord = int(blockIdx.z);
auto m_max_coord = size<0>(shape_MNKL) - size<0>(cta_tiler) * m_coord; // M - BLK_M * m_coord
auto n_max_coord = size<1>(shape_MNKL) - size<1>(cta_tiler) * n_coord; // N - BLK_N * n_coord
+3 -1
View File
@@ -48,7 +48,9 @@ __global__ void qmm_sm80_kernel(
CUTE_STATIC_ASSERT_V(congruent(select<0,1,3>(shape_MNKL), dC));
int thread_idx = int(threadIdx.x);
auto [m_coord, n_coord, l_coord] = static_cast<uint3>(blockIdx);
int m_coord = int(blockIdx.x);
int n_coord = int(blockIdx.y);
int l_coord = int(blockIdx.z);
// For gather, use index lookup for input batch slicing.
uint32_t a_batch = lhs_indices ? lhs_indices[l_coord] : l_coord;
+6 -4
View File
@@ -17,7 +17,7 @@ void QuantizedMatmul::eval_gpu(const std::vector<array>& inputs, array& out) {
auto& s = stream();
auto& encoder = cu::get_command_encoder(s);
const array& x = inputs[0];
array x = ensure_row_contiguous(inputs[0], encoder, s);
const array& w = inputs[1];
const array& scales = inputs[2];
std::optional<array> biases;
@@ -146,15 +146,17 @@ void GatherQMM::eval_gpu(const std::vector<array>& inputs, array& out) {
auto& s = stream();
auto& encoder = cu::get_command_encoder(s);
const array& x = inputs[0];
array x = ensure_row_contiguous(inputs[0], encoder, s);
const array& w = inputs[1];
const array& scales = inputs[2];
std::optional<array> biases;
if (inputs.size() == 6) {
biases = inputs[3];
}
array lhs_indices = ensure_contiguous(inputs[inputs.size() - 2], encoder, s);
array rhs_indices = ensure_contiguous(inputs[inputs.size() - 1], encoder, s);
array lhs_indices =
ensure_row_contiguous(inputs[inputs.size() - 2], encoder, s);
array rhs_indices =
ensure_row_contiguous(inputs[inputs.size() - 1], encoder, s);
int M = out.ndim() > 1 ? out.shape(-2) : 1;
int N = out.shape(-1);
+6 -1
View File
@@ -13,7 +13,12 @@ bool is_available() {
}
int device_count() {
return 1;
try {
metal::device(Device::gpu);
return 1;
} catch (...) {
return 0;
}
}
const std::unordered_map<std::string, std::variant<std::string, size_t>>&
+35 -4
View File
@@ -533,6 +533,7 @@ METAL_FUNC void fp_qvm_impl(
device T* y,
const int in_vec_size,
const int out_vec_size,
const int in_vec_stride,
uint3 tid [[threadgroup_position_in_grid]],
uint simd_gid [[simdgroup_index_in_threadgroup]],
uint simd_lid [[thread_index_in_simdgroup]]) {
@@ -563,7 +564,7 @@ METAL_FUNC void fp_qvm_impl(
int out_col = pack_factor * tn * (tid.y * num_simdgroups + simd_gid);
ws += out_col * bytes_per_pack / pack_factor + simd_lid * out_vec_size_w;
scales += out_col / group_size + simd_lid * out_vec_size_g;
x += tid.x * in_vec_size + simd_lid;
x += tid.x * in_vec_stride + simd_lid;
y += tid.x * out_vec_size + out_col;
if (out_col >= out_vec_size) {
@@ -1122,7 +1123,16 @@ template <typename T, const int group_size, int bits, bool batched>
tid);
}
fp_qvm_impl<T, group_size, bits>(
w, scales, x, y, in_vec_size, out_vec_size, tid, simd_gid, simd_lid);
w,
scales,
x,
y,
in_vec_size,
out_vec_size,
in_vec_size,
tid,
simd_gid,
simd_lid);
}
template <typename T, const int group_size, int bits, int split_k = 32>
@@ -1164,8 +1174,20 @@ template <typename T, const int group_size, int bits, int split_k = 32>
int in_vec_size_adj =
tid.z % split_k == split_k - 1 ? final_block_size : in_vec_size;
// The in_vec_stride is the full K dimension, not the partition size
int in_vec_stride = (split_k - 1) * in_vec_size + final_block_size;
fp_qvm_impl<T, group_size, bits>(
w, scales, x, y, in_vec_size_adj, out_vec_size, tid, simd_gid, simd_lid);
w,
scales,
x,
y,
in_vec_size_adj,
out_vec_size,
in_vec_stride,
tid,
simd_gid,
simd_lid);
}
template <
@@ -1423,7 +1445,16 @@ template <typename T, int group_size, int bits>
s_strides,
tid);
fp_qvm_impl<T, group_size, bits>(
w, scales, x, y, in_vec_size, out_vec_size, tid, simd_gid, simd_lid);
w,
scales,
x,
y,
in_vec_size,
out_vec_size,
in_vec_size,
tid,
simd_gid,
simd_lid);
}
template <
+3 -1
View File
@@ -80,7 +80,9 @@ template <
uint3 gsize [[threads_per_grid]]) {
Op op;
IdxT idx = IdxT(gid.z) * gsize.y + gid.y * gsize.x + gid.x * NWORK;
IdxT idx =
(IdxT(gid.z) * IdxT(gsize.y) + IdxT(gid.y)) * IdxT(gsize.x) * NWORK +
IdxT(gid.x) * NWORK;
IdxT out_idx;
IdxT update_idx;
+8 -1
View File
@@ -983,6 +983,7 @@ METAL_FUNC void qvm_impl(
device T* y,
const int in_vec_size,
const int out_vec_size,
const int in_vec_stride,
uint3 tid [[threadgroup_position_in_grid]],
uint simd_gid [[simdgroup_index_in_threadgroup]],
uint simd_lid [[thread_index_in_simdgroup]]) {
@@ -1016,7 +1017,7 @@ METAL_FUNC void qvm_impl(
ws += out_col * bytes_per_pack / pack_factor + simd_lid * out_vec_size_w;
scales += out_col / group_size + simd_lid * out_vec_size_g;
biases += out_col / group_size + simd_lid * out_vec_size_g;
x += tid.x * in_vec_size + simd_lid;
x += tid.x * in_vec_stride + simd_lid;
y += tid.x * out_vec_size + out_col;
if (out_col >= out_vec_size) {
@@ -1643,6 +1644,7 @@ template <typename T, const int group_size, const int bits, bool batched>
y,
in_vec_size,
out_vec_size,
in_vec_size,
tid,
simd_gid,
simd_lid);
@@ -1691,6 +1693,9 @@ template <typename T, const int group_size, const int bits, int split_k = 32>
int in_vec_size_adj =
tid.z % split_k == split_k - 1 ? final_block_size : in_vec_size;
// The in_vec_stride is the full K dimension, not the partition size
int in_vec_stride = (split_k - 1) * in_vec_size + final_block_size;
qvm_impl<T, group_size, bits>(
w,
scales,
@@ -1699,6 +1704,7 @@ template <typename T, const int group_size, const int bits, int split_k = 32>
y,
in_vec_size_adj,
out_vec_size,
in_vec_stride,
tid,
simd_gid,
simd_lid);
@@ -2077,6 +2083,7 @@ template <typename T, int group_size, int bits>
y,
in_vec_size,
out_vec_size,
in_vec_size,
tid,
simd_gid,
simd_lid);
+1 -10
View File
@@ -145,16 +145,7 @@ class JACCLGroup : public GroupImpl {
}
void sum_scatter(const array& input, array& output, Stream stream) override {
auto in_ptr = input.data<char>();
auto out_ptr = output.data<char>();
size_t n_bytes = input.nbytes();
int dtype = dtype_to_jaccl_dtype(output.dtype());
auto& encoder = cpu::get_command_encoder(stream);
encoder.set_input_array(input);
encoder.set_output_array(output);
encoder.dispatch([in_ptr, out_ptr, n_bytes, dtype, this]() {
group_->sum_scatter(in_ptr, out_ptr, n_bytes, dtype);
});
throw std::runtime_error("[jaccl] sum_scatter not supported.");
}
std::shared_ptr<GroupImpl> split(int color, int key = -1) override {
+5
View File
@@ -29,6 +29,8 @@ in macOS 26.2.
- **Point-to-Point Operations**:
- `send`: Send data to a specific node
- `recv`: Receive data from a specific node
- **Synchronization**:
- `barrier`: Block until all nodes in the group reach this point
- **Type Support**: Bool, Int8-64, UInt8-64, Float16, BFloat16, Float32,
Float64, Complex64
@@ -286,6 +288,9 @@ class Group {
// Simple send/recv primitives.
virtual void send(const void* input, size_t n_bytes, int dst) = 0;
virtual void recv(void* output, size_t n_bytes, int src) = 0;
// Block until every rank reaches this point.
virtual void barrier() = 0;
};
```
@@ -35,6 +35,7 @@ endfunction()
# Examples
build_example(minimal_env.cpp)
build_example(minimal_cfg.cpp)
build_example(minimal_barrier.cpp)
# Benchmarks
build_example(allreduce_bench.cpp)
@@ -0,0 +1,42 @@
// Copyright © 2026 Apple Inc.
//
// Exercises Group::barrier(). Ranks arrive at the barrier at staggered times;
// after the barrier returns we do a small all_sum to confirm the group is
// healthy and that barrier() carried the correct fence semantics.
#include <chrono>
#include <iostream>
#include <thread>
#include <jaccl/jaccl.h>
int main() {
auto group = jaccl::init();
if (!group) {
std::cerr << "Failed to initialize JACCL" << std::endl;
return 1;
}
int rank = group->rank();
int size = group->size();
std::this_thread::sleep_for(std::chrono::milliseconds(100 * rank));
std::cout << "rank " << rank << " entering barrier" << std::endl;
group->barrier();
std::cout << "rank " << rank << " exited barrier" << std::endl;
int in = rank + 1;
int out = 0;
group->all_sum(&in, &out, sizeof(in), jaccl::Int32);
int expected = size * (size + 1) / 2;
if (out != expected) {
std::cerr << "rank " << rank << ": post-barrier all_sum mismatch (got "
<< out << ", expected " << expected << ")" << std::endl;
return 1;
}
std::cout << "rank " << rank << ": post-barrier all_sum OK (" << out << ")"
<< std::endl;
return 0;
}
+1 -3
View File
@@ -28,11 +28,9 @@ class Group {
virtual void all_gather(const void* input, void* output, size_t n_bytes) = 0;
virtual void
sum_scatter(const void* input, void* output, size_t n_bytes, int dtype) = 0;
virtual void send(const void* input, size_t n_bytes, int dst) = 0;
virtual void recv(void* output, size_t n_bytes, int src) = 0;
virtual void barrier() = 0;
};
/**
+3 -21
View File
@@ -176,17 +176,6 @@ void MeshGroup::all_gather(const void* input, void* output, size_t n_bytes) {
static_cast<const char*>(input), static_cast<char*>(output), n_bytes);
}
void MeshGroup::sum_scatter(
const void* input,
void* output,
size_t n_bytes,
int dtype) {
dispatch_all_types(dtype, [&](auto type_tag) {
using T = JACCL_GET_TYPE(type_tag);
reduce_scatter<T>(input, output, n_bytes, SumOp<T>{});
});
}
void MeshGroup::send(const void* input, size_t n_bytes, int dst) {
mesh_.send(static_cast<const char*>(input), n_bytes, dst);
}
@@ -195,16 +184,9 @@ void MeshGroup::recv(void* output, size_t n_bytes, int src) {
mesh_.recv(static_cast<char*>(output), n_bytes, src);
}
template <typename T, typename ReduceOp>
void MeshGroup::reduce_scatter(
const void* input,
void* output,
size_t n_bytes,
ReduceOp reduce_op) {
auto in_ptr = static_cast<const T*>(input);
auto out_ptr = static_cast<T*>(output);
int64_t count = n_bytes / sizeof(T);
mesh_.reduce_scatter(in_ptr, out_ptr, count, reduce_op);
void MeshGroup::barrier() {
uint8_t b = 0;
all_sum(&b, &b, sizeof(b), Dtype::UInt8);
}
template <typename T, typename ReduceOp>
+2 -10
View File
@@ -44,20 +44,12 @@ class MeshGroup : public Group {
void all_gather(const void* input, void* output, size_t n_bytes) override;
void sum_scatter(const void* input, void* output, size_t n_bytes, int dtype)
override;
void send(const void* input, size_t n_bytes, int dst) override;
void recv(void* output, size_t n_bytes, int src) override;
private:
template <typename T, typename ReduceOp>
void reduce_scatter(
const void* input,
void* output,
size_t n_bytes,
ReduceOp reduce_op);
void barrier() override;
private:
template <typename T, typename ReduceOp>
void all_reduce(
const void* input,
-159
View File
@@ -29,165 +29,6 @@ class MeshImpl {
MeshImpl() : rank_(0), size_(1) {}
template <typename T, typename ReduceOp>
void reduce_scatter(const T* in, T* out, int64_t size, ReduceOp reduce_op) {
// Reduce-scatter for mesh topology.
//
// Each rank sends its entire input to all other ranks.
// Each rank reduces only its assigned chunk from all received inputs.
auto [sz, buffer_size] = buffer_size_from_message(size * sizeof(T));
int64_t N = buffer_size / sizeof(T);
constexpr int PIPELINE = 2;
constexpr int WC_NUM = PIPELINE * MESH_MAX_PEERS * 2;
int64_t total = static_cast<int64_t>(size);
int num_peers = size_ - 1;
// Calculate chunk for this rank
int64_t chunk_size = (total + size_ - 1) / size_;
int64_t my_chunk_start = rank_ * chunk_size;
int64_t my_chunk_size = std::min(chunk_size, total - my_chunk_start);
// Initialize output with our own chunk
if (my_chunk_size > 0) {
std::copy_n(in + my_chunk_start, my_chunk_size, out);
}
// A helper for convenient access to the staging buffer.
auto local_staging = [&](int buff) -> T* {
return reinterpret_cast<T*>(staging_mem_.get() + buff * MAX_BUFFER_SIZE);
};
// Counters to maintain the state of transfers
int in_flight = 0;
int64_t read_offset = 0;
int completed_send_count[PIPELINE] = {0};
int recv_end[MESH_MAX_PEERS] = {0};
int reduce_chunk = 0;
int reduce_rank = 0;
// Total number of chunks
int64_t total_chunks = (total + N - 1) / N;
// Prefill the pipeline
int buff = 0;
while (read_offset < total && buff < PIPELINE) {
post_recv_all(sz, buff);
// Copy the local data to send buffer and staging buffer
int64_t elems = std::min(N, total - read_offset);
std::copy(
in + read_offset, in + read_offset + elems, local_staging(buff));
std::copy(
in + read_offset,
in + read_offset + elems,
send_buffer(sz, buff).begin<T>());
recv_end[rank_]++;
post_send_all(sz, buff);
buff++;
in_flight += 2 * num_peers;
read_offset += N;
}
// Main loop
while (reduce_chunk < total_chunks) {
// Poll the hardware for completions.
ibv_wc wc[WC_NUM];
int n = poll(connections_, WC_NUM, wc);
for (int i = 0; i < n; i++) {
int work_type = wc[i].wr_id >> 16;
int buff = (wc[i].wr_id >> 8) & 0xff;
int rank = wc[i].wr_id & 0xff;
in_flight--;
if (work_type == SEND_WR && read_offset < total) {
completed_send_count[buff]++;
if (completed_send_count[buff] == num_peers) {
int64_t elems = std::min(N, total - read_offset);
std::copy(
in + read_offset,
in + read_offset + elems,
local_staging(buff));
std::copy(
in + read_offset,
in + read_offset + elems,
send_buffer(sz, buff).begin<T>());
recv_end[rank_]++;
post_send_all(sz, buff);
completed_send_count[buff] = 0;
in_flight += num_peers;
read_offset += N;
}
}
else if (work_type == RECV_WR) {
recv_end[rank]++;
}
}
// Process the received chunks in order, reducing only our chunk
while (reduce_chunk < total_chunks) {
int64_t w = static_cast<int64_t>(reduce_chunk) * N;
if (w >= read_offset) {
break;
}
if (recv_end[reduce_rank] <= reduce_chunk) {
break;
}
int b = reduce_chunk % PIPELINE;
int64_t elems = std::min(N, total - w);
// Check if this chunk overlaps with our output chunk
int64_t overlap_start = std::max(w, my_chunk_start);
int64_t overlap_end =
std::min(w + elems, my_chunk_start + my_chunk_size);
if (overlap_start < overlap_end) {
int64_t out_offset = overlap_start - my_chunk_start;
int64_t in_offset = overlap_start - w;
int64_t overlap_size = overlap_end - overlap_start;
// Data is read from the staging area for our own rank
if (reduce_rank == rank_) {
reduce_op(
local_staging(b) + in_offset, out + out_offset, overlap_size);
}
// Data is read from the recv buffers for other ranks
else {
reduce_op(
recv_buffer(sz, b, reduce_rank).begin<T>() + in_offset,
out + out_offset,
overlap_size);
}
}
// Check if we need to post another receive
int64_t next_chunk = static_cast<int64_t>(reduce_chunk) + PIPELINE;
if (next_chunk < total_chunks) {
recv_from(sz, reduce_rank, b);
in_flight++;
}
// Move to next rank's data for this chunk
reduce_rank++;
if (reduce_rank >= size_) {
reduce_rank = 0;
reduce_chunk++;
}
}
}
// Drain remaining in-flight completions
while (in_flight > 0) {
ibv_wc wc[WC_NUM];
int n = poll(connections_, WC_NUM, wc);
in_flight -= n;
}
}
template <typename T, typename ReduceOp>
void all_reduce(const T* in, T* out, int64_t size, ReduceOp reduce_op) {
// Fully connected all reduce with deterministic reduction order.
+3 -32
View File
@@ -166,17 +166,6 @@ void RingGroup::all_gather(const void* input, void* output, size_t n_bytes) {
n_conns_);
}
void RingGroup::sum_scatter(
const void* input,
void* output,
size_t n_bytes,
int dtype) {
dispatch_all_types(dtype, [&](auto type_tag) {
using T = JACCL_GET_TYPE(type_tag);
reduce_scatter<T>(input, output, n_bytes, SumOp<T>{});
});
}
void RingGroup::send(const void* input, size_t n_bytes, int dst) {
int right = (rank_ + 1) % size_;
int left = (rank_ + size_ - 1) % size_;
@@ -201,27 +190,9 @@ void RingGroup::recv(void* output, size_t n_bytes, int src) {
ring_.recv(static_cast<char*>(output), n_bytes, src, n_conns_);
}
template <typename T, typename ReduceOp>
void RingGroup::reduce_scatter(
const void* input,
void* output,
size_t n_bytes,
ReduceOp reduce_op) {
auto in_ptr = static_cast<const T*>(input);
auto out_ptr = static_cast<T*>(output);
int64_t count = n_bytes / sizeof(T);
if (count < size_ * 2 * n_conns_) {
ring_.reduce_scatter<1, T, ReduceOp>(in_ptr, out_ptr, count, 1, reduce_op);
return;
}
if (n_bytes <= 65536) {
ring_.reduce_scatter<2, T, ReduceOp>(in_ptr, out_ptr, count, 1, reduce_op);
return;
}
ring_.reduce_scatter<2, T, ReduceOp>(
in_ptr, out_ptr, count, n_conns_, reduce_op);
void RingGroup::barrier() {
uint8_t b = 0;
all_sum(&b, &b, sizeof(b), Dtype::UInt8);
}
template <typename T, typename ReduceOp>
+2 -10
View File
@@ -45,20 +45,12 @@ class RingGroup : public Group {
void all_gather(const void* input, void* output, size_t n_bytes) override;
void sum_scatter(const void* input, void* output, size_t n_bytes, int dtype)
override;
void send(const void* input, size_t n_bytes, int dst) override;
void recv(void* output, size_t n_bytes, int src) override;
private:
template <typename T, typename ReduceOp>
void reduce_scatter(
const void* input,
void* output,
size_t n_bytes,
ReduceOp reduce_op);
void barrier() override;
private:
template <typename T, typename ReduceOp>
void all_reduce(
const void* input,
-168
View File
@@ -45,174 +45,6 @@ class RingImpl {
RingImpl() : rank_(0), size_(1), n_conns_(0) {}
template <int MAX_DIR, typename T, typename ReduceOp>
void reduce_scatter(
const T* in_ptr,
T* out_ptr,
int64_t size,
int n_wires,
ReduceOp reduce_op) {
constexpr int PIPELINE = 2;
constexpr int WC_NUM = PIPELINE * RING_MAX_CONNS * 2 * MAX_DIR;
int64_t chunk_size = size / size_;
int64_t size_per_wire =
(chunk_size + (MAX_DIR * n_wires) - 1) / (MAX_DIR * n_wires);
auto [sz, N] = buffer_size_from_message(size_per_wire * sizeof(T));
N /= sizeof(T);
int64_t n_steps = (size_per_wire + N - 1) / N;
// Counters to maintain the state of transfers
int in_flight = 0;
int64_t send_offset[MAX_DIR];
int64_t recv_offset[MAX_DIR];
int64_t send_limits[MAX_DIR];
int64_t recv_limits[MAX_DIR];
int send_count[MAX_DIR * RING_MAX_CONNS] = {0};
int recv_count[MAX_DIR * RING_MAX_CONNS] = {0};
send_offset[0] = ((rank_ + size_ - 1) % size_) * chunk_size;
recv_offset[0] = ((rank_ + size_ - 2) % size_) * chunk_size;
if constexpr (MAX_DIR == 2) {
send_offset[1] = ((rank_ + 1) % size_) * chunk_size;
recv_offset[1] = ((rank_ + 2) % size_) * chunk_size;
send_limits[0] = std::min(
n_wires * size_per_wire, std::max<int64_t>(0, size - send_offset[0]));
send_limits[1] =
std::min(chunk_size, std::max<int64_t>(0, size - send_offset[1]));
recv_limits[0] = std::min(
n_wires * size_per_wire, std::max<int64_t>(0, size - recv_offset[0]));
recv_limits[1] =
std::min(chunk_size, std::max<int64_t>(0, size - recv_offset[1]));
} else {
send_limits[0] =
std::min(chunk_size, std::max<int64_t>(0, size - send_offset[0]));
recv_limits[0] =
std::min(chunk_size, std::max<int64_t>(0, size - recv_offset[0]));
}
for (int k = 0; k < size_ - 1; k++) {
const T* read_ptr = (k == 0) ? in_ptr : out_ptr;
int64_t read_offset[MAX_DIR];
read_offset[0] = (k == 0) ? send_offset[0] : 0;
if constexpr (MAX_DIR == 2) {
read_offset[1] = (k == 0) ? send_offset[1] : 0;
}
// Prefill the pipeline
int buff = 0;
while (buff < n_steps && buff < PIPELINE) {
// Post receives
post_recv_all<MAX_DIR>(sz, buff, n_wires);
in_flight += MAX_DIR * n_wires;
// Copy to send buffers and also to the output buffer
for (int lr = 0; lr < MAX_DIR; lr++) {
for (int lw = 0; lw < n_wires; lw++) {
// send
int64_t offset = lw * N +
send_count[lr * RING_MAX_CONNS + lw] * n_wires * N +
lr * n_wires * size_per_wire;
std::copy(
read_ptr + read_offset[lr] + offset,
read_ptr + read_offset[lr] +
std::max(offset, std::min(offset + N, send_limits[lr])),
send_buffer(sz, buff, lr, lw).begin<T>());
send_count[lr * RING_MAX_CONNS + lw]++;
in_flight++;
send_to(sz, buff, lr, lw);
// copy for in-place recv reduce
offset -= n_wires * N;
std::copy(
in_ptr + recv_offset[lr] + offset,
in_ptr + recv_offset[lr] +
std::max(offset, std::min(offset + N, recv_limits[lr])),
out_ptr + offset);
}
}
buff++;
}
while (in_flight > 0) {
ibv_wc wc[WC_NUM];
int n = poll(left_, right_, WC_NUM, wc);
for (int i = 0; i < n; i++) {
int work_type = wc[i].wr_id >> 16;
int buff = (wc[i].wr_id >> 8) & 0xff;
int wire = wc[i].wr_id & 0xff;
int lr = wire / RING_MAX_CONNS;
int lw = wire % RING_MAX_CONNS;
in_flight--;
if (work_type == SEND_WR) {
int64_t offset = lw * N + send_count[wire] * n_wires * N +
lr * n_wires * size_per_wire;
// More stuff to send
if (send_count[wire] < n_steps) {
std::copy(
read_ptr + read_offset[lr] + offset,
read_ptr + read_offset[lr] +
std::max(offset, std::min(offset + N, send_limits[lr])),
send_buffer(sz, buff, lr, lw).begin<T>());
send_to(sz, buff, lr, lw);
in_flight++;
send_count[wire]++;
}
// Copy the input chunk into the output
offset -= n_wires * N;
std::copy(
in_ptr + recv_offset[lr] + offset,
in_ptr + recv_offset[lr] +
std::max(offset, std::min(offset + N, recv_limits[lr])),
out_ptr + offset);
}
else if (work_type == RECV_WR) {
int64_t offset = lw * N + recv_count[wire] * n_wires * N +
lr * n_wires * size_per_wire;
reduce_op(
recv_buffer(sz, buff, lr, lw).begin<T>(),
out_ptr + recv_offset[lr] + offset,
std::max<int64_t>(0, std::min(N, recv_limits[lr] - offset)));
recv_count[wire]++;
if (recv_count[wire] + (PIPELINE - 1) < n_steps) {
recv_from(sz, buff, lr, lw);
in_flight++;
}
}
}
}
send_offset[0] = (send_offset[0] + size - chunk_size) % size;
recv_offset[0] = (recv_offset[0] + size - chunk_size) % size;
if constexpr (MAX_DIR == 2) {
send_offset[1] = (send_offset[1] + chunk_size) % size;
recv_offset[1] = (recv_offset[1] + chunk_size) % size;
send_limits[0] = std::min(
n_wires * size_per_wire,
std::max<int64_t>(0, size - send_offset[0]));
send_limits[1] =
std::min(chunk_size, std::max<int64_t>(0, size - send_offset[1]));
recv_limits[0] = std::min(
n_wires * size_per_wire,
std::max<int64_t>(0, size - recv_offset[0]));
recv_limits[1] =
std::min(chunk_size, std::max<int64_t>(0, size - recv_offset[1]));
} else {
send_limits[0] =
std::min(chunk_size, std::max<int64_t>(0, size - send_offset[0]));
recv_limits[0] =
std::min(chunk_size, std::max<int64_t>(0, size - recv_offset[0]));
}
for (int i = 0; i < MAX_DIR * RING_MAX_CONNS; i++) {
send_count[i] = recv_count[i] = 0;
}
}
}
template <int MAX_DIR, typename T, typename ReduceOp>
void all_reduce(
const T* in_ptr,
+5
View File
@@ -17,6 +17,11 @@ if(MLX_BUILD_GGUF)
PRIVATE $<BUILD_INTERFACE:${gguflib_SOURCE_DIR}>)
add_library(gguflib STATIC ${gguflib_SOURCE_DIR}/fp16.c
${gguflib_SOURCE_DIR}/gguflib.c)
# gguflib uses assert() to reject malformed tensor headers (e.g. ndim > 8).
# Those checks are otherwise compiled out by -DNDEBUG in release builds, which
# leaves out-of-bounds reads/writes unguarded when loading untrusted GGUF
# files. Force NDEBUG off for this target so the asserts stay live.
target_compile_options(gguflib PRIVATE -UNDEBUG)
target_link_libraries(mlx PRIVATE $<BUILD_INTERFACE:gguflib>)
target_sources(mlx PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/gguf.cpp
${CMAKE_CURRENT_SOURCE_DIR}/gguf_quants.cpp)
+3
View File
@@ -28,6 +28,7 @@ using json = nlohmann::json;
#define ST_U32 "U32"
#define ST_U64 "U64"
#define ST_F8_E4M3 "F8_E4M3"
#define ST_F8_E8M0 "F8_E8M0"
// Note: Complex numbers aren't in the spec yet so this could change -
// https://github.com/huggingface/safetensors/issues/389
@@ -97,6 +98,8 @@ Dtype dtype_from_safetensor_str(std::string_view str) {
return complex64;
} else if (str == ST_F8_E4M3) {
return uint8;
} else if (str == ST_F8_E8M0) {
return uint8;
} else {
std::ostringstream msg;
msg << "[safetensor] unsupported dtype" << str;
+167 -1
View File
@@ -705,4 +705,170 @@ array solve_triangular(
return matmul(a_inv, b, s);
}
} // namespace mlx::core::linalg
void validate_det(
const array& a,
const StreamOrDevice& stream,
const std::string& fname) {
check_cpu_stream(stream, fname);
if (issubdtype(a.dtype(), complexfloating)) {
throw std::invalid_argument(fname + " Complex inputs are not supported.");
}
if (a.ndim() < 2) {
std::ostringstream msg;
msg << fname
<< " Arrays must have >= 2 dimensions. Received array "
"with "
<< a.ndim() << " dimensions.";
throw std::invalid_argument(msg.str());
}
if (a.shape(-1) != a.shape(-2)) {
throw std::invalid_argument(fname + " Only defined for square matrices.");
}
}
array det_raw_small(const array& a, StreamOrDevice s) {
int n = a.shape(-1);
// Empty 0x0 matrix: determinant is the empty product = 1
if (n == 0) {
Shape out_shape(a.shape().begin(), a.shape().end() - 2);
return broadcast_to(array(1.0f, a.dtype()), std::move(out_shape), s);
}
// Helper to extract a[..., i, j] from the last two dims
auto elem = [&](int i, int j) {
auto starts = Shape(a.ndim(), 0);
auto stops = a.shape();
starts[a.ndim() - 2] = i;
stops[a.ndim() - 2] = i + 1;
starts[a.ndim() - 1] = j;
stops[a.ndim() - 1] = j + 1;
return squeeze(squeeze(slice(a, starts, stops, s), -1, s), -1, s);
};
if (n == 1) {
return elem(0, 0);
} else if (n == 2) {
return subtract(
multiply(elem(0, 0), elem(1, 1), s),
multiply(elem(0, 1), elem(1, 0), s),
s);
} else {
// 3x3: a00*(a11*a22 - a12*a21) - a01*(a10*a22 - a12*a20) + a02*(a10*a21 -
// a11*a20)
auto a00 = elem(0, 0), a01 = elem(0, 1), a02 = elem(0, 2);
auto a10 = elem(1, 0), a11 = elem(1, 1), a12 = elem(1, 2);
auto a20 = elem(2, 0), a21 = elem(2, 1), a22 = elem(2, 2);
return add(
subtract(
multiply(
a00,
subtract(multiply(a11, a22, s), multiply(a12, a21, s), s),
s),
multiply(
a01,
subtract(multiply(a10, a22, s), multiply(a12, a20, s), s),
s),
s),
multiply(
a02, subtract(multiply(a10, a21, s), multiply(a11, a20, s), s), s),
s);
}
}
std::pair<array, array> slogdet_impl(const array& input, StreamOrDevice s) {
int n = input.shape(-1);
auto dtype = input.dtype();
// Small-matrix fast path
if (n <= 3) {
auto raw = det_raw_small(input, s);
auto abs_raw = abs(raw, s);
auto sgn = sign(raw, s);
auto logabs = log(abs_raw, s);
return std::make_pair(sgn, logabs);
}
// General LU-based path
auto [LU, pivots] = lu_factor(input, s);
// Extract diagonal of U
auto diag = diagonal(LU, 0, -2, -1, s);
// Permutation parity: count positions where pivot[i] != i
int k = std::min(input.shape(-2), input.shape(-1));
auto iota = arange(0, k, uint32, s);
auto parity = astype(
sum(not_equal(pivots, iota, s),
/* axis = */ -1,
/* keepdims = */ false,
s),
int32,
s);
// Count negative diagonal elements
auto num_neg = astype(
sum(less(diag, array(0.0f, dtype), s),
/* axis = */ -1,
/* keepdims = */ false,
s),
int32,
s);
// sign = (-1)^(parity + num_neg)
auto total = add(parity, num_neg, s);
auto sign_val = astype(
subtract(
array(1, int32),
multiply(array(2, int32), remainder(total, array(2, int32), s), s),
s),
dtype,
s);
// logabsdet = sum(log(abs(diag)))
auto logabsdet =
sum(log(abs(diag, s), s), /* axis = */ -1, /* keepdims = */ false, s);
// Handle singular matrices: any zero on diagonal
auto is_zero =
any(equal(diag, array(0.0f, dtype), s),
/* axis = */ -1,
/* keepdims = */ false,
s);
sign_val = where(is_zero, array(0.0f, dtype), sign_val, s);
logabsdet = where(
is_zero,
array(-std::numeric_limits<float>::infinity(), dtype),
logabsdet,
s);
return std::make_pair(sign_val, logabsdet);
}
std::pair<array, array> slogdet(const array& a, StreamOrDevice s /* = {} */) {
validate_det(a, s, "[linalg::slogdet]");
auto dtype = at_least_float(a.dtype());
auto input = astype(a, dtype, s);
return slogdet_impl(input, s);
}
array det(const array& a, StreamOrDevice s /* = {} */) {
validate_det(a, s, "[linalg::det]");
auto dtype = at_least_float(a.dtype());
auto input = astype(a, dtype, s);
int n = input.shape(-1);
// Small-matrix fast path: compute directly, skip log/exp round-trip
if (n <= 3) {
return det_raw_small(input, s);
}
// General case: det = sign * exp(logabsdet)
auto [sign_val, logabsdet] = slogdet_impl(input, s);
return multiply(sign_val, exp(logabsdet, s), s);
}
} // namespace mlx::core::linalg
+4
View File
@@ -112,4 +112,8 @@ eigvalsh(const array& a, std::string UPLO = "L", StreamOrDevice s = {});
MLX_API std::pair<array, array>
eigh(const array& a, std::string UPLO = "L", StreamOrDevice s = {});
MLX_API array det(const array& a, StreamOrDevice s = {});
MLX_API std::pair<array, array> slogdet(const array& a, StreamOrDevice s = {});
} // namespace mlx::core::linalg
+2 -2
View File
@@ -5,8 +5,8 @@
#include "mlx/api.h"
#define MLX_VERSION_MAJOR 0
#define MLX_VERSION_MINOR 31
#define MLX_VERSION_PATCH 2
#define MLX_VERSION_MINOR 32
#define MLX_VERSION_PATCH 0
#define MLX_VERSION_NUMERIC \
(100000 * MLX_VERSION_MAJOR + 1000 * MLX_VERSION_MINOR + MLX_VERSION_PATCH)
+12 -3
View File
@@ -1,5 +1,8 @@
// Copyright © 2024 Apple Inc.
#include <limits>
#include <sstream>
#include <nanobind/stl/complex.h>
#include "python/src/convert.h"
@@ -15,9 +18,15 @@ enum PyScalarT {
};
int check_shape_dim(int64_t dim) {
if (dim > std::numeric_limits<int>::max()) {
throw std::invalid_argument(
"Shape dimension falls outside supported `int` range.");
if (dim > std::numeric_limits<int>::max() ||
dim < std::numeric_limits<int>::min()) {
std::ostringstream msg;
msg << "Shape dimension " << dim << " is outside the supported range ["
<< std::numeric_limits<int>::min() << ", "
<< std::numeric_limits<int>::max()
<< "]. MLX currently uses 32-bit integers for shape dimensions.";
PyErr_SetString(PyExc_OverflowError, msg.str().c_str());
nb::detail::raise_python_error();
}
return static_cast<int>(dim);
}
+4
View File
@@ -76,3 +76,7 @@ nb::object tolist(mx::array& a);
mx::array create_array(nb::object v, std::optional<mx::Dtype> t);
mx::array array_from_list(nb::list pl, std::optional<mx::Dtype> dtype);
mx::array array_from_list(nb::tuple pl, std::optional<mx::Dtype> dtype);
// Narrow a Python-side shape dimension (int64) to a C++ mx::ShapeElem (int32),
// raising a clear error if the value would overflow.
int check_shape_dim(int64_t dim);
+73
View File
@@ -660,4 +660,77 @@ void init_linalg(nb::module_& parent_module) {
Returns:
array: The unique solution to the system ``AX = B``.
)pbdoc");
m.def(
"det",
&mx::linalg::det,
"a"_a,
nb::kw_only(),
"stream"_a = nb::none(),
nb::sig(
"def det(a: array, *, stream: Union[None, Stream, Device] = None) -> array"),
R"pbdoc(
Compute the determinant of a square matrix.
This function supports arrays with at least 2 dimensions. When the
input has more than two dimensions, the determinant is computed for
each matrix in the last two dimensions.
Args:
a (array): Input array.
stream (Stream, optional): Stream or device. Defaults to ``None``
in which case the default stream of the default device is used.
Returns:
array: The determinant(s) of the input matrix (matrices).
Example:
>>> A = mx.array([[1., 2.], [3., 4.]])
>>> mx.linalg.det(A, stream=mx.cpu)
array(-2, dtype=float32)
)pbdoc");
m.def(
"slogdet",
[](const mx::array& a, mx::StreamOrDevice s) {
auto result = mx::linalg::slogdet(a, s);
return nb::make_tuple(result.first, result.second);
},
"a"_a,
nb::kw_only(),
"stream"_a = nb::none(),
nb::sig(
"def slogdet(a: array, *, stream: Union[None, Stream, Device] = None) -> Tuple[array, array]"),
R"pbdoc(
Compute the sign and natural log of the absolute value of the
determinant of a square matrix.
This function supports arrays with at least 2 dimensions. When the
input has more than two dimensions, the sign and log-absolute-determinant
are computed for each matrix in the last two dimensions.
For a singular matrix, ``sign`` is 0 and ``logabsdet`` is ``-inf``.
The determinant can be reconstructed as ``det = sign * exp(logabsdet)``.
This is more numerically stable than computing the determinant directly
for matrices with large or small determinants.
Args:
a (array): Input array.
stream (Stream, optional): Stream or device. Defaults to ``None``
in which case the default stream of the default device is used.
Returns:
tuple(array, array): The ``sign`` and ``logabsdet`` of the
determinant. ``sign`` is -1, 0, or +1. ``logabsdet`` is the
natural log of the absolute value of the determinant.
Example:
>>> A = mx.array([[1., 2.], [3., 4.]])
>>> sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
>>> sign
array(-1, dtype=float32)
>>> logabsdet
array(0.693147, dtype=float32)
)pbdoc");
}
+14 -18
View File
@@ -15,6 +15,7 @@
#include "mlx/einsum.h"
#include "mlx/ops.h"
#include "mlx/utils.h"
#include "python/src/convert.h"
#include "python/src/load.h"
#include "python/src/small_vector.h"
#include "python/src/utils.h"
@@ -45,6 +46,13 @@ double scalar_to_double(Scalar s) {
}
}
mx::Shape to_shape(const nb::object& shape) {
if (nb::isinstance<nb::int_>(shape)) {
return {check_shape_dim(nb::cast<int64_t>(shape))};
}
return nb::cast<mx::Shape>(shape);
}
void init_ops(nb::module_& m) {
m.def(
"reshape",
@@ -1702,15 +1710,11 @@ void init_ops(nb::module_& m) {
)pbdoc");
m.def(
"full",
[](const std::variant<int, mx::Shape>& shape,
[](const nb::object& shape,
const ScalarOrArray& vals,
std::optional<mx::Dtype> dtype,
mx::StreamOrDevice s) {
if (auto pv = std::get_if<int>(&shape); pv) {
return mx::full({*pv}, to_array(vals, dtype), s);
} else {
return mx::full(std::get<mx::Shape>(shape), to_array(vals, dtype), s);
}
return mx::full(to_shape(shape), to_array(vals, dtype), s);
},
"shape"_a,
"vals"_a,
@@ -1736,15 +1740,11 @@ void init_ops(nb::module_& m) {
)pbdoc");
m.def(
"zeros",
[](const std::variant<int, mx::Shape>& shape,
[](const nb::object& shape,
std::optional<mx::Dtype> dtype,
mx::StreamOrDevice s) {
auto t = dtype.value_or(mx::float32);
if (auto pv = std::get_if<int>(&shape); pv) {
return mx::zeros({*pv}, t, s);
} else {
return mx::zeros(std::get<mx::Shape>(shape), t, s);
}
return mx::zeros(to_shape(shape), t, s);
},
"shape"_a,
"dtype"_a.none() = mx::float32,
@@ -1802,15 +1802,11 @@ void init_ops(nb::module_& m) {
)pbdoc");
m.def(
"ones",
[](const std::variant<int, mx::Shape>& shape,
[](const nb::object& shape,
std::optional<mx::Dtype> dtype,
mx::StreamOrDevice s) {
auto t = dtype.value_or(mx::float32);
if (auto pv = std::get_if<int>(&shape); pv) {
return mx::ones({*pv}, t, s);
} else {
return mx::ones(std::get<mx::Shape>(shape), t, s);
}
return mx::ones(to_shape(shape), t, s);
},
"shape"_a,
"dtype"_a.none() = mx::float32,
+41 -8
View File
@@ -2,6 +2,11 @@
#pragma once
#include <cstdint>
#include <limits>
#include <sstream>
#include <type_traits>
#include "mlx/small_vector.h"
#include <nanobind/stl/detail/nb_list.h>
@@ -14,11 +19,19 @@ struct type_caster<mlx::core::SmallVector<Type, Size, Alloc>> {
using List = mlx::core::SmallVector<Type, Size, Alloc>;
using Caster = make_caster<Type>;
// For narrow integer element types we fetch each element through a wider
// integer caster so we can emit a clean OverflowError on overflow instead of
// nanobind's generic "incompatible function arguments" TypeError.
static constexpr bool kNarrowInt = std::is_integral_v<Type> &&
!std::is_same_v<Type, bool> && (sizeof(Type) < sizeof(int64_t));
NB_TYPE_CASTER(
List,
const_name("tuple[") + make_caster<Type>::Name + const_name(", ...]"))
bool from_python(handle src, uint8_t flags, cleanup_list* cleanup) noexcept {
// Not noexcept: on overflow of a narrow integer element we raise
// OverflowError so nanobind surfaces a clean error to the user.
bool from_python(handle src, uint8_t flags, cleanup_list* cleanup) {
size_t size;
PyObject* temp;
@@ -29,19 +42,39 @@ struct type_caster<mlx::core::SmallVector<Type, Size, Alloc>> {
value.clear();
value.reserve(size);
Caster caster;
bool success = o != nullptr;
flags = flags_for_local_caster<Type>(flags);
for (size_t i = 0; i < size; ++i) {
if (!caster.from_python(o[i], flags, cleanup) ||
!caster.template can_cast<Type>()) {
success = false;
break;
if constexpr (kNarrowInt) {
make_caster<int64_t> wide;
if (!wide.from_python(o[i], flags, cleanup) ||
!wide.template can_cast<int64_t>()) {
success = false;
break;
}
int64_t v = wide.operator cast_t<int64_t>();
if (v > std::numeric_limits<Type>::max() ||
v < std::numeric_limits<Type>::min()) {
std::ostringstream msg;
msg << "Integer value " << v << " is outside the supported range ["
<< static_cast<int64_t>(std::numeric_limits<Type>::min()) << ", "
<< static_cast<int64_t>(std::numeric_limits<Type>::max()) << "].";
Py_XDECREF(temp);
PyErr_SetString(PyExc_OverflowError, msg.str().c_str());
raise_python_error();
}
value.push_back(static_cast<Type>(v));
} else {
Caster caster;
if (!caster.from_python(o[i], flags, cleanup) ||
!caster.template can_cast<Type>()) {
success = false;
break;
}
value.push_back(caster.operator cast_t<Type>());
}
value.push_back(caster.operator cast_t<Type>());
}
Py_XDECREF(temp);
-17
View File
@@ -1,22 +1,5 @@
cuda_skip = {
# Lapack ops NYI
"TestLinalg.test_cholesky",
"TestLinalg.test_cholesky_inv",
"TestLinalg.test_eig",
"TestLinalg.test_eigh",
"TestLinalg.test_inverse",
"TestVmap.test_vmap_inverse",
"TestLinalg.test_lu",
"TestLinalg.test_lu_factor",
"TestLinalg.test_pseudo_inverse",
"TestLinalg.test_qr_factorization",
"TestInit.test_orthogonal",
"TestLinalg.test_svd_decomposition",
"TestVmap.test_vmap_svd",
"TestLinalg.test_tri_inverse",
# Quantization NYI
"TestQuantized.test_gather_matmul_grad",
"TestQuantized.test_gather_qmm",
"TestQuantized.test_gather_qmm_sorted",
"TestQuantized.test_gather_qmm_grad",
}
+29 -2
View File
@@ -767,10 +767,13 @@ class TestArray(mlx_tests.MLXTestCase):
def test_array_np_shape_dim_check(self):
a_npy = np.empty(2**31, dtype=np.bool_)
with self.assertRaises(ValueError) as e:
with self.assertRaises(OverflowError) as e:
mx.array(a_npy)
self.assertEqual(
str(e.exception), "Shape dimension falls outside supported `int` range."
str(e.exception),
"Shape dimension 2147483648 is outside the supported range "
"[-2147483648, 2147483647]. MLX currently uses 32-bit integers "
"for shape dimensions.",
)
def test_dtype_promotion(self):
@@ -1547,6 +1550,30 @@ class TestArray(mlx_tests.MLXTestCase):
a = a.at[:, :, :].add(update)
self.assertEqualArray(a, update)
def test_slice_update_contiguous_2d(self):
for shape in [(32, 32), (64, 64), (17, 33), (128, 128)]:
for upd_shape in [(16, 16), (8, 8), (4, 4)]:
if upd_shape[0] > shape[0] or upd_shape[1] > shape[1]:
continue
y = mx.zeros(shape)
x = mx.random.normal(upd_shape)
z = y.at[: upd_shape[0], : upd_shape[1]].add(x)
diff = z[: upd_shape[0], : upd_shape[1]] - x
self.assertTrue(mx.allclose(diff, mx.zeros_like(diff)))
# Test non-zero offset slice update
y = mx.zeros((32, 32))
x = mx.random.normal((16, 16))
z = y.at[16:, 16:].add(x)
self.assertTrue(mx.allclose(z[16:, 16:] - x, mx.zeros_like(x)))
# Test with size divisible by 4, 2, and odd
for cols in [32, 18, 15]:
y = mx.zeros((32, cols))
x = mx.random.normal((16, cols))
z = y.at[:16, :].add(x)
self.assertTrue(mx.allclose(z[:16, :] - x, mx.zeros_like(x)))
def test_slice_negative_step(self):
a_np = np.arange(20)
a_mx = mx.array(a_np)
+255
View File
@@ -520,6 +520,19 @@ class TestLinalg(mlx_tests.MLXTestCase):
P, L, U = mx.linalg.lu(a, stream=mx.cpu)
self.assertTrue(mx.allclose(L[P, :] @ U, a))
# Test singular matrix (should not throw)
a = mx.array(
[
[1.0, 2.0, 3.0, 4.0],
[2.0, 4.0, 6.0, 8.0],
[0.0, 1.0, 1.0, 0.0],
[1.0, 0.0, 0.0, 1.0],
]
)
P, L, U = mx.linalg.lu(a, stream=mx.cpu)
L_permuted = mx.take_along_axis(L, P[..., None], axis=-2)
self.assertTrue(mx.allclose(L_permuted @ U, a))
def test_lu_factor(self):
mx.random.seed(7)
@@ -616,6 +629,248 @@ class TestLinalg(mlx_tests.MLXTestCase):
expected = np.linalg.solve(a, b)
self.assertTrue(np.allclose(result, expected))
def test_det(self):
# 1x1 fast path
A = mx.array([[5.0]])
self.assertTrue(np.allclose(mx.linalg.det(A, stream=mx.cpu), 5.0))
# 2x2 fast path
A = mx.array([[1.0, 2.0], [3.0, 4.0]])
d = mx.linalg.det(A, stream=mx.cpu)
self.assertTrue(np.allclose(d, -2.0))
# 3x3 fast path
A = mx.array([[1.0, 2.0, 3.0], [0.0, 1.0, 4.0], [5.0, 6.0, 0.0]])
d = mx.linalg.det(A, stream=mx.cpu)
expected = np.linalg.det(np.array(A))
self.assertTrue(np.allclose(d, expected, atol=1e-5))
# 4x4 LU path: compare with numpy
np.random.seed(42)
A_np = np.random.randn(4, 4).astype(np.float32)
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-4))
# 5x5 LU path
A_np = np.random.randn(5, 5).astype(np.float32)
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-4))
# Identity matrix
A = mx.eye(5)
self.assertTrue(np.allclose(mx.linalg.det(A, stream=mx.cpu), 1.0))
# Batched: (3, 4, 4)
A_np = np.random.randn(3, 4, 4).astype(np.float32)
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-4))
# Multi-batch: (2, 3, 3, 3)
A_np = np.random.randn(2, 3, 3, 3).astype(np.float32)
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-4))
# Integer input auto-promotes to float
A = mx.array([[1, 2], [3, 4]])
d = mx.linalg.det(A, stream=mx.cpu)
self.assertTrue(np.allclose(d, -2.0))
# float64
A_np = np.random.randn(4, 4).astype(np.float64)
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-10))
# Singular 4x4 matrix (LU path): det should be 0
A = mx.array(
[
[1.0, 2.0, 3.0, 4.0],
[2.0, 4.0, 6.0, 8.0],
[0.0, 1.0, 1.0, 0.0],
[1.0, 0.0, 0.0, 1.0],
]
)
d = mx.linalg.det(A, stream=mx.cpu)
self.assertTrue(np.allclose(d, 0.0, atol=1e-5))
# Singular 5x5 matrix (LU path)
A_np = np.ones((5, 5), dtype=np.float32)
A_mx = mx.array(A_np)
d = mx.linalg.det(A_mx, stream=mx.cpu)
self.assertTrue(np.allclose(d, 0.0, atol=1e-5))
# Batched singular matrices (LU path)
A_np = np.array([np.diag([1.0, 2.0, 0.0, 3.0]), np.eye(4, dtype=np.float32)])
A_mx = mx.array(A_np)
d_mx = mx.linalg.det(A_mx, stream=mx.cpu)
d_np = np.linalg.det(A_np)
self.assertTrue(np.allclose(d_mx, d_np, atol=1e-5))
# Empty 0x0 matrix: det is the empty product = 1
d = mx.linalg.det(mx.zeros((0, 0)), stream=mx.cpu)
self.assertEqual(d.shape, ())
self.assertEqual(float(d), 1.0)
# Batched empty matrices: shape preserves batch dims
d = mx.linalg.det(mx.zeros((3, 0, 0)), stream=mx.cpu)
self.assertTrue(np.allclose(d, np.linalg.det(np.zeros((3, 0, 0)))))
# Error: non-square
with self.assertRaises(ValueError):
mx.linalg.det(mx.array([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]]), stream=mx.cpu)
# Error: 1D
with self.assertRaises(ValueError):
mx.linalg.det(mx.array([1.0, 2.0]), stream=mx.cpu)
# Error: complex unsupported (small-matrix path)
with self.assertRaises(ValueError):
mx.linalg.det(mx.array([[1.0 + 1j, 2.0], [3.0, 4.0]]), stream=mx.cpu)
# Error: complex unsupported (LU path)
with self.assertRaises(ValueError):
mx.linalg.det(mx.eye(4).astype(mx.complex64), stream=mx.cpu)
def test_slogdet(self):
# 2x2: det = -2 => sign = -1, logabsdet = log(2)
A = mx.array([[1.0, 2.0], [3.0, 4.0]])
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertTrue(np.allclose(sign, -1.0))
self.assertTrue(np.allclose(logabsdet, np.log(2.0), atol=1e-5))
# Identity: sign = 1, logabsdet = 0
A = mx.eye(4)
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertTrue(np.allclose(sign, 1.0))
self.assertTrue(np.allclose(logabsdet, 0.0, atol=1e-6))
# Compare with numpy for random matrices
np.random.seed(42)
for n in [1, 2, 3, 4, 5]:
A_np = np.random.randn(n, n).astype(np.float32)
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(A_np)
with self.subTest(n=n):
self.assertTrue(np.allclose(sign_mx, sign_np, atol=1e-5))
self.assertTrue(np.allclose(logabs_mx, logabs_np, atol=1e-4))
# Singular matrix 2x2 (fast path): sign = 0, logabsdet = -inf
A = mx.array([[1.0, 2.0], [2.0, 4.0]])
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertEqual(float(sign), 0.0)
self.assertEqual(float(logabsdet), float("-inf"))
# Singular 4x4 matrix (LU path): sign = 0, logabsdet = -inf
A = mx.array(
[
[1.0, 2.0, 3.0, 4.0],
[2.0, 4.0, 6.0, 8.0],
[0.0, 1.0, 1.0, 0.0],
[1.0, 0.0, 0.0, 1.0],
]
)
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertEqual(float(sign), 0.0)
self.assertEqual(float(logabsdet), float("-inf"))
# Singular 5x5 matrix (LU path): all-ones matrix
A = mx.array(np.ones((5, 5), dtype=np.float32))
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertEqual(float(sign), 0.0)
self.assertEqual(float(logabsdet), float("-inf"))
# Batched with mix of singular and non-singular (LU path)
A_np = np.array([np.diag([1.0, 2.0, 0.0, 3.0]), np.eye(4, dtype=np.float32)])
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(A_np)
self.assertTrue(np.allclose(sign_mx, sign_np, atol=1e-5))
# Check -inf for singular, 0.0 for identity
self.assertEqual(float(logabs_mx[0]), float("-inf"))
self.assertTrue(np.allclose(logabs_mx[1], 0.0, atol=1e-6))
# Batched
A_np = np.random.randn(3, 4, 4).astype(np.float32)
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(A_np)
self.assertTrue(np.allclose(sign_mx, sign_np, atol=1e-5))
self.assertTrue(np.allclose(logabs_mx, logabs_np, atol=1e-4))
# Multi-batch
A_np = np.random.randn(2, 3, 3, 3).astype(np.float32)
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(A_np)
self.assertTrue(np.allclose(sign_mx, sign_np, atol=1e-5))
self.assertTrue(np.allclose(logabs_mx, logabs_np, atol=1e-4))
# Numerical stability: large matrix where det overflows
# 0.1 * I_100 has det = 0.1^100 which underflows in float32
# but slogdet should give sign=1, logabsdet = 100*log(0.1)
n = 100
A = mx.array(0.1) * mx.eye(n)
sign, logabsdet = mx.linalg.slogdet(A, stream=mx.cpu)
self.assertTrue(np.allclose(sign, 1.0))
self.assertTrue(np.allclose(logabsdet, n * np.log(0.1), atol=1e-3))
# Verify det = sign * exp(logabsdet) for non-singular cases
A_np = np.random.randn(5, 5).astype(np.float32)
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
det_mx = mx.linalg.det(A_mx, stream=mx.cpu)
reconstructed = float(sign_mx) * np.exp(float(logabs_mx))
self.assertTrue(np.allclose(float(det_mx), reconstructed, rtol=1e-4))
# float64
A_np = np.random.randn(4, 4).astype(np.float64)
A_mx = mx.array(A_np)
sign_mx, logabs_mx = mx.linalg.slogdet(A_mx, stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(A_np)
self.assertTrue(np.allclose(sign_mx, sign_np))
self.assertTrue(np.allclose(logabs_mx, logabs_np, atol=1e-10))
# Empty 0x0 matrix: sign = 1, logabsdet = 0 (empty product)
sign, logabsdet = mx.linalg.slogdet(mx.zeros((0, 0)), stream=mx.cpu)
self.assertEqual(sign.shape, ())
self.assertEqual(logabsdet.shape, ())
self.assertEqual(float(sign), 1.0)
self.assertEqual(float(logabsdet), 0.0)
# Batched empty matrices
sign, logabsdet = mx.linalg.slogdet(mx.zeros((3, 0, 0)), stream=mx.cpu)
sign_np, logabs_np = np.linalg.slogdet(np.zeros((3, 0, 0)))
self.assertTrue(np.allclose(sign, sign_np))
self.assertTrue(np.allclose(logabsdet, logabs_np))
# Error: non-square
with self.assertRaises(ValueError):
mx.linalg.slogdet(
mx.array([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]]), stream=mx.cpu
)
# Error: 1D
with self.assertRaises(ValueError):
mx.linalg.slogdet(mx.array([1.0, 2.0]), stream=mx.cpu)
# Error: complex unsupported (small-matrix path)
with self.assertRaises(ValueError):
mx.linalg.slogdet(mx.array([[1.0 + 1j, 2.0], [3.0, 4.0]]), stream=mx.cpu)
# Error: complex unsupported (LU path)
with self.assertRaises(ValueError):
mx.linalg.slogdet(mx.eye(4).astype(mx.complex64), stream=mx.cpu)
if __name__ == "__main__":
mlx_tests.MLXTestRunner()
+46
View File
@@ -92,6 +92,52 @@ class TestOps(mlx_tests.MLXTestCase):
self.assertEqual(y.dtype, t)
self.assertTrue(mx.array_equal(y, x))
def test_shape_overflow_error(self):
# Shape dimensions that don't fit in int32 should raise a clear
# OverflowError that names the offending value, rather than a generic
# "incompatible function arguments" TypeError. The overflow check
# lives in the mx::Shape type caster, so it applies to every op that
# takes a shape. See issue #2681.
too_big = 2**31
# Array creation ops — also exercise the scalar shape path.
for ctor in (mx.zeros, mx.ones):
with self.assertRaises(OverflowError) as cm:
ctor(too_big)
self.assertIn(str(too_big), str(cm.exception))
with self.assertRaises(OverflowError) as cm:
ctor([too_big])
self.assertIn(str(too_big), str(cm.exception))
with self.assertRaises(OverflowError) as cm:
mx.full(too_big, 0.0)
self.assertIn(str(too_big), str(cm.exception))
with self.assertRaises(OverflowError) as cm:
mx.full([too_big], 0.0)
self.assertIn(str(too_big), str(cm.exception))
# Other shape-taking ops should surface the same clean error.
a = mx.zeros(4)
with self.assertRaises(OverflowError) as cm:
mx.reshape(a, [too_big])
self.assertIn(str(too_big), str(cm.exception))
with self.assertRaises(OverflowError) as cm:
mx.broadcast_to(a, [too_big, 1])
self.assertIn(str(too_big), str(cm.exception))
# Negative overflow (< int32 min) is caught too.
too_negative = -(2**31) - 1
with self.assertRaises(OverflowError) as cm:
mx.zeros([too_negative])
self.assertIn(str(too_negative), str(cm.exception))
# Shapes that fit in int32 still go through unchanged.
self.assertEqual(mx.zeros(4).shape, (4,))
self.assertEqual(mx.zeros((2, 3)).shape, (2, 3))
self.assertEqual(mx.ones([2, 3]).shape, (2, 3))
self.assertEqual(mx.full((2, 3), 1.5).tolist(), [[1.5] * 3] * 2)
def test_scalar_inputs(self):
# Check combinations of python types
a = mx.add(False, True)
+25 -1
View File
@@ -470,6 +470,30 @@ class TestQuantized(mlx_tests.MLXTestCase):
self.assertEqual(y_q.shape, y_hat.shape)
self.assertLess((y_q - y_hat).abs().max(), 2e-3)
def test_qvm_splitk_multi_row(self):
# Test qvm split_k with M > 1 to ensure the x row stride is correct
key = mx.random.key(0)
k1, k2 = mx.random.split(key)
tests = product(
[64, 32], # group_size
[4, 8], # bits
[128], # out dim (N)
[2048, 4096], # in dim (K) >= 1024 to trigger split_k
[2, 3], # M (multiple rows)
)
for group_size, bits, N, K, M in tests:
with self.subTest(M=M, K=K, N=N, group_size=group_size, bits=bits):
x = 1e-1 * mx.random.normal(shape=(M, K), key=k1)
w = 1e-1 * mx.random.normal(shape=(K, N), key=k2)
w_q, scales, biases = mx.quantize(w, group_size, bits)
w_hat = mx.dequantize(w_q, scales, biases, group_size, bits)
y_q = mx.quantized_matmul(
x, w_q, scales, biases, False, group_size, bits
)
y_hat = x @ w_hat
self.assertEqual(y_q.shape, y_hat.shape)
self.assertLess((y_q - y_hat).abs().max(), 2e-3)
def test_fp_qvm(self):
key = mx.random.key(0)
k1, k2 = mx.random.split(key)
@@ -1046,7 +1070,7 @@ class TestQuantized(mlx_tests.MLXTestCase):
y3 = scatter_unsort(y3, inv_order, indices.shape)
y4 = scatter_unsort(y4, inv_order, indices.shape)
tol = 1.5e-5 if (dtype == mx.float32) else 2.5e-4
tol = 1.5e-5 if (dtype == mx.float32) else 1e-3
self.assertLess((y1 - y2).abs().max(), tol)
self.assertLess((y1 - y3).abs().max(), tol)
+65
View File
@@ -637,3 +637,68 @@ TEST_CASE("test solve_triangluar") {
expected = array({-3., 2., 3.});
CHECK(allclose(expected, result).item<bool>());
}
TEST_CASE("test det") {
// 1x1 fast path
{
array a = array({5.0f}, {1, 1});
auto d = det(a, Device::cpu);
CHECK_EQ(d.item<float>(), doctest::Approx(5.0f));
}
// 2x2 fast path: det([[1,2],[3,4]]) = -2
{
array a = array({1.0f, 2.0f, 3.0f, 4.0f}, {2, 2});
auto d = det(a, Device::cpu);
CHECK_EQ(d.item<float>(), doctest::Approx(-2.0f));
}
// 3x3 fast path: det([[1,2,3],[0,1,4],[5,6,0]]) = 1
{
array a =
array({1.0f, 2.0f, 3.0f, 0.0f, 1.0f, 4.0f, 5.0f, 6.0f, 0.0f}, {3, 3});
auto d = det(a, Device::cpu);
CHECK_EQ(d.item<float>(), doctest::Approx(1.0f));
}
// 4x4 LU path: identity matrix det = 1
{
array a = eye(4);
auto d = det(a, Device::cpu);
CHECK_EQ(d.item<float>(), doctest::Approx(1.0f));
}
// Non-square should throw
CHECK_THROWS(
det(array({1.0f, 2.0f, 3.0f, 4.0f, 5.0f, 6.0f}, {2, 3}), Device::cpu));
// 1D should throw
CHECK_THROWS(det(array({1.0f, 2.0f}), Device::cpu));
}
TEST_CASE("test slogdet") {
// 2x2: det = -2, so sign = -1, logabsdet = log(2)
{
array a = array({1.0f, 2.0f, 3.0f, 4.0f}, {2, 2});
auto [s, logabs] = slogdet(a, Device::cpu);
CHECK_EQ(s.item<float>(), doctest::Approx(-1.0f));
CHECK_EQ(logabs.item<float>(), doctest::Approx(std::log(2.0f)));
}
// Identity: sign = 1, logabsdet = 0
{
array a = eye(4);
auto [s, logabs] = slogdet(a, Device::cpu);
CHECK_EQ(s.item<float>(), doctest::Approx(1.0f));
CHECK_EQ(logabs.item<float>(), doctest::Approx(0.0f));
}
// Singular: sign = 0, logabsdet = -inf
{
array a = array({1.0f, 2.0f, 2.0f, 4.0f}, {2, 2});
auto [s, logabs] = slogdet(a, Device::cpu);
CHECK_EQ(s.item<float>(), 0.0f);
CHECK(std::isinf(logabs.item<float>()));
CHECK(logabs.item<float>() < 0);
}
}