PaperlayerView exact revision on arXiv

Source-linked AI summary

Double-Precision Matrix Multiplication Emulation via Ozaki-II Scheme with FP8 Quantization

Yuki Uchino, Katsuhisa Ozaki, Toshiyuki Imamura

arXiv:2603.10634v2cs.DC

TL;DR

Ozaki-II’s original algorithm cannot directly exploit FP8 MMA units, although some emerging architectures prioritize FP8 and reduce INT8 resources. The paper combines a Karatsuba-based extension with modular reduction for selected moduli to reduce FP8 matrix multiplications. It presents the approach as a practical option for architectures with limited INT8 capability or FP8-dominated computation.

  • Problem

    Ozaki-II’s original algorithm cannot directly use FP8 MMA units, despite emerging architectures prioritizing FP8 and reducing INT8 resources.

    Source: B. Problem Statement and Related Work · p. 2

  • Method

    The method combines a Karatsuba-based extension with modular reduction for selected moduli to emulate FP64 DGEMM on FP8 MMA units.

    Source: §VI · p. 10

  • Results

    The hybrid construction reduces the number of required FP8 matrix multiplications.

    Source: §VI · p. 10

  • Takeaways & Limitations

    The FP8-based emulation provides a practical option for architectures dominated by low-precision floating-point computation or with limited INT8 capability.

    Source: §VI · p. 11

    Direct FP8 E4M3 quantization of Ozaki-II restricts the CRT range, making accurate FP32 and FP64 emulation difficult.

    Source: A. Limitations of Direct FP8… · p. 3

Abstract

from arXiv · show

In this paper, we propose a method for emulating double-precision general matrix--matrix multiplication (DGEMM), a fundamental and performance-critical kernel in many high-performance computing applications. Ozaki-I and Ozaki-II are established DGEMM emulation schemes via low-precision matrix multiply-accumulate (MMA) units. For the Ozaki-I scheme, INT8-, FP8-, and FP16-based implementations have been proposed, all of which can be realized based on the same underlying algorithmic structure. In contrast, although INT8-based implementations of the Ozaki-II scheme have been reported, the original algorithm cannot be directly adapted to exploit FP8 MMA units. In several recent architectures, such as NVIDIA Blackwell Ultra and NVIDIA Rubin, INT8 performance has been reduced, making reliance on INT8 alone insufficient. Therefore, we introduce a novel technique to demonstrate DGEMM emulation based on the Ozaki-II scheme that operates on FP8 MMA units. Compared to the FP8-based Ozaki-I scheme, our method significantly reduces the computational cost and enables efficient FP64 emulation.

I. INTRODUCTION

Low-precision MMA units offer a route to accelerating FP64 DGEMM, but Ozaki-II’s original structure does not directly support FP8. The paper introduces an FP8-based Ozaki-II approach and compares it with existing emulation methods.

B. Problem Statement and Related Work

  • Ozaki-I and Ozaki-II decompose FP64 matrices into low-precision matrices and combine their products, and existing Ozaki-based DGEMM methods exploit INT8 MMA units.The fixed-point representation of the decomposed matrices is particularly suited to INT8 MMA units.Source: B. Problem Statement and Related Work · p. 2
  • FP8-based Ozaki-II addresses architectures that prioritize FP8 and reduce INT8 resources, while INT8 remains relevant on systems with comparable throughput or continued INT8 support.The paper describes this hardware balance as heterogeneous and treats both emulation paths as practically important.Source: A. Background · p. 1, B. Problem Statement and Related Work · p. 2

C. Our Contributions

  • The paper identifies why a direct INT8-to-FP8 substitution cannot preserve exactness and presents FP8 as a design choice shaped by exactness under FP32 accumulation and practical performance.The stated assumption is INT8 inputs with INT32 accumulation and FP8 E4M3 inputs with FP32 accumulation.Source: C. Our Contributions · p. 2
  • The proposed FP8-based Ozaki-II approach combines a Karatsuba-based extension with modular reduction for selected moduli to reduce FP8 matrix multiplications while meeting target precision.The paper also compares methods by effective precision, performance models, working memory, and numerical experiments.Source: C. Our Contributions · p. 2

II. OZAKI-II SCHEME

  • For direct Ozaki-II CRT emulation, the paper states that at least 14 moduli are required in principle to emulate FP64 arithmetic.With N = 14, the selected modulus product satisfies P/2 > 2^109 > 2^(53+53).Source: §II · p. 3

