Do W8A16 for FP8 in GEMV case

(base) osimons@ub-osimons:~/lcpp_fp8.git$ ./build-x64-linux-gcc-reldbg/bin/llama-bench -m /mnt/share/gguf/nvidia/Qwen3.6-27B-NVFP4/Qwen3.6-27B-NVFP4_fp8.gguf -p 8196 -ub 512,2048 -d 0,32000 -n 256 --load-mode dio
ggml_cuda_init: found 1 CUDA devices (Total VRAM: 32165 MiB):
  Device 0: NVIDIA RTX PRO 4500 Blackwell, compute capability 12.0, VMM: yes, VRAM: 32165 MiB
| model                          |       size |     params | backend    | ngl | n_ubatch |         lm |            test |                  t/s |
| ------------------------------ | ---------: | ---------: | ---------- | --: | -------: | ---------: | --------------: | -------------------: |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |      512 |        dio |          pp8196 |       2781.53 ± 8.48 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |      512 |        dio |           tg256 |         37.56 ± 0.02 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |      512 |        dio | pp8196 @ d32000 |       1616.31 ± 2.61 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |      512 |        dio |  tg256 @ d32000 |         34.11 ± 0.03 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |     2048 |        dio |          pp8196 |       2542.00 ± 0.51 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |     2048 |        dio |           tg256 |         37.36 ± 0.02 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |     2048 |        dio | pp8196 @ d32000 |       1762.86 ± 1.33 |
| qwen35 27B NVFP4               |  19.56 GiB |    27.32 B | CUDA       |  -1 |     2048 |        dio |  tg256 @ d32000 |         34.09 ± 0.02 |

ggml/src/ggml-cuda/common.cuh

@@ -1107,13 +1107,6 @@ struct ggml_cuda_type_traits<GGML_TYPE_NVFP4> {
     static constexpr int qi = QI_NVFP4;
 };
 
