Implement MatmulInteger on GPU (#3070)

* Implement MatmulInteger

onnxruntime/core/providers/cuda/cuda_execution_provider.cc

@@ -258,6 +258,7 @@ class ONNX_OPERATOR_VERSIONED_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kO
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, float, MatMul);
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, double, MatMul);
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, MLFloat16, MatMul);
+class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 10, int8_t, MatMulInteger);
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, float, Tile);
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, double, Tile);
 class ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, MLFloat16, Tile);
@@ -724,6 +725,7 @@ static void RegisterCudaKernels(KernelRegistry& kernel_registry) {
       BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, float, MatMul)>,
       BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, double, MatMul)>,
       BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 9, MLFloat16, MatMul)>,
+      BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 10, int8_t, MatMulInteger)>,
       BuildKernelCreateInfo<ONNX_OPERATOR_VERSIONED_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, 10, float, Clip)>,
       BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, float, Tile)>,
       BuildKernelCreateInfo<ONNX_OPERATOR_TYPED_KERNEL_CLASS_NAME(kCudaExecutionProvider, kOnnxDomain, 6, double, Tile)>,

onnxruntime/core/providers/cuda/math/matmul_integer.cc

@@ -0,0 +1,156 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+#include "matmul_integer.h"
+#include "matmul_integer.cuh"
+#include "core/providers/cpu/math/matmul_helper.h"
+#include "core/providers/cuda/shared_inc/fpgeneric.h"
+#include "core/providers/cuda/cuda_allocator.h"
+#include "core/providers/common.h"
+
+namespace onnxruntime {
+namespace cuda {
+
+ONNX_OPERATOR_TYPED_KERNEL_EX(
+    MatMulInteger,
+    kOnnxDomain,
+    10,
+    int8_t,
+    kCudaExecutionProvider,
+    KernelDefBuilder()
+        .TypeConstraint("T1", DataTypeImpl::GetTensorType<int8_t>())
+        .TypeConstraint("T2", DataTypeImpl::GetTensorType<int8_t>())
+        .TypeConstraint("T3", DataTypeImpl::GetTensorType<int32_t>()),
+    MatMulInteger<int8_t, int8_t>);
+
+template <>
+Status MatMulInteger<int8_t, int8_t>::PadMatrix(
+    int64_t row,
+    int64_t col,
+    int64_t align_size,
+    const int8_t*& src,
+    int64_t& pad_size,
+    IAllocatorUniquePtr<int8_t>& temp_mem_holder) const {
+  pad_size = align_size - col % align_size;
+  if (pad_size != align_size) {
+    temp_mem_holder = GetScratchBuffer<int8_t>(row * (col + pad_size));
+    ORT_RETURN_IF_ERROR(PadMatrixInLeadingDimension(src, temp_mem_holder.get(), row, col, pad_size));
+    src = temp_mem_holder.get();
+  } else {
+    pad_size = 0;
+  }
+
+  return Status::OK();
+}
+
+template <>
+Status MatMulInteger<int8_t, int8_t>::ComputeInternal(OpKernelContext* ctx) const {
+  auto a = ctx->Input<Tensor>(0);
+  auto b = ctx->Input<Tensor>(1);
+  ORT_ENFORCE(a != nullptr && b != nullptr);
+
+  MatMulComputeHelper helper;
+  ORT_RETURN_IF_ERROR(helper.Compute(a->Shape(), b->Shape()));
+  Tensor* Y = ctx->Output(0, helper.OutputShape());
+
+  const int8_t* a_ptr = a->template Data<int8_t>();
+  const int8_t* b_ptr = b->template Data<int8_t>();
+  int32_t* output_ptr = Y->template MutableData<int32_t>();
+
+  // validate zero points
+  const int8_t* a_offset = nullptr;
+  const int8_t* b_offset = nullptr;
+  if (has_a_zero_point_) {
+    auto a_zero_point = ctx->Input<Tensor>(2);
+    ORT_ENFORCE(IsScalarOr1ElementVector(a_zero_point),
+                "MatmulInteger : input1 zero point must be a scalar or 1D tensor of size 1");
+    a_offset = a_zero_point->template Data<int8_t>();
+  }
+  if (has_b_zero_point_) {
+    auto b_zero_point = ctx->Input<Tensor>(3);
+    ORT_ENFORCE(IsScalarOr1ElementVector(b_zero_point),
+                "MatmulInteger : input2 zero point must be a scalar or 1D tensor of size 1");
+    b_offset = b_zero_point->template Data<int8_t>();
+  }
+
+  // offset output c[i,j] to
+  // k*a_offset*b_offset -
+  // b_offset * (a[i,0] + a[i,1] ...+a[i,k]) -
+  // a_offset * (b[0,j] + b[1,j] ... + b[k,j])
+  // ReduceRowSumOnMatrixA computes the b_offset * (a[i,0] + a[i,1] ...+a[i,k]) part
+  // ReduceColSumOnMatrixB computes the a_offset * (b[0,j] + b[1,j] ... + b[k,j]) part
+  // OffsetOutput computes gets the final result
+  IAllocatorUniquePtr<int32_t> a_row_buf;
+  if (has_b_zero_point_) {
+    a_row_buf = GetScratchBuffer<int32_t>(helper.OutputShape().Size() / helper.N());
+    ORT_RETURN_IF_ERROR(ReduceRowSumOnMatrixA(a_ptr, a_row_buf.get(), b_offset, helper));
+  }
+
+  IAllocatorUniquePtr<int32_t> b_col_buf;
+  if (has_a_zero_point_) {
+    b_col_buf = GetScratchBuffer<int32_t>(helper.OutputShape().Size() / helper.M());
+    ORT_RETURN_IF_ERROR(ReduceColSumOnMatrixB(b_ptr, b_col_buf.get(), a_offset, helper));
+  }
+
+  int alpha = 1;
+  int beta = 0;
+  if (has_a_zero_point_ || has_b_zero_point_) {
+    OffsetOutput(a_row_buf.get(),
+                 b_col_buf.get(),
+                 output_ptr,
+                 a_offset,
+                 b_offset,
+                 helper);
+    beta = 1;
+  }
+
+  // pad A and B to make their leading dimension be multiples of 32
+  // because cublasGemmEx requires:
+  // 1. leading dimension is multiples of 4
+  // 2. A, B is 32-bit aligned
+  const int64_t align_size = 32;
+  int64_t a_pad_size = 0;
+  int64_t b_pad_size = 0;
+  IAllocatorUniquePtr<int8_t> a_padded;
+  IAllocatorUniquePtr<int8_t> b_padded;
+  ORT_RETURN_IF_ERROR(PadMatrix(a->Shape().Size() / helper.K(),
+                                helper.K(),
+                                align_size,
+                                a_ptr,
+                                a_pad_size,
+                                a_padded));
+  ORT_RETURN_IF_ERROR(PadMatrix(b->Shape().Size() / helper.N(),
+                                helper.N(),
+                                align_size,
+                                b_ptr,
+                                b_pad_size,
+                                b_padded));
+
+  for (int batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+    CUBLAS_RETURN_IF_ERROR(cublasGemmEx(
+        Base::CublasHandle(),
+        CUBLAS_OP_N,
+        CUBLAS_OP_N,
+        static_cast<int>(helper.N()),
+        static_cast<int>(helper.M()),
+        static_cast<int>(helper.K()),
+        &alpha,
+        b_ptr + helper.RightOffsets()[batch] + helper.RightOffsets()[batch] / helper.N() * b_pad_size,
+        CUDA_R_8I,
+        static_cast<int>(helper.N() + b_pad_size),
+        a_ptr + helper.LeftOffsets()[batch] + helper.LeftOffsets()[batch] / helper.K() * a_pad_size,
+        CUDA_R_8I,
+        static_cast<int>(helper.K() + a_pad_size),
+        &beta,
+        output_ptr + helper.OutputOffsets()[batch],
+        CUDA_R_32I,
+        static_cast<int>(helper.N()),
+        CUDA_R_32I,
+        CUBLAS_GEMM_DFALT));
+  }
+
+  return Status::OK();
+}
+
+}  // namespace cuda
+}  // namespace onnxruntime

