cuda: fix sparse F16 compilation and temporary buffer leaks

ggml-cuda.cu

@@ -4510,7 +4510,7 @@ static __global__ void dequantize_mul_mat_axpy_sparse_pro(const void * __restric
 
 #pragma unroll 8
     for(int gpu_row = wid; gpu_row < nrows; gpu_row += AXPY_BLOCK_Y) {        
-        if(y_tmp[gpu_row] == 0.0) continue;
+        if(float(y_tmp[gpu_row]) == 0.0f) continue;
 
         const int ib = ((gpu_row + row_offset)*ncols + col) / qk; // x block index
 
@@ -4550,7 +4550,7 @@ static __global__ void dequantize_mul_mat_axpy_sparse(const void * __restrict__
     const int tid = threadIdx.x;
     short *d = (short *)((char *)vx + ncols * gpu_row * 2);
 
-    if (y[row] == 0)
+    if (float(y[row]) == 0.0f)
         return;
     if (idx[row] < dev_sparse_threshold) {
         return;
@@ -4755,7 +4755,7 @@ static __global__ void dequantize_mul_mat_vec_sparse(const void * __restrict__ v
 }
 
 template <int qk, int qr, dequantize_kernel_t dequantize_kernel>
-static __global__ void dequantize_mul_mat_batch_sparse(const void * __restrict__ vx, const dfloat * __restrict__ y, float * __restrict__ dst, const int ncols, const int nrows, int src1_cols, int dst_ne0,int * lst, float * idx) {
+static __global__ void dequantize_mul_mat_batch_sparse(const void * __restrict__ vx, const float * __restrict__ y, float * __restrict__ dst, const int ncols, const int nrows, int src1_cols, int dst_ne0,int * lst, float * idx) {
     // qk = quantized weights per x block
     // qr = number of quantized weights per data value in x block
     const int gpu_row = blockIdx.y*blockDim.y + threadIdx.y;
@@ -4772,7 +4772,7 @@ static __global__ void dequantize_mul_mat_batch_sparse(const void * __restrict__
     const int vals_per_iter = iter_stride / WARP_SIZE; // num quantized vals per thread and i iter
     const int y_offset = qr == 1 ? 1 : qk/2;
     float * loop_idx = idx;;
-    dfloat * loop_y = (dfloat *)y;
+    const float * loop_y = y;
     float * loop_dst = dst;
 
 
@@ -4811,8 +4811,8 @@ static __global__ void dequantize_mul_mat_batch_sparse(const void * __restrict__
 
                 // matrix multiplication
 
-                tmp += v.x * loop_y[iybs + iqs + j / qr + 0];
-                tmp += v.y * loop_y[iybs + iqs + j / qr + y_offset];
+                tmp += float(v.x) * loop_y[iybs + iqs + j / qr + 0];
+                tmp += float(v.y) * loop_y[iybs + iqs + j / qr + y_offset];
                 // #endif
             }
         }
@@ -5387,7 +5387,7 @@ static void dequantize_mul_mat_vec_q4_0_cuda_sparse(const void * vx, const dfloa
     dequantize_mul_mat_vec_sparse<QK4_0, QR4_0, dequantize_q4_0>
         <<<block_nums, block_dims, 0, stream>>>(vx, y, dst, ncols, nrows, lst, idx);
 }
-static void dequantize_mul_mat_batch_q4_0_cuda_sparse(const void * vx, const dfloat * y, float * dst, const int ncols, const int nrows, int src1_ncols, int dst_ne0, cudaStream_t stream, int *lst, float *idx) {
+static void dequantize_mul_mat_batch_q4_0_cuda_sparse(const void * vx, const float * y, float * dst, const int ncols, const int nrows, int src1_ncols, int dst_ne0, cudaStream_t stream, int *lst, float *idx) {
     GGML_ASSERT(ncols % GGML_CUDA_DMMV_X == 0);
     const int block_num_y = (nrows + GGML_CUDA_MMV_Y - 1) / GGML_CUDA_MMV_Y;
     const dim3 block_nums(1, block_num_y, 1);
@@ -5644,7 +5644,7 @@ static void convert_mul_mat_vec_f16_cuda_sparse(const void * vx, const dfloat *
         <<<block_nums, block_dims, 0, stream>>>(vx, y, dst, ncols, nrows, lst, idx);
 
 }
-static void convert_mul_mat_batch_f16_cuda_sparse(const void * vx, const dfloat * y, float * dst, const int ncols, const int nrows, int src1_ncols, int dst_ne0, cudaStream_t stream, int *lst, float *idx) {
+static void convert_mul_mat_batch_f16_cuda_sparse(const void * vx, const float * y, float * dst, const int ncols, const int nrows, int src1_ncols, int dst_ne0, cudaStream_t stream, int *lst, float *idx) {
     GGML_ASSERT(ncols % GGML_CUDA_DMMV_X == 0);
     const int block_num_y = (nrows + GGML_CUDA_MMV_Y - 1) / GGML_CUDA_MMV_Y;
     const dim3 block_nums(1, block_num_y, 1);
@@ -7235,31 +7235,12 @@ inline void ggml_cuda_op_mul_mat_vec_sparse_q(
     const int64_t src1_padded_row_size, const cudaStream_t & stream) {
 
     const int64_t ne00 = src0->ne[0];
-    const int64_t ne10 = src1->ne[1];
     const int64_t row_diff = row_high - row_low;
 
     float * sparse_idx = static_cast<float *>(ggml_cuda_get_tensor_data(dst->src[2]));
     int32_t * row_lookup = dst->src[3] != NULL ? static_cast<int32_t *>(ggml_cuda_get_tensor_data(dst->src[3])) : NULL;
 
-    // on some GPUs it is faster to convert src1 to half and to use half precision intrinsics
-#ifdef GGML_CUDA_F16
-    size_t ash;
-    dfloat * src1_dfloat = nullptr; // dfloat == half
-
-    bool src1_convert_f16 = src0->type == GGML_TYPE_Q4_0 || src0->type == GGML_TYPE_Q4_1 ||
-        src0->type == GGML_TYPE_Q5_0 || src0->type == GGML_TYPE_Q5_1 ||
-        src0->type == GGML_TYPE_Q8_0 || src0->type == GGML_TYPE_F16;
-
-    if (src1_convert_f16) {
-        src1_dfloat = (half *) ggml_cuda_pool_malloc(ne00*sizeof(half), &ash);
-        ggml_cpy_f32_f16_cuda((const char *) src1_ddf_i, (char *) src1_dfloat, ne00,
-                                ne00, 1, sizeof(float), 0, 0,
-                                ne00, 1, sizeof(half),  0, 0, stream);
-    }
-#else
-    const dfloat * src1_dfloat = (const dfloat *) src1_ddf_i; // dfloat == float, no conversion
-#endif // GGML_CUDA_F16
-
+    // This kernel consumes Q8_1 activations, so no F16 temporary is needed.
     cudaMemsetAsync((void *)dst_dd_i, 0, ggml_nbytes(dst), stream);
     switch (src0->type) {
         case GGML_TYPE_Q4_0:
@@ -7321,6 +7302,12 @@ inline void ggml_cuda_op_mul_mat_vec_sparse_dequantized(
             break;
     }
 
+#ifdef GGML_CUDA_F16
+    if (src1_convert_f16) {
+        ggml_cuda_pool_free(src1_dfloat, ash);
+    }
+#endif // GGML_CUDA_F16
+
     (void) src1;
     (void) dst;
     (void) src1_ddf_i;

tests/CMakeLists.txt

@@ -49,3 +49,17 @@ llama_build_and_test_executable(test-rope.cpp)
 get_filename_component(TEST_TARGET test-c.c NAME_WE)
 add_executable(${TEST_TARGET} test-c.c)
 target_link_libraries(${TEST_TARGET} PRIVATE llama)
+
+if (LLAMA_CUBLAS AND CUDAToolkit_FOUND)
+    add_executable(test-cuda-sparse test-cuda-sparse.cpp)
+    target_link_libraries(test-cuda-sparse PRIVATE ggml CUDA::cudart)
+    foreach(SPARSE_TYPE q4_0 f16)
+        foreach(SPARSE_COLUMNS 1 3)
+            foreach(SPARSE_MAPPED 0 1)
+                set(SPARSE_TEST test-cuda-sparse-${SPARSE_TYPE}-${SPARSE_COLUMNS}-${SPARSE_MAPPED})
+                add_test(NAME ${SPARSE_TEST} COMMAND test-cuda-sparse ${SPARSE_TYPE} ${SPARSE_COLUMNS} ${SPARSE_MAPPED})
+                set_tests_properties(${SPARSE_TEST} PROPERTIES SKIP_RETURN_CODE 77 RUN_SERIAL TRUE LABELS "cuda;sparse" TIMEOUT 120)
+            endforeach()
+        endforeach()
+    endforeach()
+endif()

tests/README-cuda-sparse.md

@@ -0,0 +1,45 @@
+# CUDA sparse operation regression
+
+`test-cuda-sparse` exercises the production `ggml_cuda_compute_forward` dispatcher
+using synthetic tensors. No model weights or Python dependencies are required.
+
+Build with a CUDA toolkit and a supported host C++ compiler. The example targets
+an RTX 5080; change `CMAKE_CUDA_ARCHITECTURES` for a different GPU.
+
+```sh
+cmake -S . -B build-sparse-on \
+  -DLLAMA_CUBLAS=ON -DLLAMA_CUDA_F16=ON \
+  -DCMAKE_CUDA_ARCHITECTURES=120 -DCMAKE_CUDA_STANDARD=17 \
+  -DLLAMA_BUILD_TESTS=ON -DLLAMA_BUILD_EXAMPLES=OFF \
+  -DLLAMA_BUILD_SERVER=OFF -DBUILD_SHARED_LIBS=OFF
+cmake --build build-sparse-on --config Release --target test-cuda-sparse
+ctest --test-dir build-sparse-on -C Release -R '^test-cuda-sparse-' --output-on-failure
+```
+
+Repeat in a separate build directory with `LLAMA_CUDA_F16=OFF` to cover the
+default configuration. The tests are registered only for the CUDA backend,
+run serially, and return the CTest skip code 77 if no usable CUDA device exists.
+
+The eight cases cover Q4_0 and F16 weights, one and three input vectors, and
+direct and mapped output rows. Each checks active, inactive, and mixed masks
+against a CPU dot-product reference using the decoded stored weights. Inputs
+are chosen to avoid conflating Q8_1 input quantization with indexing errors.
+
+After warmup, each case repeats the operation 2048 times with input width 16384
+and 64 stored weight rows. Retained CUDA memory may grow by at most 16 MiB.
+Run on an otherwise idle GPU: `cudaMemGetInfo` can include other processes,
+and the tolerance accommodates allocator granularity and device bookkeeping.
+The original F16-enabled vector paths exceed this tolerance by losing a
+temporary activation buffer on every call. Batch paths are controls.
+
+For a short sanitizer check, invoke the executable directly:
+
+```sh
+compute-sanitizer --tool memcheck --error-exitcode 1 \
+  build-sparse-on/bin/test-cuda-sparse q4_0 1 1 16
+```
+
+Arguments are `q4_0|f16`, input-vector count, mapped-row flag, and optionally
+iteration count. Use `bin/Release/test-cuda-sparse.exe` with a Windows
+multi-configuration generator. Sixteen iterations are a sanitizer smoke test;
+use the default 2048 iterations for the retained-memory regression.

tests/test-cuda-sparse.cpp

@@ -0,0 +1,193 @@
+#include "ggml.h"
+#include "ggml-cuda.h"
+
+#include <cuda_runtime.h>
+
+#include <algorithm>
+#include <cmath>
+#include <cstdio>
+#include <cstdlib>
+#include <cstring>
+#include <vector>
+
+// Exercise the actual CUDA dispatcher without loading a model. Run on an idle GPU:
+// cudaMemGetInfo includes allocations made by other processes on some platforms.
+#define CHECK_CUDA(call) do { \
+    const cudaError_t err = (call); \
+    if (err != cudaSuccess) { \
+        std::fprintf(stderr, "%s: %s\n", #call, cudaGetErrorString(err)); \
+        std::exit(1); \
+    } \
+} while (0)
+
+static void upload(ggml_tensor * tensor) {
+    tensor->backend = GGML_BACKEND_GPU;
+    ggml_cuda_transform_tensor(tensor->data, tensor);
+}
+
+static bool compute(ggml_tensor * output) {
+    ggml_compute_params params = {};
+    params.type = GGML_TASK_COMPUTE;
+    params.nth = 1;
+    if (!ggml_cuda_compute_forward(&params, output)) {
+        std::fprintf(stderr, "Sparse operation did not run on CUDA\n");
+        return false;
+    }
+    return true;
+}
+
+static bool check_output(ggml_tensor * output, const std::vector<float> & weights,
+                         const float * input, const float * mask,
+                         int k, int rows, int columns, bool mapped) {
+    CHECK_CUDA(cudaDeviceSynchronize());
+    ggml_cuda_copy_to_host(output);
+    const float * actual = static_cast<const float *>(output->data);
+    const int out_rows = static_cast<int>(output->ne[0]);
+    std::vector<double> expected(out_rows * columns, 0.0);
+    for (int col = 0; col < columns; ++col) {
+        for (int row = 0; row < rows; ++row) {
+            const int out_row = mapped ? 2 * row + 1 : row;
+            if (mask[col * out_rows + out_row] < 0.5f) {
+                continue;
+            }
+            double sum = 0.0;
+            for (int i = 0; i < k; ++i) {
+                sum += double(weights[row * k + i]) * input[col * k + i];
+            }
+            expected[col * out_rows + out_row] = sum;
+        }
+    }
+    double max_error = 0.0;
+    for (size_t i = 0; i < expected.size(); ++i) {
+        const double error = std::fabs(actual[i] - expected[i]);
+        max_error = std::max(max_error, error);
+        // The F16 vector kernel accumulates in half2; Q4_0 vectors quantize input.
+        if (!std::isfinite(actual[i]) || error > 0.01 + 0.01 * std::fabs(expected[i])) {
+            std::fprintf(stderr, "Output mismatch at %zu: actual=%g expected=%.9g\n",
+                         i, actual[i], expected[i]);
+            return false;
+        }
+        if (expected[i] == 0.0 && mask[i] < 0.5f && actual[i] != 0.0f) {
+            std::fprintf(stderr, "Inactive row was not cleared at %zu\n", i);
+            return false;
+        }
+    }
+    std::printf("max_abs_error=%.9g\n", max_error);
+    return true;
+}
+
+int main(int argc, char ** argv) {
+    if (argc < 4 || argc > 5 ||
+        (std::strcmp(argv[1], "q4_0") != 0 && std::strcmp(argv[1], "f16") != 0)) {
+        std::fprintf(stderr, "Usage: %s q4_0|f16 columns mapped [iterations]\n", argv[0]);
+        return 1;
+    }
+    const ggml_type type = std::strcmp(argv[1], "q4_0") == 0 ? GGML_TYPE_Q4_0 : GGML_TYPE_F16;
+    const int columns = std::atoi(argv[2]);
+    const bool mapped = std::atoi(argv[3]) != 0;
+    const int iterations = argc == 5 ? std::atoi(argv[4]) : 2048;
+    if (columns < 1 || columns > 4 || iterations < 1 || iterations > 16384) {
+        return 1;
+    }
+    int devices = 0;
+    const cudaError_t device_status = cudaGetDeviceCount(&devices);
+    if (device_status != cudaSuccess || devices == 0) {
+        std::fprintf(stderr, "No usable CUDA device: skipping\n");
+        return 77;
+    }
+    CHECK_CUDA(cudaSetDevice(0));
+    const int k = 16384;
+    const int rows = 64;
+    const int out_rows = mapped ? 2 * rows : rows;
+    ggml_init_params init = { 32 * 1024 * 1024, nullptr, false };
+    ggml_context * ctx = ggml_init(init);
+    if (!ctx || !ggml_cublas_loaded()) {
+        std::fprintf(stderr, "Could not initialize ggml CUDA\n");
+        return 1;
+    }
+    ggml_cuda_set_device_constants(0.5f);
+    ggml_tensor * weight = ggml_new_tensor_2d(ctx, type, k, rows);
+    ggml_tensor * input = ggml_new_tensor_2d(ctx, GGML_TYPE_F32, k, columns);
+    ggml_tensor * mask = ggml_new_tensor_2d(ctx, GGML_TYPE_F32, out_rows, columns);
+    ggml_tensor * lookup = mapped ? ggml_new_tensor_1d(ctx, GGML_TYPE_I32, rows) : nullptr;
+    ggml_tensor * output = mapped
+        ? ggml_mul_mat_idx_upscale(ctx, weight, input, mask, lookup, out_rows)
+        : ggml_mul_mat_idx(ctx, weight, input, mask, nullptr);
+
+    std::vector<float> reference_weights(k * rows);
+    for (int row = 0; row < rows; ++row) {
+        for (int i = 0; i < k; ++i) {
+            reference_weights[row * k + i] = float((i * 13 + row * 7) % 17 - 8) / 128.0f;
+        }
+    }
+    const ggml_type_traits_t traits = ggml_internal_get_type_traits(type);
+    traits.from_float_reference(reference_weights.data(), weight->data, k * rows);
+    traits.to_float(weight->data, reference_weights.data(), k * rows);
+    float * x = static_cast<float *>(input->data);
+    for (int col = 0; col < columns; ++col) {
+        for (int i = 0; i < k; ++i) {
+            // Exactly representable inputs make the reference insensitive to Q8_1 rounding.
+            x[col * k + i] = ((i / 32 + col) % 2 ? -1.0f : 1.0f) * (col + 1) / 8.0f;
+        }
+    }
+    if (mapped) {
+        int32_t * indices = static_cast<int32_t *>(lookup->data);
+        for (int row = 0; row < rows; ++row) {
+            indices[row] = 2 * row + 1;
+        }
+    }
+    std::memset(mask->data, 0, ggml_nbytes(mask));
+    std::memset(output->data, 0, ggml_nbytes(output));
+    upload(weight);
+    upload(input);
+    upload(mask);
+    upload(output);
+    if (lookup) {
+        upload(lookup);
+    }
+
+    bool ok = true;
+    float * mask_data = static_cast<float *>(mask->data);
+    // Active -> inactive -> mixed catches stale outputs as well as indexing errors.
+    for (int pattern = 0; pattern < 3; ++pattern) {
+        for (int col = 0; col < columns; ++col) {
+            for (int row = 0; row < out_rows; ++row) {
+                mask_data[col * out_rows + row] = pattern == 0 ? 1.0f :
+                    pattern == 1 ? 0.0f : float((row + col) % 3 != 0);
+            }
+        }
+        ggml_cuda_copy_to_device(mask);
+        ok = compute(output) && ok;
+        ok = check_output(output, reference_weights, x, mask_data, k, rows, columns, mapped) && ok;
+    }
+    // Warm kernels and reusable allocator buffers before measuring retained memory.
+    for (int i = 0; i < 32; ++i) {
+        ok = compute(output) && ok;
+    }
+    CHECK_CUDA(cudaDeviceSynchronize());
+    size_t free_before, free_after, total;
+    CHECK_CUDA(cudaMemGetInfo(&free_before, &total));
+    for (int i = 0; i < iterations; ++i) {
+        ok = compute(output) && ok;
+    }
+    CHECK_CUDA(cudaDeviceSynchronize());
+    CHECK_CUDA(cudaMemGetInfo(&free_after, &total));
+    ok = check_output(output, reference_weights, x, mask_data, k, rows, columns, mapped) && ok;
+    const size_t growth = free_before > free_after ? free_before - free_after : 0;
+    const size_t allowance = 16 * 1024 * 1024;
+    std::printf("type=%s columns=%d mapped=%d iterations=%d free_before=%zu free_after=%zu growth_bytes=%zu\n",
+                argv[1], columns, int(mapped), iterations, free_before, free_after, growth);
+    if (growth > allowance) {
+        std::fprintf(stderr, "CUDA memory grew by %.2f MiB after warmup (limit 16 MiB)\n",
+                     double(growth) / (1024 * 1024));
+        ok = false;
+    }
+    for (ggml_tensor * tensor : { weight, input, mask, output, lookup }) {
+        if (tensor) {
+            ggml_cuda_free_data(tensor);
+        }
+    }
+    ggml_free(ctx);
+    std::printf("%s\n", ok ? "PASS" : "FAIL");
+    return ok ? 0 : 1;
+}