-template<>
-struct ggml_cuda_type_traits<GGML_TYPE_F8_E4M3> {
-    static constexpr int qk = QK8_1;
-    static constexpr int qr = QR8_1;
-    static constexpr int qi = QI8_1;
-};
-
 template<>
 struct ggml_cuda_type_traits<GGML_TYPE_Q2_K> {
     static constexpr int qk = QK_K;

ggml/src/ggml-cuda/ggml-cuda.cu

@@ -1773,7 +1773,8 @@ static bool ggml_cuda_should_fuse_mul_mat_vec_f(const ggml_tensor * tensor) {
     const bool is_mul_mat_id = tensor->op == GGML_OP_MUL_MAT_ID;
 
     bool use_mul_mat_vec_f =
-        (src0->type == GGML_TYPE_F32 || src0->type == GGML_TYPE_F16 || src0->type == GGML_TYPE_BF16) &&
+        (src0->type == GGML_TYPE_F32 || src0->type == GGML_TYPE_F16 || src0->type == GGML_TYPE_BF16 ||
+         src0->type == GGML_TYPE_F8_E4M3) &&
         src1->type == GGML_TYPE_F32 && dst->type == GGML_TYPE_F32;
 
     const int cc      = ggml_cuda_info().devices[ggml_cuda_get_device()].cc;
@@ -1804,7 +1805,6 @@ static bool ggml_cuda_should_fuse_mul_mat_vec_q(const ggml_tensor * tensor) {
     const int cc = ggml_cuda_info().devices[ggml_cuda_get_device()].cc;
     bool use_mul_mat_vec_q = ggml_cuda_should_use_mmvq(src0->type, cc, src1->ne[1]) && !bad_padding_clear &&
                              src1->type == GGML_TYPE_F32 && dst->type == GGML_TYPE_F32;
-    use_mul_mat_vec_q = use_mul_mat_vec_q && (src0->type != GGML_TYPE_F8_E4M3 || src0->ne[0] % QK8_1 == 0);
 
     // fusion is not universally faster on Pascal
     if (cc <= GGML_CUDA_CC_PASCAL) {
@@ -1844,12 +1844,12 @@ static void ggml_cuda_mul_mat(ggml_backend_cuda_context & ctx, const ggml_tensor
     const int warp_size = ggml_cuda_info().devices[ctx.device].warp_size;
 
     if (ggml_cuda_should_use_mmvf(src0->type, cc, src0->ne, src0->nb, ne11)) {
-        // The custom F16 vector kernel can be used over batched cuBLAS GEMM.
+        // The custom vector kernel can be used over batched cuBLAS GEMM.
         // But this is only faster for GPUs without tensor cores or with a thin src0 matrix (particularly KQV in attention)
         ggml_cuda_mul_mat_vec_f(ctx, src0, src1, nullptr, dst);
         return;
     }
-    // A transposed vector can still use MMVQ (i.e. ne01 == 1)
+    // A transposed vector can still use MMVF (i.e. ne01 == 1)
     if (ne01 == 1 && ne11 > MMVF_MAX_BATCH_SIZE && ne2 == 1 && ne3 == 1
             && src0->type == GGML_TYPE_F32
             && ggml_is_contiguous(src0) && ggml_is_contiguous(src1) && ggml_is_contiguous(dst)
@@ -1867,8 +1867,7 @@ static void ggml_cuda_mul_mat(ggml_backend_cuda_context & ctx, const ggml_tensor
         ggml_cuda_mul_mat_f(ctx, src0, src1, nullptr, dst);
         return;
     }
-    if (ggml_cuda_should_use_mmvq(src0->type, cc, ne11) &&
-            (src0->type != GGML_TYPE_F8_E4M3 || ne00 % QK8_1 == 0)) {
+    if (ggml_cuda_should_use_mmvq(src0->type, cc, ne11)) {
         ggml_cuda_mul_mat_vec_q(ctx, src0, src1, nullptr, dst);
         return;
     }
@@ -1896,12 +1895,14 @@ static bool ggml_cuda_mul_mat_id_needs_sync(const ggml_tensor * dst, const int c
     }
 
     if (dst->ne[2] <= MMVQ_MAX_BATCH_SIZE) {
-        if (ggml_is_quantized(src0->type) || src0->type == GGML_TYPE_F8_E4M3) {
-            if (dst->ne[2] <= get_mmvq_mmid_max_batch(src0->type, cc) &&
-                    (src0->type != GGML_TYPE_F8_E4M3 || src0->ne[0] % QK8_1 == 0)) {
+        if (ggml_is_quantized(src0->type)) {
+            if (dst->ne[2] <= get_mmvq_mmid_max_batch(src0->type, cc)) {
                 return false;
             }
-        } else if (GGML_CUDA_CC_IS_AMD(cc)) {
+        } else if (src0->type == GGML_TYPE_F8_E4M3 &&
+                ggml_cuda_should_use_mmvf(src0->type, cc, src0->ne, src0->nb, dst->ne[2])) {
+            return false;
+        } else if (src0->type != GGML_TYPE_F8_E4M3 && GGML_CUDA_CC_IS_AMD(cc)) {
             return false;
         }
     }
@@ -1933,15 +1934,18 @@ static void ggml_cuda_mul_mat_id(ggml_backend_cuda_context & ctx, ggml_tensor *
     if (src1->type == GGML_TYPE_F32 && dst->type == GGML_TYPE_F32) {
         static_assert(MMVQ_MAX_BATCH_SIZE == MMVF_MAX_BATCH_SIZE);
         if (ne2 <= MMVQ_MAX_BATCH_SIZE) {
-            if (ggml_is_quantized(src0->type) || src0->type == GGML_TYPE_F8_E4M3) {
+            if (ggml_is_quantized(src0->type)) {
                 const int mmvq_mmid_max = get_mmvq_mmid_max_batch(src0->type, cc);
-                if (ne2 <= mmvq_mmid_max &&
-                        (src0->type != GGML_TYPE_F8_E4M3 || ne00 % QK8_1 == 0)) {
+                if (ne2 <= mmvq_mmid_max) {
                     ggml_cuda_mul_mat_vec_q(ctx, src0, src1, ids, dst);
                     return;
                 }
+            } else if (src0->type == GGML_TYPE_F8_E4M3 &&
+                    ggml_cuda_should_use_mmvf(src0->type, cc, src0->ne, src0->nb, ne2)) {
+                ggml_cuda_mul_mat_vec_f(ctx, src0, src1, ids, dst);
+                return;
             } else {
-                if (GGML_CUDA_CC_IS_AMD(cc)) {
+                if (src0->type != GGML_TYPE_F8_E4M3 && GGML_CUDA_CC_IS_AMD(cc)) {
                     ggml_cuda_mul_mat_vec_f(ctx, src0, src1, ids, dst);
                     return;
                 }
@@ -3775,6 +3779,13 @@ static int ggml_cuda_try_fuse(ggml_backend_cuda_context * cuda_ctx, ggml_cgraph
                 fusion_data.glu_op     = ggml_get_glu_op(glu);
                 fusion_data.glu_limit  = ggml_get_op_params_f32(glu, 3);
 
+                if (ggml_cuda_should_fuse_mul_mat_vec_f(up_n)) {
+                    ggml_cuda_mul_mat_vec_f(*cuda_ctx, src0, src1, ids, cgraph->nodes[glu_idx], &fusion_data);
+                    fused_mul_mat_vec = true;
+                    fused_node_count  = n_ops;
+                    break;
+                }
+
                 if (ggml_cuda_should_fuse_mul_mat_vec_q(up_n)) {
                     ggml_cuda_mul_mat_vec_q(*cuda_ctx, src0, src1, ids, cgraph->nodes[glu_idx], &fusion_data);
                     fused_mul_mat_vec = true;
@@ -3869,6 +3880,13 @@ static int ggml_cuda_try_fuse(ggml_backend_cuda_context * cuda_ctx, ggml_cgraph
                 fusion_data.glu_op     = ggml_get_glu_op(glu);
                 fusion_data.glu_limit  = ggml_get_op_params_f32(glu, 3);
 
+                if (ggml_cuda_should_fuse_mul_mat_vec_f(up_n)) {
+                    ggml_cuda_mul_mat_vec_f(*cuda_ctx, src0, src1, ids, cgraph->nodes[glu_idx], &fusion_data);
+                    fused_mul_mat_vec = true;
+                    fused_node_count  = n_ops;
+                    break;
+                }
+
                 if (ggml_cuda_should_fuse_mul_mat_vec_q(up_n)) {
                     ggml_cuda_mul_mat_vec_q(*cuda_ctx, src0, src1, ids, cgraph->nodes[glu_idx], &fusion_data);
                     fused_mul_mat_vec = true;
@@ -4066,6 +4084,13 @@ static int ggml_cuda_try_fuse(ggml_backend_cuda_context * cuda_ctx, ggml_cgraph
             fusion_data.x_bias  = bias;
             fusion_data.x_scale = scale;
 
+            if (ggml_cuda_should_fuse_mul_mat_vec_f(mm_node)) {
+                ggml_cuda_mul_mat_vec_f(*cuda_ctx, src0, src1, ids, out_node, &fusion_data);
+                fused_mul_mat_vec = true;
+                fused_node_count  = n_ops;
+                break;
+            }
+
             if (ggml_cuda_should_fuse_mul_mat_vec_q(mm_node)) {
                 ggml_cuda_mul_mat_vec_q(*cuda_ctx, src0, src1, ids, out_node, &fusion_data);
                 fused_mul_mat_vec = true;

ggml/src/ggml-cuda/mmvf.cu

@@ -4,15 +4,30 @@
 #include "mmvf.cuh"
 #include "convert.cuh"
 
+template <typename T>
+using mmvf_y_t = std::conditional_t<std::is_same_v<T, ggml_fp8_e4m3_t>, nv_bfloat16, float>;
+
+#if !defined(GGML_USE_HIP) && !defined(GGML_USE_MUSA)
+static __device__ __forceinline__ nv_bfloat16 mmvf_f8_e4m3_to_bf16(uint8_t bits) {
+#if defined(FP8_AVAILABLE)
+    __nv_fp8_e4m3 value;
+    value.__x = bits;
+    return static_cast<nv_bfloat16>(value);
+#else
+    return static_cast<nv_bfloat16>(ggml_cuda_f8_e4m3_to_fp32(bits));
+#endif
+}
+#endif
+
 template <typename T, typename type_acc, int ncols_dst, int block_size, bool has_fusion = false, bool is_multi_token_id = false>
 static __global__ void mul_mat_vec_f(
-        const T * x_ptr, const float * y_ptr, const int32_t * ids_ptr, const ggml_cuda_mm_fusion_args_device fusion, float * dst_ptr,
+        const T * x_ptr, const mmvf_y_t<T> * y_ptr, const int32_t * ids_ptr, const ggml_cuda_mm_fusion_args_device fusion, float * dst_ptr,
         const int ncols2, const uint3 nchannels_y, const int stride_row, const int stride_col_y2, const int stride_col_dst,
         const uint3 channel_ratio, const int stride_channel_x, const int stride_channel_y, const int stride_channel_dst,
         const uint3 sample_ratio, const int stride_sample_x, const int stride_sample_y, const int stride_sample_dst,
         const int ids_stride) {
     const T       * GGML_CUDA_RESTRICT x   = x_ptr;
-    const float   * GGML_CUDA_RESTRICT y   = y_ptr;
+    const mmvf_y_t<T> * GGML_CUDA_RESTRICT y = y_ptr;
     const int32_t * GGML_CUDA_RESTRICT ids = ids_ptr;
     float         * GGML_CUDA_RESTRICT dst = dst_ptr;
     const int row         = blockIdx.x;
@@ -55,16 +70,24 @@ static __global__ void mul_mat_vec_f(
     bool use_gate = false;
     bool use_bias = false;
     bool use_gate_bias = false;
+    bool use_scale = false;
+    bool use_gate_scale = false;
     ggml_glu_op glu_op = ggml_glu_op::GGML_GLU_OP_SWIGLU;
     float glu_limit = 0.0f;
     const T * gate_x = nullptr;
     const float * x_bias = nullptr;
     const float * gate_bias = nullptr;
+    const float * x_scale = nullptr;
+    const float * gate_scale = nullptr;
 
     if constexpr (has_fusion) {
         use_gate = fusion.gate != nullptr;
         use_bias = fusion.x_bias != nullptr;
         use_gate_bias = fusion.gate_bias != nullptr;
+        if constexpr (std::is_same_v<T, ggml_fp8_e4m3_t>) {
+            use_scale = fusion.x_scale != nullptr;
+            use_gate_scale = fusion.gate_scale != nullptr && use_gate;
+        }
         glu_op = fusion.glu_op;
         glu_limit = fusion.glu_limit;
 
@@ -80,6 +103,12 @@ static __global__ void mul_mat_vec_f(
         } else {
             use_gate_bias = false;
         }
+        if (use_scale) {
+            x_scale = static_cast<const float *>(fusion.x_scale);
+        }
+        if (use_gate_scale) {
+            gate_scale = static_cast<const float *>(fusion.gate_scale);
+        }
     }
 
     if (use_gate) {
@@ -96,8 +125,6 @@ static __global__ void mul_mat_vec_f(
         }
     }
 
-    const float2 * y2 = (const float2 *) y;
-
     extern __shared__ char data_mmv[];
     float * buf_iw = (float *) data_mmv;
     [[maybe_unused]] float * buf_iw_gate = nullptr;
@@ -127,6 +154,7 @@ static __global__ void mul_mat_vec_f(
     }
 
     if constexpr (std::is_same_v<T, float>) {
+        const float2 * y2 = (const float2 *) y;
         const float2 * x2 = (const float2 *) x;
         [[maybe_unused]] const float2 * gate_x2 = nullptr;
         if constexpr (has_fusion) {
@@ -159,6 +187,7 @@ static __global__ void mul_mat_vec_f(
             }
         }
     } else if constexpr (std::is_same_v<T, half>) {
+        const float2 * y2 = (const float2 *) y;
         const half2 * x2 = (const half2 *) x;
         [[maybe_unused]] const half2 * gate_x2 = nullptr;
         if constexpr (has_fusion) {
@@ -234,6 +263,7 @@ static __global__ void mul_mat_vec_f(
 #endif // FP16_AVAILABLE
         }
     } else if constexpr (std::is_same_v<T, nv_bfloat16>) {
+        const float2 * y2 = (const float2 *) y;
 //TODO: add support for ggml_cuda_mad for hip_bfloat162
 #if defined(GGML_USE_HIP)
         const int * x2 = (const int *) x;
@@ -299,6 +329,73 @@ static __global__ void mul_mat_vec_f(
                 }
             }
         }
+#endif
+    } else if constexpr (std::is_same_v<T, ggml_fp8_e4m3_t>) {
+        const nv_bfloat162 * y2 = (const nv_bfloat162 *) y;
+#if !defined(GGML_USE_HIP) && !defined(GGML_USE_MUSA) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 800
+        nv_bfloat162 sum_bf[ncols_dst] = {};
+        nv_bfloat162 sum_bf_gate[ncols_dst] = {};
+
+        for (int col2 = tid; col2 < ncols2; col2 += block_size) {
+            const nv_bfloat162 tmpx = make_bfloat162(
+                mmvf_f8_e4m3_to_bf16(x[2*col2 + 0].bits),
+                mmvf_f8_e4m3_to_bf16(x[2*col2 + 1].bits));
+            nv_bfloat162 tmpx_gate = {};
+            if constexpr (has_fusion) {
+                if (use_gate) {
+                    tmpx_gate = make_bfloat162(
+                        mmvf_f8_e4m3_to_bf16(gate_x[2*col2 + 0].bits),
+                        mmvf_f8_e4m3_to_bf16(gate_x[2*col2 + 1].bits));
+                }
+            }
+#pragma unroll
+            for (int j = 0; j < ncols_dst; ++j) {
+                const nv_bfloat162 tmpy = y2[j*stride_col_y2 + col2];
+                sum_bf[j] = __hfma2(tmpx, tmpy, sum_bf[j]);
+
+                if constexpr (has_fusion) {
+                    if (use_gate) {
+                        sum_bf_gate[j] = __hfma2(tmpx_gate, tmpy, sum_bf_gate[j]);
+                    }
+                }
+            }
+        }
+
+#pragma unroll
+        for (int j = 0; j < ncols_dst; ++j) {
+            sumf[j] = __low2float(sum_bf[j]) + __high2float(sum_bf[j]);
+            if constexpr (has_fusion) {
+                if (use_gate) {
+                    sumf_gate[j] = __low2float(sum_bf_gate[j]) + __high2float(sum_bf_gate[j]);
+                }
+            }
+        }
+#else
+        for (int col2 = tid; col2 < ncols2; col2 += block_size) {
+            const float x0 = ggml_cuda_f8_e4m3_to_fp32(x[2*col2 + 0].bits);
+            const float x1 = ggml_cuda_f8_e4m3_to_fp32(x[2*col2 + 1].bits);
+            float gate_x0 = 0.0f;
+            float gate_x1 = 0.0f;
+            if constexpr (has_fusion) {
+                if (use_gate) {
+                    gate_x0 = ggml_cuda_f8_e4m3_to_fp32(gate_x[2*col2 + 0].bits);
+                    gate_x1 = ggml_cuda_f8_e4m3_to_fp32(gate_x[2*col2 + 1].bits);
+                }
+            }
+#pragma unroll
+            for (int j = 0; j < ncols_dst; ++j) {
+                const float2 tmpy = ggml_cuda_cast<float2>(y2[j*stride_col_y2 + col2]);
+                ggml_cuda_mad(sumf[j], x0, tmpy.x);
+                ggml_cuda_mad(sumf[j], x1, tmpy.y);
+
+                if constexpr (has_fusion) {
+                    if (use_gate) {
+                        ggml_cuda_mad(sumf_gate[j], gate_x0, tmpy.x);
+                        ggml_cuda_mad(sumf_gate[j], gate_x1, tmpy.y);
+                    }
+                }
+            }
+        }
 #endif
     } else {
         static_assert(std::is_same_v<T, void>, "unsupported type");
@@ -347,12 +444,18 @@ static __global__ void mul_mat_vec_f(
     float value = sumf[tid];
 
     if constexpr (has_fusion) {
+        if (use_scale) {
+            value *= x_scale[ids ? channel_x : 0];
+        }
         if (use_bias) {
             value += x_bias[tid*stride_col_dst + row];
         }
 
         if (use_gate) {
             float gate_value = sumf_gate[tid];
+            if (use_gate_scale) {
+                gate_value *= gate_scale[ids ? channel_x : 0];
+            }
             if (use_gate_bias) {
                 gate_value += gate_bias[tid*stride_col_dst + row];
             }
@@ -379,13 +482,17 @@ static __global__ void mul_mat_vec_f(
     dst[tid*stride_col_dst + row] = value;
 
     if constexpr (!has_fusion) {
-        GGML_UNUSED_VARS(use_gate, use_bias, use_gate_bias, glu_op, glu_limit, gate_x, x_bias, gate_bias, sumf_gate);
+        GGML_UNUSED_VARS(use_gate, use_bias, use_gate_bias, use_scale, use_gate_scale, glu_op, glu_limit, gate_x, x_bias, gate_bias,
+            x_scale, gate_scale, sumf_gate);
+    }
+    if constexpr (!std::is_same_v<T, ggml_fp8_e4m3_t>) {
+        GGML_UNUSED_VARS(use_scale, use_gate_scale, x_scale, gate_scale);
     }
 }
 
 template<typename T, typename type_acc, int ncols_dst, int block_size, bool is_multi_token_id = false>
 static void mul_mat_vec_f_switch_fusion(
-        const T * x, const float * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
+        const T * x, const mmvf_y_t<T> * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
         const int64_t ncols, const uint3 nchannels_y,
         const int64_t stride_row, const int64_t stride_col_y, const int64_t stride_col_dst,
         const uint3 channel_ratio, const int stride_channel_x, const int stride_channel_y, const int stride_channel_dst,
@@ -394,7 +501,8 @@ static void mul_mat_vec_f_switch_fusion(
 
     const ggml_cuda_kernel_launch_params launch_params = {block_nums, block_dims, nbytes_shared, stream};
 
-    const bool has_fusion = fusion.gate != nullptr || fusion.x_bias != nullptr || fusion.gate_bias != nullptr;
+    const bool has_fusion = fusion.gate != nullptr || fusion.x_bias != nullptr || fusion.gate_bias != nullptr ||
+        fusion.x_scale != nullptr || fusion.gate_scale != nullptr;
     if constexpr (ncols_dst == 1) {
         if (has_fusion) {
             ggml_cuda_kernel_launch(mul_mat_vec_f<T, type_acc, ncols_dst, block_size, true, is_multi_token_id>, launch_params,
@@ -416,7 +524,7 @@ static void mul_mat_vec_f_switch_fusion(
 
 template <typename T, typename type_acc, int ncols_dst, bool is_multi_token_id = false>
 void launch_mul_mat_vec_f_cuda(
-        const T * x, const float * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
+        const T * x, const mmvf_y_t<T> * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
         const int64_t ncols, const int64_t nrows,
         const int64_t stride_row, const int64_t stride_col_y, const int64_t stride_col_dst,
         const int64_t nchannels_x, const int64_t nchannels_y, const int64_t nchannels_dst,
@@ -511,7 +619,7 @@ void launch_mul_mat_vec_f_cuda(
 
 template <typename T, typename type_acc>
 static void mul_mat_vec_f_cuda_switch_ncols_dst(
-        const T * x, const float * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
+        const T * x, const mmvf_y_t<T> * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
         const int64_t ncols, const int64_t nrows, const int64_t ncols_dst,
         const int64_t stride_row, const int64_t stride_col_y, const int64_t stride_col_dst,
         const int64_t nchannels_x, const int64_t nchannels_y, const int64_t nchannels_dst,
@@ -608,7 +716,7 @@ static void mul_mat_vec_f_cuda_switch_ncols_dst(
 
 template<typename T>
 static void mul_mat_vec_f_cuda(
-        const T * x, const float * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
+        const T * x, const mmvf_y_t<T> * y, const int32_t * ids, const ggml_cuda_mm_fusion_args_device fusion, float * dst,
         const int64_t ncols, const int64_t nrows, const int64_t ncols_dst,
         const int64_t stride_row, const int64_t stride_col_y, const int stride_col_dst,
         const int64_t nchannels_x, const int64_t nchannels_y, const int64_t nchannels_dst,
@@ -663,6 +771,7 @@ void ggml_cuda_mul_mat_vec_f(ggml_backend_cuda_context & ctx, const ggml_tensor
     if (fusion) {
         GGML_ASSERT( !ids || dst->ne[2] == 1);
         GGML_ASSERT(  ids || dst->ne[1] == 1);
+        GGML_ASSERT((fusion->x_scale == nullptr && fusion->gate_scale == nullptr) || src0->type == GGML_TYPE_F8_E4M3);
         if (fusion->x_bias) {
             GGML_ASSERT(fusion->x_bias->type == GGML_TYPE_F32);
             GGML_ASSERT(fusion->x_bias->ne[0] == dst->ne[0]);
@@ -679,6 +788,18 @@ void ggml_cuda_mul_mat_vec_f(ggml_backend_cuda_context & ctx, const ggml_tensor
             GGML_ASSERT(!ids || fusion->gate_bias->ne[1] == src0->ne[2]);
             fusion_local.gate_bias = fusion->gate_bias->data;
         }
+        if (fusion->x_scale) {
+            GGML_ASSERT(fusion->x_scale->type == GGML_TYPE_F32);
+            GGML_ASSERT(ggml_is_contiguous(fusion->x_scale));
+            GGML_ASSERT(ggml_nelements(fusion->x_scale) == (ids ? src0->ne[2] : 1));
+            fusion_local.x_scale = fusion->x_scale->data;
+        }
+        if (fusion->gate_scale) {
+            GGML_ASSERT(fusion->gate_scale->type == GGML_TYPE_F32);
+            GGML_ASSERT(ggml_is_contiguous(fusion->gate_scale));
+            GGML_ASSERT(ggml_nelements(fusion->gate_scale) == (ids ? src0->ne[2] : 1));
+            fusion_local.gate_scale = fusion->gate_scale->data;
+        }
         fusion_local.glu_op = fusion->glu_op;
         fusion_local.glu_limit = fusion->glu_limit;
     }
@@ -723,6 +844,22 @@ void ggml_cuda_mul_mat_vec_f(ggml_backend_cuda_context & ctx, const ggml_tensor
                 ne02, nchannels_y, nchannels_dst, s02, stride_channel_y, stride_channel_dst,
                 ne03,              ne3,           s03, s13,              s3,                 ids_stride, prec, ctx.stream());
         } break;
+        case GGML_TYPE_F8_E4M3: {
+            ggml_cuda_pool_alloc<nv_bfloat16> src1_bf16(ctx.pool(), ggml_nelements(src1));
+            to_bf16_nc_cuda_t to_bf16 = ggml_get_to_bf16_nc_cuda(GGML_TYPE_F32);
+            GGML_ASSERT(to_bf16 != nullptr);
+            to_bf16(src1->data, src1_bf16.get(), ne10, ne11, ne12, ne13, s11, s12, s13, ctx.stream());
+
+            const int64_t bs11 = ne10;
+            const int64_t bs12 = ne11*bs11;
+            const int64_t bs13 = ne12*bs12;
+            const int64_t stride_col_y_bf16 = ids ? bs12 : bs11;
+            const int64_t stride_channel_y_bf16 = ids ? bs11 : bs12;
+            const ggml_fp8_e4m3_t * src0_d = (const ggml_fp8_e4m3_t *) src0->data;
+            mul_mat_vec_f_cuda(src0_d, src1_bf16.get(), ids_d, fusion_local, dst_d, ne00, ne01, ncols_dst, s01,
+                stride_col_y_bf16, stride_col_dst, ne02, nchannels_y, nchannels_dst, s02, stride_channel_y_bf16,
+                stride_channel_dst, ne03, ne3, s03, bs13, s3, ids_stride, prec, ctx.stream());
+        } break;
         default:
             GGML_ABORT("unsupported type: %s", ggml_type_name(src0->type));
     }
@@ -869,6 +1006,8 @@ bool ggml_cuda_should_use_mmvf(enum ggml_type type, int cc, const int64_t * src0
                 return ne11 <= 8;
             }
             return ne11 <= 8;
+        case GGML_TYPE_F8_E4M3:
+            return ne11 <= MMVF_MAX_BATCH_SIZE;
         default:
             return false;
     }

ggml/src/ggml-cuda/mmvq.cu

@@ -9,7 +9,7 @@
 typedef float (*vec_dot_q_cuda_t)(const void * __restrict__ vbq, const block_q8_1 * __restrict__ bq8_1, const int & kbx, const int & iqs);
 
 static constexpr __host__ __device__ bool is_scaled_low_precision_type(ggml_type type) {
-    return type == GGML_TYPE_NVFP4 || type == GGML_TYPE_F8_E4M3;
+    return type == GGML_TYPE_NVFP4;
 }
 
 static constexpr __device__ vec_dot_q_cuda_t get_vec_dot_q_cuda(ggml_type type) {
@@ -23,7 +23,6 @@ static constexpr __device__ vec_dot_q_cuda_t get_vec_dot_q_cuda(ggml_type type)
         case GGML_TYPE_Q8_0:    return vec_dot_q8_0_q8_1;
         case GGML_TYPE_MXFP4:   return vec_dot_mxfp4_q8_1;
         case GGML_TYPE_NVFP4:   return vec_dot_nvfp4_q8_1;
-        case GGML_TYPE_F8_E4M3: return vec_dot_f8_e4m3_q8_1;
         case GGML_TYPE_Q2_K:    return vec_dot_q2_K_q8_1;
         case GGML_TYPE_Q3_K:    return vec_dot_q3_K_q8_1;
         case GGML_TYPE_Q4_K:    return vec_dot_q4_K_q8_1;
@@ -53,7 +52,6 @@ static constexpr __host__ __device__ int get_vdr_mmvq(ggml_type type) {
         case GGML_TYPE_Q8_0:    return VDR_Q8_0_Q8_1_MMVQ;
         case GGML_TYPE_MXFP4:   return VDR_MXFP4_Q8_1_MMVQ;
         case GGML_TYPE_NVFP4:   return VDR_NVFP4_Q8_1_MMVQ;
-        case GGML_TYPE_F8_E4M3: return VDR_F8_E4M3_Q8_1_MMVQ;
         case GGML_TYPE_Q2_K:    return VDR_Q2_K_Q8_1_MMVQ;
         case GGML_TYPE_Q3_K:    return VDR_Q3_K_Q8_1_MMVQ;
         case GGML_TYPE_Q4_K:    return VDR_Q4_K_Q8_1_MMVQ;
@@ -293,7 +291,7 @@ int get_mmvq_mmid_max_batch(ggml_type type, int cc) {
 }
 
 bool ggml_cuda_should_use_mmvq(enum ggml_type type, int cc, int64_t ne11) {
-    if (!ggml_is_quantized(type) && type != GGML_TYPE_F8_E4M3) {
+    if (!ggml_is_quantized(type)) {
         return false;
     }
     // k-quants cost more to decode and mvq redoes that per column, so MMQ wins sooner.
@@ -576,7 +574,6 @@ static __global__ void mul_mat_vec_q(
     constexpr int qk  = ggml_cuda_type_traits<type>::qk;
     constexpr int qi  = ggml_cuda_type_traits<type>::qi;
     constexpr int vdr = get_vdr_mmvq(type);
-    constexpr int kbx_stride = type == GGML_TYPE_F8_E4M3 ? QK8_1 : 1;
     constexpr mmvq_parameter_table_id table_id = get_device_table_id();
     constexpr int nwarps = calc_nwarps(type, ncols_dst, table_id, small_k, halve_iters);
     constexpr int rows_per_cuda_block = calc_rows_per_block(ncols_dst, table_id, small_k, nwarps);
@@ -687,11 +684,11 @@ static __global__ void mul_mat_vec_q(
 #pragma unroll
             for (int i = 0; i < rows_per_cuda_block; ++i) {
                 tmp[j][i] += vec_dot_q_cuda(
-                    vx, &y[j*stride_col_y + kby], kbx_offset + i*stride_row_x + kbx_stride*kbx, kqs);
+                    vx, &y[j*stride_col_y + kby], kbx_offset + i*stride_row_x + kbx, kqs);
                 if constexpr (has_fusion) {
                     if (use_gate) {
                         tmp_gate[j][i] += vec_dot_q_cuda(
-                            vgate, &y[j*stride_col_y + kby], kbx_offset + i*stride_row_x + kbx_stride*kbx, kqs);
+                            vgate, &y[j*stride_col_y + kby], kbx_offset + i*stride_row_x + kbx, kqs);
                     }
                 }
             }
@@ -809,7 +806,6 @@ static __global__ void mul_mat_vec_q_moe(
     constexpr int qk  = ggml_cuda_type_traits<type>::qk;
     constexpr int qi  = ggml_cuda_type_traits<type>::qi;
     constexpr int vdr = get_vdr_mmvq(type);
-    constexpr int kbx_stride = type == GGML_TYPE_F8_E4M3 ? QK8_1 : 1;
     constexpr int warp_size = ggml_cuda_get_physical_warp_size();
 
     constexpr vec_dot_q_cuda_t vec_dot_q_cuda = get_vec_dot_q_cuda(type);
@@ -865,10 +861,10 @@ static __global__ void mul_mat_vec_q_moe(
 
 #pragma unroll
         for (int i = 0; i < c_rows_per_block; ++i) {
-            tmp[i] += vec_dot_q_cuda(vx, &y[kby], kbx_offset + i*stride_row_x + kbx_stride*kbx, kqs);
+            tmp[i] += vec_dot_q_cuda(vx, &y[kby], kbx_offset + i*stride_row_x + kbx, kqs);
             if constexpr (has_fusion) {
                 if (use_gate) {
-                    tmp_gate[i] += vec_dot_q_cuda(vgate, &y[kby], kbx_offset + i*stride_row_x + kbx_stride*kbx, kqs);
+                    tmp_gate[i] += vec_dot_q_cuda(vgate, &y[kby], kbx_offset + i*stride_row_x + kbx, kqs);
                 }
             }
         }
@@ -1273,12 +1269,6 @@ static void mul_mat_vec_q_switch_type(
                  nchannels_x, nchannels_y, nchannels_dst, stride_channel_x, stride_channel_y, stride_channel_dst,
                  nsamples_x, nsamples_dst, stride_sample_x, stride_sample_y, stride_sample_dst, ids_stride, stream);
             break;
-        case GGML_TYPE_F8_E4M3:
-            mul_mat_vec_q_switch_ncols_dst<GGML_TYPE_F8_E4M3>
-                (vx, vy, ids, fusion, dst, ncols_x, nrows_x, ncols_dst, stride_row_x, stride_col_y, stride_col_dst,
-                 nchannels_x, nchannels_y, nchannels_dst, stride_channel_x, stride_channel_y, stride_channel_dst,
-                 nsamples_x, nsamples_dst, stride_sample_x, stride_sample_y, stride_sample_dst, ids_stride, stream);
-            break;
         case GGML_TYPE_Q2_K:
             mul_mat_vec_q_switch_ncols_dst<GGML_TYPE_Q2_K>
                 (vx, vy, ids, fusion, dst, ncols_x, nrows_x, ncols_dst, stride_row_x, stride_col_y, stride_col_dst,

ggml/src/ggml-cuda/vecdotq.cuh

@@ -361,24 +361,6 @@ static __device__ __forceinline__ float vec_dot_nvfp4_q8_1(
     return sum;
 }
 
-#define VDR_F8_E4M3_Q8_1_MMVQ 1
-
-static __device__ __forceinline__ float vec_dot_f8_e4m3_q8_1(
-        const void * __restrict__ vbq,
-        const block_q8_1 * __restrict__ bq8_1,
-        const int32_t & kbx,
-        const int32_t & iqs) {
-    const ggml_fp8_e4m3_t * bq8 = (const ggml_fp8_e4m3_t *) vbq + kbx;
-    const int8_t * q8 = bq8_1->qs + 4*iqs;
-
-    float sum = 0.0f;
-#pragma unroll
-    for (int i = 0; i < 4; ++i) {
-        sum += ggml_cuda_f8_e4m3_to_fp32(bq8[4*iqs + i].bits) * q8[i];
-    }
-    return __low2float(bq8_1->ds) * sum;
-}
-
 #define VDR_Q2_K_Q8_1_MMVQ 1
 #define VDR_Q2_K_Q8_1_MMQ  4