onnxruntime/core/providers/cuda/math/matmul_integer.cu

@@ -0,0 +1,152 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+#include "matmul_integer.cuh"
+
+#include <cub/cub.cuh>
+#include "core/providers/cuda/cu_inc/common.cuh"
+
+namespace onnxruntime {
+namespace cuda {
+
+template <int TPB>
+__global__ void ReduceRowSumOnMatrixAKernel(const int8_t* matrix, int32_t* row_sum, const int8_t* offset, int32_t K) {
+  int32_t thread_data = 0;
+  const int8_t* row_ptr = matrix + blockIdx.x * K;
+  for (int i = threadIdx.x; i < K; i += TPB) {
+    thread_data += *(row_ptr + i);
+  }
+
+  using BlockReduce = cub::BlockReduce<int32_t, TPB>;
+  __shared__ typename BlockReduce::TempStorage temp_storage;
+  int32_t sum = BlockReduce(temp_storage).Sum(thread_data);
+
+  if (threadIdx.x == 0) {
+    row_sum[blockIdx.x] = (*offset) * sum;
+  }
+}
+
+Status ReduceRowSumOnMatrixA(const int8_t* matrix, int32_t* row_sum, const int8_t* offset, const MatMulComputeHelper& helper) {
+  for (size_t batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+    ReduceRowSumOnMatrixAKernel<static_cast<int>(GridDim::maxThreadsPerBlock)>
+        <<<static_cast<int>(helper.M()), GridDim::maxThreadsPerBlock, 0>>>(matrix + helper.LeftOffsets()[batch],
+                                                                           row_sum + batch * helper.M(),
+                                                                           offset,
+                                                                           static_cast<int>(helper.K()));
+  }
+
+  return CUDA_CALL( cudaPeekAtLastError() ) ? Status::OK() : Status( common::ONNXRUNTIME, common::FAIL );
+}
+
+template <int TPB>
+__global__ void ReduceColSumOnMatrixBKernel(const int8_t* matrix, int32_t* col_sum, const int8_t* offset, int32_t row, int32_t col) {
+  int32_t thread_data = 0;
+  const int8_t* col_ptr = matrix + blockIdx.x;
+  for (int i = threadIdx.x; i < row; i += TPB) {
+    thread_data += *(col_ptr + i * col);
+  }
+
+  using BlockReduce = cub::BlockReduce<int32_t, TPB>;
+  __shared__ typename BlockReduce::TempStorage temp_storage;
+  int32_t sum = BlockReduce(temp_storage).Sum(thread_data);
+
+  if (threadIdx.x == 0) {
+    col_sum[blockIdx.x] = (*offset) * sum;
+  }
+}
+
+Status ReduceColSumOnMatrixB(const int8_t* matrix, int32_t* col_sum, const int8_t* offset, const MatMulComputeHelper& helper) {
+  for (size_t batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+    ReduceColSumOnMatrixBKernel<static_cast<int>(GridDim::maxThreadsPerBlock)>
+        <<<static_cast<int>(helper.N()), GridDim::maxThreadsPerBlock, 0>>>(matrix + helper.RightOffsets()[batch],
+                                                                           col_sum + batch * helper.N(),
+                                                                           offset,
+                                                                           static_cast<int32_t>(helper.K()),
+                                                                           static_cast<int32_t>(helper.N()));
+  }
+
+  return CUDA_CALL( cudaPeekAtLastError() ) ? Status::OK() : Status( common::ONNXRUNTIME, common::FAIL );
+}
+
+__global__ void ComputeOffsetOfMatrixAB(const int32_t* row_sum,
+                                        const int32_t* col_sum,
+                                        int32_t* output,
+                                        const int8_t* a_offset,
+                                        const int8_t* b_offset,
+                                        int32_t K,
+                                        int32_t N) {
+  for (int32_t i = threadIdx.x; i < N; i += blockDim.x) {
+    *(output + blockIdx.x * N + i) = K * (*a_offset) * (*b_offset) - row_sum[blockIdx.x] - col_sum[i];
+  }
+}
+
+__global__ void ComputeOffsetOfMatrixA(const int32_t* col_sum,
+                                       int32_t* output,
+                                       int32_t N) {
+  for (int32_t i = threadIdx.x; i < N; i += blockDim.x) {
+    *(output + blockIdx.x * N + i) = -col_sum[i];
+  }
+}
+
+__global__ void ComputeOffsetOfMatrixB(const int32_t* row_sum,
+                                       int32_t* output,
+                                       int32_t N) {
+  for (int32_t i = threadIdx.x; i < N; i += blockDim.x) {
+    *(output + blockIdx.x * N + i) = -row_sum[blockIdx.x];
+  }
+}
+
+Status OffsetOutput(const int32_t* row_sum,
+                    const int32_t* col_sum,
+                    int32_t* output,
+                    const int8_t* a_offset,
+                    const int8_t* b_offset,
+                    const MatMulComputeHelper& helper) {
+  if (a_offset && b_offset) {
+    for (size_t batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+      ComputeOffsetOfMatrixAB<<<static_cast<int>(helper.M()), GridDim::maxThreadsPerBlock, 0>>>(
+          row_sum + batch * helper.M(),
+          col_sum + batch * helper.N(),
+          output + helper.OutputOffsets()[batch],
+          a_offset,
+          b_offset,
+          static_cast<int32_t>(helper.K()),
+          static_cast<int32_t>(helper.N()));
+    }
+  } else if (a_offset) {
+    for (size_t batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+      ComputeOffsetOfMatrixA<<<static_cast<int>(helper.M()), GridDim::maxThreadsPerBlock, 0>>>(
+          col_sum + batch * helper.N(),
+          output + helper.OutputOffsets()[batch],
+          static_cast<int32_t>(helper.N()));
+    }
+  } else if (b_offset) {
+    for (size_t batch = 0; batch < helper.OutputOffsets().size(); batch++) {
+      ComputeOffsetOfMatrixB<<<static_cast<int>(helper.M()), GridDim::maxThreadsPerBlock, 0>>>(
+          row_sum + batch * helper.M(),
+          output + helper.OutputOffsets()[batch],
+          static_cast<int32_t>(helper.N()));
+    }
+  }
+
+  return CUDA_CALL( cudaPeekAtLastError() ) ? Status::OK() : Status( common::ONNXRUNTIME, common::FAIL );
+}
+
+__global__ void PadMatrixInLeadingDimensionKernel(const int8_t* src, int8_t* dst, int col_src, int col_dst) {
+  for (int32_t i = threadIdx.x; i < col_src; i += blockDim.x) {
+    *(dst + blockIdx.x * col_dst + i) = *(src + blockIdx.x * col_src + i);
+  }
+}
+
+Status PadMatrixInLeadingDimension(const int8_t* src, int8_t* dst, int64_t row, int64_t col, int64_t pad_size) {
+  PadMatrixInLeadingDimensionKernel<<<static_cast<int>(row), GridDim::maxThreadsPerBlock, 0>>>(
+      src,
+      dst,
+      static_cast<int>(col),
+      static_cast<int>(col + pad_size));
+
+  return CUDA_CALL(cudaPeekAtLastError()) ? Status::OK() : Status(common::ONNXRUNTIME, common::FAIL);
+  ;
+}
+}  // namespace cuda
+}  // namespace onnxruntime