III. PROPOSED METHOD USING FP8 E4M3

Directly applying FP8 E4M3 quantization to Ozaki-II restricts the usable moduli and the CRT reconstruction range. The resulting range is insufficient for accurate FP32 emulation.

A. Limitations of Direct FP8 Quantization

  • Direct FP8 E4M3 quantization assumes k ≤ 2^16 so matrix multiplication has no rounding error, and requires pℓ ≤ 32 for the residues to be represented.FP8 E4M3 exactly represents consecutive integers from −16 to 16.Source: A. Limitations of Direct FP8… · p. 3
  • With the available FP8 E4M3 moduli, P/2 < 2^47, limiting CRT reconstruction enough to make accurate FP32 emulation difficult.The passage states that this limitation also makes accurate FP64 emulation difficult.Source: A. Limitations of Direct FP8… · p. 3
  • Even using all listed candidate moduli does not raise the modulus product enough for the required dynamic range.The candidate set is {32, 31, 29, 27, 25, 23, 19, 17, 13, 11, 7}.Source: A. Limitations of Direct FP8… · p. 3

B. Karatsuba-Based Extension

The Karatsuba-based extension represents the decomposed integer matrices as sums of FP8 E4M3 matrices and computes their products through three error-free FP8 matrix multiplications. It requires at least 13 moduli to reach precision comparable to INT8-based emulation with 14 moduli.

B. Karatsuba-Based Extension

  • The Karatsuba extension represents each decomposed matrix as a sum of two FP8 E4M3 matrices and applies Karatsuba reconstruction to their product.The resulting computation uses three error-free FP8 matrix multiplications.Source: B. Karatsuba-Based Extension · p. 3
  • The moduli are selected as pairwise coprime integers no greater than 513, with a greedily constructed descending sequence.The passage gives an example sequence beginning at 481.Source: B. Karatsuba-Based Extension · p. 3
  • N ≥ 13 is required to achieve precision comparable to INT8-based emulation with 14 moduli.The cited bound is P/2 > 2^115 > 2^(53+53).Source: B. Karatsuba-Based Extension · p. 3

C. Modular Reduction without Karatsuba

Modular reduction uses square moduli to avoid Karatsuba reconstruction for selected products, while the hybrid method combines this with Karatsuba for other moduli. The hybrid reduces the required moduli and FP8 matrix multiplications relative to the Karatsuba-based method.

C. Modular Reduction without Karatsuba

  • Karatsuba requires sums to be represented rigorously in FP8 E4M3, restricting the effective representable range.The modular-reduction formulation is introduced to relax this constraint.Source: C. Modular Reduction without Karatsuba · p. 3
  • For square moduli pℓ = s^2, modular arithmetic allows the product to be evaluated with three FP8 matrix multiplications without Karatsuba reconstruction.The possible square-modulus values are limited, so the hybrid method uses Karatsuba for remaining moduli.Source: C. Modular Reduction without Karatsuba · p. 3

D. Hybrid Method

  • The hybrid construction prioritizes pairwise coprime square moduli in descending order from 1089, then adds pairwise coprime nonsquares.It applies modular reduction to the resulting matrix products and Karatsuba to the remaining moduli.Source: D. Hybrid Method · p. 4
  • N ≥ 12 suffices for precision comparable to INT8-based emulation with N = 14, while the hybrid method reduces moduli and FP8 matrix multiplications versus the Karatsuba-based method.The passage reports that the reduced multiplication count improves computational efficiency.Source: D. Hybrid Method · p. 4

E. Method for Conversion to Integer Matrices

The conversion method turns the input floating-point matrices into integer matrices and chooses scaling vectors using a computable upper bound. In accurate mode, it focuses on FP8 matrix multiplication, with FP32 products exact and rounding error arising during accumulation.

