Randomized Singular Value Decomposition (RSVD)

• RSVD is a class of probabilistic algorithms that efficiently approximates the low-rank SVD via random projections, oversampling, and power iterations.
• It leverages techniques like QR decomposition and controlled error bounds to offer substantial computational and memory advantages over classical SVD.
• Widely used in tensor network contraction, inverse problems, and data analytics, RSVD can achieve significant speedups, e.g., over 100-fold in TRG applications.

Randomized Singular Value Decomposition (RSVD) is a class of probabilistic algorithms for approximating the low-rank or truncated SVD of large matrices, leveraging random projections and sketching to accelerate the core computations underlying dimensionality reduction, model order reduction, inverse problems, and high-dimensional data analysis. RSVD has established itself as a standard tool in numerical linear algebra for problems where classical deterministic SVD algorithms are computationally prohibitive or memory-bound, particularly in applications such as tensor network contraction, regularized inversion, and large-scale data analytics.

1. Algorithmic Foundations of RSVD

The canonical RSVD algorithm proceeds as follows: Given a matrix $A \in \mathbb{R}^{m \times n}$, target rank $k\ll\min(m,n)$, oversampling parameter $p$ (typically $p\simeq k$), and number of power iterations $q\ge0$, RSVD computes an approximate rank-$k$ SVD $A\approx U \Sigma V^T$ using the following steps (Morita et al., 2017):

1. Random projection: Draw a test matrix $\Omega\in\mathbb{R}^{n\times(k+p)}$ with entries $\Omega_{ij}\sim\mathcal{N}(0,1)$.
2. Sampling (“sketch”): Compute $Y_0 = A\Omega \in \mathbb{R}^{m\times(k+p)}$.
3. (Power iteration, optional): For $k\ll\min(m,n)$0, alternate:
   • $k\ll\min(m,n)$1,
   • Orthonormalize $k\ll\min(m,n)$2 to get $k\ll\min(m,n)$3,
   • $k\ll\min(m,n)$4,
   • Orthonormalize $k\ll\min(m,n)$5 to get $k\ll\min(m,n)$6, and finally set $k\ll\min(m,n)$7 (or $k\ll\min(m,n)$8 if $k\ll\min(m,n)$9).
4. Orthonormalization: $p$0 via thin QR ($p$1).
5. Small projection: Form $p$2.
6. Truncated SVD: Compute SVD $p$3, retain top $p$4 components.
7. Recovery of singular vectors: $p$5 ($p$6, $p$7, $p$8 diagonal $p$9).

This construction guarantees that $p\simeq k$0, with the approximation error controlled by the singular spectrum of $p\simeq k$1, $p\simeq k$2, and $p\simeq k$3.

2. Complexity, Memory, and Practical Acceleration

RSVD offers significant practical computational and memory advantages over classical deterministic SVD methods, especially in large-scale applications (Morita et al., 2017):

• Computational scaling: Each matrix–matrix multiply (e.g., $p\simeq k$4 or $p\simeq k$5) costs $p\simeq k$6 flops. For example, in coarse-graining steps in TRG with bond dimension $p\simeq k$7, $p\simeq k$8, $p\simeq k$9, $q\ge0$0, yielding overall $q\ge0$1 scaling—compared to $q\ge0$2 for full SVD on the same matrix.
• Memory usage: RSVD-based workflows hold only third-order tensors and $q\ge0$3 intermediate blocks in memory (using “loop blocking”, i.e., small blockwise contractions), as opposed to $q\ge0$4 for explicitly constructing and storing fourth-order tensors required in standard SVD-based TRG frameworks (Morita et al., 2017).
• Empirical speedup: For TRG applied to the 2D Ising model at $q\ge0$5, RSVD-based contraction is over 100-fold faster than full SVD, while achieving machine-precision accuracy (Morita et al., 2017).

Such scaling benefits extend beyond tensor networks to regularized inversion, imaging, and matrix completion settings (Li et al., 2023, Feng et al., 2018).

3. Error Bounds and Probabilistic Guarantees

RSVD delivers controlled approximation errors with high probability, which admit sharp theoretical analyses:

• Frobenius-norm error, $q\ge0$6:

$q\ge0$7

where the $q\ge0$8 are singular values of $q\ge0$9.

• With power iterations ($k$0): Each singular value $k$1 is replaced by $k$2, so the tail sum decays much faster. Empirically, $k$3 suffices at the critical point of the 2D Ising model for $k$4; errors decay as

$k$5

with $k$6 (Morita et al., 2017).

• Spectral norm bounds: With $k$7 power iterations,

$k$8

matching deterministic SVD up to the constant for reasonable $k$9 and $A\approx U \Sigma V^T$0 (Morita et al., 2017, Li, 2024).

High-fidelity matching to full SVD is empirically validated even in challenging spectral regimes (e.g., at phase transitions in statistical mechanics; see Ising model at criticality in (Morita et al., 2017)).

4. Integration into Scientific and Data Analysis Workflows

RSVD is now a standard drop-in for classical SVD in numerous large-scale computations:

• Tensor Renormalization Group (TRG):
  • Coarse-graining of a square lattice contracts rank-four tensors $A\approx U \Sigma V^T$1 (size $A\approx U \Sigma V^T$2) by reshaping to $A\approx U \Sigma V^T$3, then performing SVD for truncation.
  • RSVD-based TRG avoids explicit construction of $A\approx U \Sigma V^T$4 or its dense SVD, sequentially contracting random third-order tensors, with all intermediates $A\approx U \Sigma V^T$5 in memory (Morita et al., 2017).
