⚛️ quantum physics

Approximation and composition of functions in quantized tensor trains via orthogonal polynomial expansions

This paper presents a constructive algorithm that utilizes orthogonal polynomial expansions and Clenshaw evaluations to efficiently represent analytic functions as quantized tensor trains (QTT), enabling stable, rapidly convergent function composition in high-dimensional settings.

Original authors: Juan José Rodríguez-Aldavero, Paula García-Molina, Luca Tagliacozzo, Juan José García-Ripoll

Published 2026-09-18

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: Juan José Rodríguez-Aldavero, Paula García-Molina, Luca Tagliacozzo, Juan José García-Ripoll

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

In the modern world of science and engineering, researchers often face a daunting problem: how to describe a system with hundreds or thousands of moving parts without drowning in data. Imagine trying to map every single grain of sand on a beach; the sheer volume of information would quickly overwhelm any computer. To solve this, mathematicians and physicists have developed ways to compress this information, stripping away the unnecessary details while keeping the essential shape of the problem intact. One powerful method for doing this is called a tensor train, a technique that breaks a massive, complex object into a chain of smaller, manageable pieces. When these pieces are arranged in a specific, layered way, they form what is known as a quantized tensor train. This structure is incredibly efficient, allowing computers to handle problems that would otherwise be impossible, such as simulating the behavior of quantum particles or solving complex equations in high-dimensional spaces. However, a persistent challenge remains: how do you take a smooth, continuous function—a mathematical description of a curve or a surface—and translate it into this compressed format without losing accuracy or stability?

A team of researchers at the Institute of Fundamental Physics in Madrid has developed a new way to answer this question. They created a constructive algorithm that translates smooth, continuous functions into these compressed tensor formats by using a specific type of mathematical building block called orthogonal polynomials. Think of these polynomials as a set of standard, well-behaved curves that can be mixed together to recreate almost any smooth shape. The researchers found that by expanding a function into a sum of these curves and then carefully translating that sum into the tensor format, they could create highly accurate approximations. Their method is particularly effective for functions that are smooth and do not have sharp, jagged edges. It works by building the solution step-by-step, using a stable mathematical recipe that prevents errors from piling up, even when the calculation involves thousands of variables.

The team tested their approach on a variety of mathematical functions, ranging from simple bell-shaped curves to complex, oscillating waves. They found that for smooth functions, their method converged rapidly, meaning it reached a high level of accuracy with relatively few computational steps. In tests involving univariate functions—those with a single variable—their technique required far fewer data points to achieve the same precision as other popular methods. While other techniques often rely on randomly sampling points from a function to guess its shape, which can be inefficient and unpredictable, this new method uses the known mathematical structure of the function to build the solution directly. This deterministic approach ensures that the result is stable and reproducible. The researchers also demonstrated that their method could handle multivariate functions, which involve many variables at once, by chaining together simpler, single-variable approximations. This allowed them to tackle problems with up to 200 variables, representing a system with more than a trillion possible states, a scale that is far beyond the reach of traditional, uncompressed methods.

One of the key strengths of this new algorithm is its ability to maintain stability even as the complexity of the problem grows. In many numerical methods, increasing the number of variables or the precision of the calculation can lead to a breakdown in accuracy, where tiny errors multiply and ruin the result. The researchers showed that their use of orthogonal polynomials, combined with a specific evaluation technique known as a Clenshaw recurrence, keeps these errors in check. They observed that the method scales efficiently, meaning the time and memory required to solve the problem grow at a manageable rate rather than exploding exponentially. This is crucial for applications in quantum-inspired computing, where the goal is to simulate complex physical systems that are too large for standard computers. The team compared their results with existing state-of-the-art techniques, such as tensor cross-interpolation, and found that while their method might not always be the fastest for every single type of problem, it offers a robust and reliable alternative, especially when dealing with smooth, highly differentiable functions.

The work also highlights the importance of how the data is organized within the computer's memory. The researchers explored different ways of ordering the variables in their calculations, finding that a specific arrangement, which they called a serial order, often performed better than a more scrambled, interleaved arrangement for certain types of complex, nonlinear models. This discovery suggests that the way we structure our mathematical models can be just as important as the algorithms we use to solve them. By carefully choosing the order of operations and the type of polynomial expansion, the researchers were able to push the boundaries of what is computationally feasible, handling systems with dense interactions and strong correlations that would typically cause other methods to fail.

Ultimately, this research provides a general framework for composing functions within these compressed formats. It allows scientists to take a known function and apply it to another function that is already in a compressed state, enabling the construction of complex, layered models without ever needing to expand them into their full, unwieldy form. This capability opens the door to solving nonlinear equations and simulating intricate physical processes with a level of efficiency that was previously out of reach. The algorithms developed in this study are now available as open-source software, allowing other researchers to apply these techniques to their own problems. By turning the abstract challenge of high-dimensional data into a concrete, solvable process, this work offers a new tool for navigating the vast and complex landscapes of modern scientific computation.

Technical Summary: Approximation and Composition of Functions in Quantized Tensor Trains via Orthogonal Polynomial Expansions