E. Method for Conversion to Integer Matrices

  • The method first converts A and B into integer matrices A′ and B′ satisfying condition (3), using scaling vectors μ and ν to estimate an upper bound.The passage distinguishes fast and accurate computing modes.Source: E. Method for Conversion to Integer… · p. 4
  • In accurate mode, the upper bound is estimated by matrix multiplication using INT8 MMA units, and scaling vectors are chosen so that 2^μ_i w̄_ij 2^ν_j < P for all (i, j).Source: E. Method for Conversion to Integer… · p. 4
  • The scaled matrices are cast in round-up mode to FP8 E4M3, and each matrix-product term is formed in FP32 without rounding error.The passage states that rounding errors arise only in accumulation.Source: E. Method for Conversion to Integer… · p. 4
  • For accurate accumulation, the computed result is multiplied by an upper-bound factor replaced with 1 + k2^-24, using round-up mode.This replacement avoids division.Source: E. Method for Conversion to Integer… · p. 4

F. Rationale for Choosing FP8 over FP16, BF16, and FP4

The rationale for FP8 is that its E4M3 significand can support error-free FP32 accumulation over a wider dot-product length than comparable FP16/BF16 precision without decomposition. The passage also contrasts FP8's stated peak-throughput envelope with FP16/BF16 and identifies conditions that limit recursive FP4 approaches.

F. Rationale for Choosing FP8 over FP16, BF16, and FP4

  • FP32 accumulation requires up to 2β + ⌈log2 k⌉ bits in the worst case, against FP32's 24-bit significand.Here β is the input format's significand-bit count, including the implicit leading bit.Source: F. Rationale for Choosing FP8 over… · p. 5
  • FP8 E4M3 uses β = 4 and permits error-free FP32 accumulation for k ≤ 2^16, whereas comparable FP16/BF16 precision requires k-blocking.The passage notes that small k-blocks can reduce throughput and make the kernel memory-bound.Source: F. Rationale for Choosing FP8 over… · p. 5
  • NVIDIA Rubin's specified dense throughput is 17.5 PFLOP/s for FP8 versus 4.0 PFLOP/s for FP16/BF16, a peak ratio of 4.375.The comparison applies when FP8 matrix multiplication preserves the required exactness.Source: F. Rationale for Choosing FP8 over… · p. 5
  • Exact intermediate sums make straightforward recursive FP4 Karatsuba or Toom–Cook inapplicable, though a single-level Karatsuba implementation can use three FP4 multiplications per FP8 multiplication.The passage says this construction may become advantageous if future FP4 throughput exceeds three times FP8 throughput.Source: F. Rationale for Choosing FP8 over… · p. 5

IV. COMPARISON WITH PREVIOUS WORK

The proposed FP8-based Ozaki-II scheme reaches FP64-level precision with far fewer matrix multiplications than FP8-based Ozaki-I, though it requires more multiplications than INT8-based Ozaki-II.

A. Number of Low-Precision Matrix Multiplications

  • INT8-based Ozaki-I reduces the count to S(S + 1)/2 by neglecting smaller terms, but this modification introduces accuracy loss.The full Ozaki-I scheme uses S^2 matrix multiplications.Source: A. Number of Low-Precision Matrix… · p. 5
  • At least 53 bits of precision require 36 FP8 matrix multiplications in proposed Ozaki-II, versus 121 in FP8-based Ozaki-I.The proposed scheme uses N ≥12; accurate mode requires one additional FP8 multiplication.Source: A. Number of Low-Precision Matrix… · p. 5
  • The proposed approach achieves comparable or higher precision than FP8-based Ozaki-I while significantly reducing multiplication count.Table II summarizes effective precision and multiplication counts for the compared methods.Source: A. Number of Low-Precision Matrix… · p. 5, Table II · p. 6
  • The proposed FP8-based method requires approximately 2.5× the matrix multiplications of INT8-based Ozaki-II.INT8-based Ozaki-II requires N ≥14 for at least 53 bits; its multiplication count is N in fast mode and N +1 in accurate mode.Source: A. Number of Low-Precision Matrix… · pp. 5–6

B. Performance Modeling and Throughput Comparison

The paper develops analytic performance models for real-valued INT8-based and proposed FP8-based Ozaki-II emulation. Its predictions indicate that FP8 emulation could exceed a 200 TFLOP/s reference under Rubin-like specifications, while INT8 may remain faster when FP8 throughput is only about twice INT8 throughput.