• Inverse Problems and Imaging:
  • Embedded in full-waveform inversion, regularized linear inversion, and imaging pipelines, RSVD provides efficient dimension reduction of velocity increment or system matrices, typically within an augmented Lagrangian or regularization scheme (Li et al., 2023, Ito et al., 2019).
  • Inversion accuracy is preserved while convergence and noise suppression are improved; e.g., RSVD-WTNNR-iALM outperforms Tikhonov-regularized FWI in both convergence and noise resilience (Li et al., 2023).
• Matrix Completion and Machine Learning:
  • Accelerated SVT algorithms replace Lanczos or Krylov-projected SVDs with RSVD-BKI or fast RSVD-PI variants, yielding 6–15× speedup in image inpainting and recommender systems, at no loss in solution accuracy (Feng et al., 2018).

5. Parameter Selection, Algorithmic Variants, and Implementation

Optimal use of RSVD involves selection and tuning of key parameters:

Parameter	Standard Value/Role	Effect
Target rank $A\approx U \Sigma V^T$6	Target subspace dim (e.g., $A\approx U \Sigma V^T$7)	Sets SVD truncation.
Oversampling $A\approx U \Sigma V^T$8	$A\approx U \Sigma V^T$9 recommended	Increases capture probability for $\Omega\in\mathbb{R}^{n\times(k+p)}$0.
Power iterations $\Omega\in\mathbb{R}^{n\times(k+p)}$1	$\Omega\in\mathbb{R}^{n\times(k+p)}$2–$\Omega\in\mathbb{R}^{n\times(k+p)}$3 (depends on spectrum)	Amplifies decay in slowly-varying spectra.
Loop-block size	$\Omega\in\mathbb{R}^{n\times(k+p)}$4 in TRG	Ensures no intermediates exceed $\Omega\in\mathbb{R}^{n\times(k+p)}$5.

Advanced variants include:

• Adaptive rank estimation and stopping criteria based on estimated energy capture, as in R3SVD (Ji et al., 2016).
• Memory-aware block-wise implementations to ensure constant working memory independent of final rank (Ji et al., 2016).
• Specialized numerical kernels for GPU/CPU acceleration (e.g., RSVDPACK (Voronin et al., 2015)).

The algorithm is robust to implementation details; e.g., choice of Gaussian vs. SRFT test matrices, block size, and reorthogonalization frequency. Empirical studies confirm the theoretical scaling and error guarantees across diverse applications (Morita et al., 2017, Li et al., 2023, Feng et al., 2018).

6. Application to Tensor and Higher-Order Decompositions

RSVD generalizes naturally to tensor-network algorithms by exploiting sequential random projections along tensor cores:

• Tensor contraction: In TRG, random tensors are contracted into sequence of third-order S-tensors, allowing on-the-fly column space sampling without constructing high-rank tensors (Morita et al., 2017).
• Complexity for higher-order: For a $\Omega\in\mathbb{R}^{n\times(k+p)}$6-dimensional tensor network, each RSVD step reduces the contraction cost by a full factor of $\Omega\in\mathbb{R}^{n\times(k+p)}$7 (the mode size), saving both memory and compute (Morita et al., 2017, Huber et al., 2017).

The integration of RSVD in tensor decompositions is particularly beneficial at challenging points such as criticality, where slow singular value decay would severely degrade deterministic algorithms due to the unfavorable scaling of traditional SVD (Morita et al., 2017).

7. Outlook: Strengths, Limitations, and Empirical Best Practices

RSVD is an indispensable tool for large-scale approximate linear algebra for several reasons:

• Scaling: Drastic reduction in arithmetic cost (e.g., $\Omega\in\mathbb{R}^{n\times(k+p)}$8 vs. $\Omega\in\mathbb{R}^{n\times(k+p)}$9 in TRG) and memory consumption ($\Omega_{ij}\sim\mathcal{N}(0,1)$0 instead of $\Omega_{ij}\sim\mathcal{N}(0,1)$1) (Morita et al., 2017).
• Accuracy: Controlled approximation error, even at physically difficult points (critical Ising model), provided oversampling $\Omega_{ij}\sim\mathcal{N}(0,1)$2 and moderate $\Omega_{ij}\sim\mathcal{N}(0,1)$3 (Morita et al., 2017).
• Universality: Effective in a wide spectrum of settings, from geometric inverse problems and machine learning to high-dimensional tensor contractions (Li et al., 2023, Feng et al., 2018).
• Practicality: Straightforward to implement, compatible with most high-level languages and optimized BLAS/LAPACK environments (Voronin et al., 2015).

However, selection of $\Omega_{ij}\sim\mathcal{N}(0,1)$4 and $\Omega_{ij}\sim\mathcal{N}(0,1)$5 remains spectrum-dependent. Insufficient oversampling or power iterations can compromise the capture of slow-decay singular directions, particularly in physically critical or ill-posed problems. Loop-blocking and attention to contraction order are crucial for avoiding hidden $\Omega_{ij}\sim\mathcal{N}(0,1)$6 intermediates in tensor workflows (Morita et al., 2017). Nonetheless, extensive empirical benchmarks confirm RSVD-based pipelines yield numerically stable, reproducible, and highly efficient solutions for large-scale matrix and tensor computations.