Problem Statement
The central challenge addressed in this work is the construction of accurate, compressed tensorized representations of functions within the Quantized Tensor Train (QTT) format. While QTT representations offer exponential reductions in storage and computational complexity for high-dimensional objects compared to explicit grid discretizations, generating these representations for arbitrary functions remains difficult. Existing methods generally fall into three categories: adaptive sampling (e.g., Tensor Cross-Interpolation or TCI), randomized projection, and methods exploiting explicit algebraic structure. The authors identify a gap for a deterministic, structure-preserving method that leverages classical polynomial approximation theory to construct QTT representations, particularly for analytic and highly differentiable functions, and extends efficiently to multivariate settings through function composition.

Methodology
The proposed approach constructs QTT approximants by expanding functions in orthogonal polynomial bases (specifically Chebyshev polynomials in the numerical experiments) and evaluating these expansions directly within the QTT algebra.

Finite-Precision QTT Algebra: The algorithm relies on a closed algebra of QTT operations, including:
- Truncation: A variational two-site scheme that reduces the rank of a QTT to a prescribed tolerance $\epsilon$ while maintaining orthonormality.
- Hadamard Product: Element-wise multiplication of two QTTs, which increases ranks and necessitates subsequent truncation.
- Outer Product: Concatenation of QTT cores for serial ordering or embedding via identity tensors for interleaved ordering.
Orthogonal Polynomial Expansion: The method approximates a function $f(x)$ as a truncated series $f_d(x) = \sum c_k P_k(x)$ . To avoid the numerical instability and rank explosion associated with explicitly constructing basis polynomials $P_k(x)$ via recurrence, the authors employ a Clenshaw-type recurrence. This allows the evaluation of the polynomial sum directly on the QTT argument $G$ (which encodes the grid or a multivariate function) without intermediate full-rank tensor construction.
Function Composition: The framework naturally handles composite functions $f(g(x))$ . The argument $g(x)$ is first represented as a QTT tensor (potentially constructed from lower-dimensional blocks via outer products). The polynomial expansion of $f$ is then applied to this QTT argument. This allows for the approximation of complex multivariate functions with compositional algebraic structure (e.g., $f(\sum x_i^2)$ ) without constructing the full high-dimensional coefficient tensor, which would suffer from the curse of dimensionality.

Key Contributions

Constructive Algorithm: The paper introduces a deterministic algorithm that combines orthogonal polynomial expansions with Clenshaw recurrences and variational TT compression to generate QTT representations.
Generalization to Multivariate Settings: Unlike methods that directly approximate full high-dimensional tensors, this approach exploits compositional structure to scale efficiently to multivariate functions, handling up to 200 QTT cores (corresponding to $2^{200}$ grid points).
Numerical Stability: The use of orthonormal QTT forms during truncation and the Clenshaw recurrence ensures numerical stability even at high polynomial orders and for highly differentiable functions.
Benchmarking: The work provides a systematic comparison against Tensor Cross-Interpolation (TCI) and multiscale interpolative constructions, analyzing error, rank, and computational time across univariate and multivariate benchmarks.

Results

Univariate Performance: For smooth, analytic functions (e.g., Gaussian, oscillatory functions), the method exhibits exponential convergence consistent with classical Chebyshev theory. The computational cost scales linearly with the number of QTT cores ( $O(n)$ ) and sub-quadratically with polynomial degree ( $O(d^{1.48})$ ), significantly better than worst-case estimates.
Comparison with TCI: In univariate settings, the method does not universally outperform TCI, particularly for non-smooth or discontinuous functions where TCI's adaptive sampling excels. However, for smooth functions, the proposed method requires orders of magnitude fewer function evaluations (approximately two orders of magnitude in benchmarks) while maintaining competitive accuracy and wall-clock time.
Multivariate Performance: The method successfully approximates multivariate Gaussian-like functions and nonlinear quartic models. It demonstrates that functions with compositional structure can be represented with sub-exponential complexity. The choice of QTT core ordering (serial vs. interleaved) significantly impacts compressibility; for example, "squeezed" Gaussian functions remain efficiently representable in both orderings, while product Gaussians show rank growth sensitivity to ordering.
Scalability: The algorithm scales to dimensions where direct dense representations are infeasible (e.g., $10^{60}$ tensor elements), with computational complexity growing polynomially with dimension for low-rank scenarios.

Significance and Claims
The authors position this work as a general framework for function composition in QTT form that generalizes efficiently to multidimensional settings. They claim that while the method is not a universal replacement for sampling-based techniques like TCI (especially for non-smooth targets), it offers a deterministic, structure-preserving alternative grounded in classical approximation theory.

The primary significance lies in its ability to:

Provide accurate approximations for highly differentiable functions with minimal sampling requirements.
Enable the efficient manipulation of composite functions within the QTT format, facilitating the solution of nonlinear problems in quantum-inspired numerical analysis.
Offer a tunable trade-off between accuracy and computational cost via the truncation tolerance parameter.

The paper concludes that orthogonal polynomial expansions provide a robust mechanism for loading smooth functions into QTT formats, with favorable asymptotic scaling that complements existing methods, particularly in regimes where high density and smoothness are present. The implementation is made available via the open-source SeeMPS library.

Technical Summary: Approximation and Composition of Functions in Quantized Tensor Trains via Orthogonal Polynomial Expansions

More like this