B. Performance Modeling and Throughput Comparison

  • The authors develop an analytic performance model for real-valued INT8-based Ozaki-II GEMM emulation, addressing the absence of an explicit model in prior real-valued studies.They also construct performance models for the proposed FP8-based method.Source: B. Performance Modeling and Throughput… · p. 6
  • The INT8 fast-mode model treats sustained memory bandwidth and BLAS GEMM throughput as platform inputs, with correction parameter c summarizing FP64- and control-heavy overheads.The model's correction parameter can reflect overhead components invoked once per modulus.Source: B. Performance Modeling and Throughput… · p. 6
  • The FP8 workflow forms scaled FP8 E4M3 matrices, computes residue products and modular reductions, then reconstructs the approximate product.The workflow assumes N < 34, with square moduli p1 through p6.Source: B. Performance Modeling and Throughput… · p. 6
  • Under Rubin-like hardware specifications, predicted FP8-based Ozaki-II throughput could substantially exceed the 200 TFLOP/s FP64 DGEMM emulation reference.If FP8 matrix multiplication is only about twice as fast as INT8, the model indicates INT8-based emulation will likely remain faster.Source: B. Performance Modeling and Throughput… · p. 7, Figure 2 · p. 7, Figure 1 · p. 7
  • Figures 1 and 2 show predicted-throughput heatmaps for INT8- and FP8-based Ozaki-II, respectively, across sustained GEMM throughput and memory bandwidth.Each figure compares fast and accurate modes at 16384 × 16384 × 16384 dimensions.Source: Figure 1 · p. 7, Figure 2 · p. 7

C. Working Memory Footprint Comparison

FP8-based Ozaki-II has a larger working-memory footprint than INT8-based Ozaki-II, although blocking can reduce workspace for both methods. Keeping k large while blocking m and n is presented as a practical way to limit performance loss.

C. Working Memory Footprint Comparison

  • The INT8 working-memory model is Wi8(m, n, k, N) = (mk + kn + 5mn)N + 2(m + n), while the FP8 model includes an additional 2Nmn term.The FP8 model also scales (mk + kn + 4mn) by MN.Source: C. Working Memory Footprint Comparison · p. 7
  • At m = n = k = 16384, INT8-based Ozaki-II with N = 14 requires 27 GB, while FP8-based Ozaki-II with N = 12 requires 55 GB.Both methods require substantial working memory, making workspace reduction a practical challenge in large HPC applications.Source: C. Working Memory Footprint Comparison · p. 7
  • The FP8-based method requires more working memory, mainly because residues use multiple FP8 matrices and reduced residues use INT16 rather than INT8.The FP8 exponent bits are effectively unused when storing fixed-point residues.Source: C. Working Memory Footprint Comparison · p. 7
  • Blocking emulation can reduce working memory for both INT8- and FP8-based methods.The working-memory discussion defines workspace as temporary computation buffers, excluding input and output matrices.Source: C. Working Memory Footprint Comparison · p. 7
  • Blocking only m and n while keeping k large can reduce workspace while suppressing performance degradation.Reducing k can lower arithmetic intensity and underutilize MMA units; the models omit edge effects and additional accumulation overhead across k-blocks.Source: C. Working Memory Footprint Comparison · p. 8

V. NUMERICAL RESULTS

Experiments use an open-source library supporting INT8- and FP8-based Ozaki-II with bitwise reproducibility under a fixed toolchain. The evaluation spans several GPU and accelerator platforms and reports medians from timed runs.

V. NUMERICAL RESULTS

  • The open-source library supports both INT8- and FP8-based Ozaki-II and yields bitwise reproducible results under a fixed toolchain.All experiments in the section use this library.Source: §V · p. 8
  • The evaluation platforms include RTX 4090 Laptop, RTX 5080, Radeon RX 9070 XT, GH200, GB10, and B200 systems.The B200 experiments used SAKURA internet Inc.'s managed HPC cluster service, SAKURAONE.Source: §V · p. 8
  • Throughput and time-breakdown measurements use 30 warm-up runs followed by 30 timed runs, with the median reported.Native FP64 DGEMM uses cuBLAS cublasDgemm; FP8 emulation GEMMs use cuBLASLt and INT8 GEMMs use cuBLAS.Source: §V · p. 8

A. Accuracy

DGEMM emulation accuracy is compared across modes and baselines, while throughput varies sharply by platform and matrix shape. INT8-based Ozaki-II is generally faster, but FP8-based Ozaki-II offers an alternative where FP8 hardware is favored.