onnxruntime/core/providers/cuda/math/matmul_integer.cuh

@@ -0,0 +1,26 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+#pragma once
+
+#include "matmul_integer.h"
+#include "core/providers/cpu/math/matmul_helper.h"
+#include "core/providers/cuda/cuda_common.h"
+#include "core/providers/cuda/shared_inc/cuda_utils.h"
+
+namespace onnxruntime {
+namespace cuda {
+
+Status ReduceRowSumOnMatrixA(const int8_t* matrix, int32_t* row_sum, const int8_t* offset, const MatMulComputeHelper& helper);
+Status ReduceColSumOnMatrixB(const int8_t* matrix, int32_t* col_sum, const int8_t* offset, const MatMulComputeHelper& helper);
+Status OffsetOutput(const int32_t* row_sum,
+                    const int32_t* col_sum,
+                    int32_t* output,
+                    const int8_t* a_offset,
+                    const int8_t* b_offset,
+                    const MatMulComputeHelper& helper);
+
+Status PadMatrixInLeadingDimension(const int8_t* src, int8_t* dst, int64_t row, int64_t col, int64_t pad_size);
+
+}  // namespace cuda
+}  // namespace onnxruntime

onnxruntime/core/providers/cuda/math/matmul_integer.h

@@ -0,0 +1,44 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+#pragma once
+
+#include "core/providers/cuda/cuda_common.h"
+
+namespace onnxruntime {
+namespace cuda {
+template <typename T1, typename T2>
+class MatMulInteger final : public CudaKernel {
+  using Base = CudaKernel;
+
+ public:
+  MatMulInteger(const OpKernelInfo& info) : CudaKernel(info) {
+    has_a_zero_point_ = false;
+    has_b_zero_point_ = false;
+    if (info.GetInputCount() > 2) {
+      has_a_zero_point_ = true;
+    }
+    if (info.GetInputCount() > 3) {
+      has_b_zero_point_ = true;
+    }
+  }
+
+  Status ComputeInternal(OpKernelContext* context) const override;
+
+ private:
+  // pad matrix and B to make their leading dimension be multiples of *align_size*
+  Status PadMatrix(
+      int64_t row,
+      int64_t col,
+      int64_t align_size,
+      const int8_t*& src,
+      int64_t& pad_size,
+      IAllocatorUniquePtr<int8_t>& temp_mem_holder) const;
+
+ private:
+  bool has_a_zero_point_;
+  bool has_b_zero_point_;
+};
+
+}  // namespace cuda
+}  // namespace onnxruntime