A. Accuracy

  • For normally distributed test matrices, both Ozaki-II variants achieve accuracy close to the cuBLAS Ozaki-I baseline using seven slices.The compared configurations use 15 INT8 moduli in fast mode and 16 in accurate mode, versus 13 FP8 moduli in fast mode and 12 in accurate mode.Source: A. Accuracy · p. 8, Figure 3 · p. 8
  • In accurate mode, FP8-based emulation has a smaller accuracy gap between k = 1024 and k = 65536 than INT8-based emulation.Round-up conversion to FP8 E4M3 inflates small values less than the INT8 ceiling conversion, enabling a tighter bound estimate; accurate mode typically yields higher accuracy than fast mode.Source: A. Accuracy · p. 8

B. Throughput

  • On low-FP64-throughput GPUs, emulation outperforms native FP64 DGEMM even at 1024-square sizes, while its advantage narrows on GH200 and B200.For sufficiently large problems, the observed speed order on all platforms except GH200 is INT8-based Ozaki-II, INT8-based Ozaki-I, FP8-based Ozaki-II, then native FP64.Source: B. Throughput · p. 9, Figure 4 · p. 9
  • On RTX 5080, INT8-based Ozaki-II reaches 4.8–15× and 7.4–24× native-DGEMM speedups at square sizes 1024 and 8192, respectively.FP8-based Ozaki-II reaches 3.7–5.0× and 4.3–9.4× at those respective sizes; INT8-based Ozaki-II is 1.3–3.1× faster across evaluated ranges.Source: B. Throughput · p. 9, Figure 5 · p. 9
  • On B200, both Ozaki-II schemes reach only 0.1–0.7× native throughput at size 1024, while at size 16384 INT8 reaches 1.3–3.9× and FP8 reaches 0.7–1.9×.INT8-based Ozaki-II is 0.9–2.3× faster than FP8-based Ozaki-II across evaluated ranges.Source: B. Throughput · p. 9, Figure 6 · p. 9
  • On B200 at 16384-square size, INT8-based Ozaki-II measures 137 TFLOP/s in fast mode and 138 TFLOP/s in accurate mode, versus 61 and 65 TFLOP/s for FP8-based Ozaki-II.The analytic model predicts 140 TFLOP/s for INT8 in both modes and 69 and 73 TFLOP/s for FP8 in fast and accurate modes, respectively.Source: B. Throughput · pp. 9–10, Figure 6 · p. 9

C. Time Breakdown

Runtime shifts from memory-bound behavior at small k toward compute-bound GEMMs as k grows, with phase shares also depending on matrix shape and hardware. The FP8 approach broadens Ozaki-II emulation to FP8-oriented platforms, although INT8 remains preferable in throughput and memory footprint where its MMA capacity is sufficient.

C. Time Breakdown

  • For large m and n, GEMMs take a growing runtime share as k increases and eventually dominate; smaller m and n make non-GEMM phases more significant.At small k, requantization and dequantization take larger shares, while quantization grows with k because it scales with input sizes mk and kn rather than output size mn.Source: C. Time Breakdown · p. 10
  • FP8-based emulation spends a larger runtime fraction in GEMMs because each modulus uses three FP8 matrix multiplications, versus one INT8 multiplication per modulus.The workflow partitions time into quantization, GEMMs, requantization, dequantization, and other overhead.Source: C. Time Breakdown · p. 10
  • RTX 5080 has a relatively large dequantization share, reflecting limited FP64 throughput; B200 has a smaller GEMM share on small problems, consistent with lower MMA-unit utilization.Similar dequantization behavior was observed on RTX 4090 Laptop, RX 9070 XT, and GB10.Source: C. Time Breakdown · p. 10, Figure 7 · p. 10, Figure 8 · p. 10

VI. CONCLUSION

  • The hybrid Ozaki-II construction combines a Karatsuba-based extension with modular reduction for selected square moduli to reduce the required FP8 matrix multiplications.This enables FP64-targeted emulation on FP8 MMA units.Source: §VI · p. 10
  • Where INT8 MMA capacity is at least about half FP8 throughput, INT8 approaches are preferable in both throughput and working-memory footprint.FP8-based emulation nevertheless extends applicability to platforms dominated by low-precision floating-point computation or with limited INT8 capability, including B300 and Rubin.Source: §VI · pp. 10–11

Related papers

Loading 2603.10634v2…