onnxruntime/test/providers/cpu/math/matmul_integer_test.cc

@@ -2,6 +2,7 @@
 // Licensed under the MIT License.
 
 #include "gtest/gtest.h"
+#include "test/common/cuda_op_test_utils.h"
 #include "test/providers/provider_test_utils.h"
 
 #include "core/common/common.h"
@@ -32,6 +33,209 @@ TEST(MatmulIntegerOpTest, MatMulInteger) {
   test.AddOutput<int32_t>("T3", {1, 1}, {-1});
   test.Run();
 }
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 4},
+                        {-3, 7, 5, -6,
+                         4, -5, 8, 7});
+  test.AddInput<int8_t>("T2",
+                        {4, 4},
+                        {5, -3, 7, 8,
+                         -6, -8, -3, 6,
+                         7, 9, 9, -5,
+                         8, 7, -6, 7});
+  test.AddInput<int8_t>("a_zero_point", {}, {5});
+  test.AddInput<int8_t>("b_zero_point", {}, {5});
+  test.AddOutput<int32_t>("T3",
+                          {2, 4},
+                          {-55, 16, 89, -44,
+                           122, 154, 68, -39});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t_A_ND) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 2, 4},
+                        {-3, 7, 5, -6,
+                         4, -5, 8, 7,
+
+                         7, -4, 3, 6,
+                         -4, -5, 5, 7});
+
+  test.AddInput<int8_t>("T2",
+                        {4, 3},
+                        {5, -3, 7,
+                         8, -6, -8,
+                         -3, 6, 7,
+                         9, 9, -5});
+
+  test.AddInput<int8_t>("a_zero_point", {}, {3});
+  test.AddInput<int8_t>("b_zero_point", {}, {4});
+  test.AddOutput<int32_t>("T3",
+                          {2, 2, 3},
+                          {-49, -39, 21,
+                           -46, 103, 78,
+
+                           -9, 57, 69,
+                           -33, 153, 45});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t_B_ND) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 4},
+                        {-3, 7, 5, -6,
+                         4, -5, 8, 7});
+  test.AddInput<int8_t>("T2",
+                        {2, 4, 3},
+                        {5, -3, 7,
+                         8, -6, -8,
+                         -3, 6, 7,
+                         9, 9, -5,
+
+                         5, -3, 7,
+                         8, -6, -8,
+                         -3, 6, 7,
+                         9, 9, -5});
+  test.AddInput<int8_t>("a_zero_point", {}, {1});
+  test.AddInput<int8_t>("b_zero_point", {}, {2});
+  test.AddOutput<int32_t>("T3",
+                          {2, 2, 3},
+                          {-45, -61, -11,
+                           -20, 103, 68,
+
+                           -45, -61, -11,
+                           -20, 103, 68});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t_A_ND_B_ND) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 2, 4},
+                        {-3, 7, 5, -6,
+                         4, -5, 8, 7,
+
+                         -3, 7, 5, -6,
+                         4, -5, 8, 7});
+  test.AddInput<int8_t>("T2",
+                        {2, 4, 4},
+                        {5, -3, 7, 8,
+                         -6, -8, -3, 6,
+                         7, 9, 9, -5,
+                         8, 7, -6, 7,
+
+                         5, -3, 7, 8,
+                         -6, -8, -3, 6,
+                         7, 9, 9, -5,
+                         8, 7, -6, 7});
+  test.AddInput<int8_t>("a_zero_point", {}, {5});
+  test.AddInput<int8_t>("b_zero_point", {}, {5});
+  test.AddOutput<int32_t>("T3",
+                          {2, 2, 4},
+                          {-55, 16, 89, -44,
+                           122, 154, 68, -39,
+
+                           -55, 16, 89, -44,
+                           122, 154, 68, -39});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t_A_Has_Zero_Point) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 2, 4},
+                        {-3, 7, 5, -6,
+                         4, -5, 8, 7,
+
+                         -3, 7, 5, -6,
+                         4, -5, 8, 7});
+  test.AddInput<int8_t>("T2",
+                        {2, 4, 4},
+                        {0, -8, 2, 3,
+                         -11, -13, -8, 1,
+                         2, 4, 4, -10,
+                         3, 2, -11, 2,
+
+                         0, -8, 2, 3,
+                         -11, -13, -8, 1,
+                         2, 4, 4, -10,
+                         3, 2, -11, 2});
+  test.AddInput<int8_t>("a_zero_point", {}, {5});
+  test.AddOutput<int32_t>("T3",
+                          {2, 2, 4},
+                          {-55, 16, 89, -44,
+                           122, 154, 68, -39,
+
+                           -55, 16, 89, -44,
+                           122, 154, 68, -39});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
+TEST(MatmulIntegerOpTest, MatMulInteger_int8_t_No_Zero_Point) {
+  if (!DefaultCudaExecutionProvider() || !HasCudaEnvironment(530 /*min_cuda_architecture*/)) return;
+
+  OpTester test("MatMulInteger", 10);
+  test.AddInput<int8_t>("T1",
+                        {2, 2, 4},
+                        {-8, 2, 0, -11,
+                         -1, -10, 3, 2,
+
+                         -8, 2, 0, -11,
+                         -1, -10, 3, 2});
+  test.AddInput<int8_t>("T2",
+                        {2, 4, 4},
+                        {0, -8, 2, 3,
+                         -11, -13, -8, 1,
+                         2, 4, 4, -10,
+                         3, 2, -11, 2,
+
+                         0, -8, 2, 3,
+                         -11, -13, -8, 1,
+                         2, 4, 4, -10,
+                         3, 2, -11, 2});
+  test.AddOutput<int32_t>("T3",
+                          {2, 2, 4},
+                          {-55, 16, 89, -44,
+                           122, 154, 68, -39,
+
+                           -55, 16, 89, -44,
+                           122, 154, 68, -39});
+
+  std::vector<std::unique_ptr<IExecutionProvider>> execution_providers;
+  execution_providers.push_back(DefaultCudaExecutionProvider());
+  test.Run(OpTester::ExpectResult::kExpectSuccess, "", {}, nullptr, &execution_providers);
+}
+
 TEST(MatmulIntegerOpTest, MatMulInteger_WithZero_ZeroPoint) {
   OpTester test("MatMulInteger", 10);
   test.AddInput<uint8_t>("T1", {4, 3}, {11, 7, 3, 10, 6, 2, 9, 5, 1, 8, 4, 0});