Direct Iteration Inversion Methods

• Inversion by Direct Iteration is a family of methods that decompose inverse problems into sequences of tractable subproblems solved via explicit update maps.
• These techniques span diverse applications, from nonlinear PDE eigenvalue problems and non-uniform difference operators to stochastic simulation and image restoration.
• The methods prioritize iterative accumulation of inverse information, achieving high accuracy and efficiency by replacing hard inversion with structured, repeatable computations.

Searching arXiv for relevant papers on “Inversion by Direct Iteration” and closely related direct inversion formulations. In the cited arXiv literature, “inversion by direct iteration” denotes a family of inversion strategies in which the inverse action is realized by an explicitly specified update map applied repeatedly to the object of interest, rather than by forming a closed-form inverse or by embedding a new forward solve inside every refinement step. The phrase covers several distinct but related practices: inverse iteration for the Monge–Ampère eigenvalue problem, orbit-sum inversion of non-uniform shift operators, diffusion-based realization of $(I-P)^{-1}$, direct inverse-CDF sampling for stochastic processes, one-shot inverse NFFT operators obtained after precomputation, structural inversion of pair potentials from fixed particle configurations, and supervised image restoration by learned multi-step refinement (Abedin et al., 2020, Temple et al., 2018, Hong, 2012, Kircheis et al., 2018, Kircheis et al., 2023, Malham et al., 2024, Malham et al., 2020, Coquand et al., 12 Mar 2026, Delbracio et al., 2023).

1. Terminological scope

The adjective “direct” does not have a single meaning across these works. In stochastic simulation, it refers to evaluating an explicit approximation of an inverse distribution function, so that there is no per-sample Newton iteration to solve $F(x)=u$. In inverse NFFT, it means that after a one-time precomputation, each inversion costs about as much as one adjoint NFFT. In direct Boltzmann inversion, it means that the same fixed ensemble of configurations is reused and that no new Monte Carlo or molecular-dynamics simulation is run at each update. In supervised image restoration, it means that the restoration process is learned directly from paired low-quality and high-quality examples, without requiring knowledge of any analytic form of the degradation process, and that inference starts from the degraded observation rather than from pure noise (Malham et al., 2024, Malham et al., 2020, Kircheis et al., 2023, Coquand et al., 12 Mar 2026, Delbracio et al., 2023).

The noun “iteration” is equally heterogeneous. It may denote repeated inversion of a nonlinear PDE operator with Dirichlet boundary condition, repeated composition with a non-uniform shift map $\Phi$, repeated diffusion of residual “fluid” through a graph, or repeated application of a learned denoising-style update. What remains common is that the inverse problem is decomposed into a sequence of smaller or more structured subproblems, each of which is simpler than the original global inversion (Abedin et al., 2020, Temple et al., 2018, Hong, 2012, Delbracio et al., 2023).

2. Nonlinear operator inversion in analysis

A canonical analytical instance is the Monge–Ampère eigenvalue problem on a bounded, convex domain $\Omega\subset\mathbb{R}^n$,

$$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$

with $u$ convex and nontrivial. Because the equation is homogeneous, the eigenfunction is unique only up to positive multiplicative constants. The inverse-iteration scheme defines $u_{k+1}$ by

$$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$

where

$$R(u) := \frac{\displaystyle\int_{\Omega} (-u)\, dMu}{\displaystyle\int_{\Omega} (-u)^{n+1}}.$$

The normalized iterates $\hat u_k=u_k/\|u_k\|_{L^\infty(\Omega)}$ converge uniformly to the unique eigenfunction of unit height, and $F(x)=u$0. The method is explicitly presented as a nonlinear analogue of inverse iteration for matrices, with convergence established in the uniform topology for convex functions and with an explicit admissible initial class, including the paraboloid $F(x)=u$1 restricted to $F(x)=u$2 (Abedin et al., 2020).

A different analytical use of direct iteration appears for the non-uniform difference operator

$F(x)=u$3

Here Fourier diagonalization fails because the shift is non-uniform. The inversion is carried out in physical space by iterating the functional equation

$F(x)=u$4

Forward and backward compositions yield

$F(x)=u$5

and, after passing to the limit along orbits converging to the fixed points $F(x)=u$6 and $F(x)=u$7,

$F(x)=u$8

Under the solvability condition

$F(x)=u$9

the operator has a bounded inverse on its range in $\Phi$0, uniquely determined up to an additive constant. In this setting, direct iteration means explicit summation along the discrete characteristics generated by $\Phi$1, and it is introduced precisely because non-uniform periodicity destroys the translation invariance needed for Fourier methods (Temple et al., 2018).

3. Direct inversion in stochastic simulation

For the squared Bessel process $\Phi$2 of dimension $\Phi$3,

$\Phi$4

the exact transition law satisfies

$\Phi$5

The direct inversion method samples this transition by inverting the central chi-square CDF through a two-dimensional Chebyshev expansion in $\Phi$6, while the non-central case is reduced to the central one by standard Poisson-mixture decompositions. The $\Phi$7-axis is split into three regions, region-specific nonlinear reparametrizations $\Phi$8 are used, and evaluation is performed by two nested Chebyshev sums via Clenshaw’s recurrence. The method is explicitly non-iterative at sampling time: there is no per-sample root-finding loop. Numerical results state that the inverse-CDF approximation can achieve accuracy on the order of $\Phi$9, and that to reach approximately $\Omega\subset\mathbb{R}^n$0 relative error in option prices, direct inversion can be markedly faster than QE-M and full truncation Euler, with one reported CPU-time comparison of $\Omega\subset\mathbb{R}^n$1 versus $\Omega\subset\mathbb{R}^n$2 versus $\Omega\subset\mathbb{R}^n$3 for a CIR-priced put option (Malham et al., 2024).

An analogous but more elaborate inversion architecture is developed for the Heston model. The time-integrated variance conditional on its endpoints is represented, after a change of measure to a squared Bessel setting, by a random series

$\Omega\subset\mathbb{R}^n$4

The components are built from infinite weighted sums of parameter-free variables such as

$\Omega\subset\mathbb{R}^n$5

and the corresponding inverse CDFs are approximated offline by piecewise Chebyshev polynomials on carefully chosen probability regimes. Truncation errors decay exponentially in the truncation level $\Omega\subset\mathbb{R}^n$6, gamma tails are used for the outer-series remainder, and the resulting direct inversion routines are reused under any market conditions because the building-block distributions are independent of the Heston parameters. The final target law is obtained from the auxiliary law by an acceptance–rejection step with acceptance probability $\Omega\subset\mathbb{R}^n$7 (Malham et al., 2020).

4. Linear algebra, graph diffusion, and Fourier inversion

In linear systems, a direct-iteration viewpoint appears in the D-iteration method. Starting from a fixed-point formulation

$\Omega\subset\mathbb{R}^n$8

the method maintains a fluid vector $\Omega\subset\mathbb{R}^n$9 and a history vector $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$0. At each step, all fluid at one chosen index is removed, accumulated in $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$1, and diffused through the outgoing links encoded in the corresponding column of $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$2. The scheme is presented as a Gauss–Seidel–like method recast in diffusion terms, with a stronger order-independence and a graph-based interpretation. For PageRank-type problems,

$$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$3

the method realizes the inverse action of $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$4 without forming it explicitly, and $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$5 provides the exact distance to the limit when $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$6 has no zero column, otherwise an upper bound (Hong, 2012).

For nonequispaced Fourier inversion, “direct” means something different: the inverse operator is precomputed and then applied in one NFFT-like pass. In the square case $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$7, a direct method of complexity $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$8 uses Lagrange interpolation together with fast summation and a final FFT. In the rectangular cases, the matrix representation $$\begin{cases} \det D^2 u = \lambda_{MA} (-u)^n & \text{in } \Omega,\ u=0 & \text{on } \partial\Omega, \end{cases}$$9 is modified so that a weighted adjoint NFFT or a modified adjoint NFFT approximates an inverse. One route chooses density-compensation weights $u$0 so that

$u$1

and, when the double-bandwidth system

$u$2

is satisfied, reconstruction is exact for trigonometric polynomials of degree $u$3. A second route optimizes the sparse NFFT matrix in Frobenius norm so that

$u$4

In both cases, the precomputation may itself use conjugate-gradient iterations or local least-squares solves, but each subsequent inversion requires only one adjoint NFFT-like application, with cost $u$5 (Kircheis et al., 2018, Kircheis et al., 2023).

5. Data-driven and statistical-mechanical inversion

In supervised image restoration, “Inversion by Direct Iteration” is the name of a specific learned framework. Given paired data $u$6, where $u$7 is a high-quality target image and $u$8 is the degraded observation, the method introduces the interpolation path

$u$9

A network $u_{k+1}$0 is trained to predict the clean image $u_{k+1}$1, and inference proceeds by the deterministic update

$u_{k+1}$2

optionally with an added noise term under the Brownian-noise formulation. The method is explicitly positioned against single-step regression, which under an $u_{k+1}$3 loss yields the conditional expectation $u_{k+1}$4 and therefore suffers from “regression to the mean.” It is also distinguished from conditional denoising diffusion, because it does not require a pre-specified forward diffusion/noising process and because inference starts from the degraded image rather than pure noise. Reported applications include motion and out-of-focus deblurring, super-resolution, compression artifact removal, and denoising (Delbracio et al., 2023).

In statistical mechanics, direct Boltzmann inversion addresses the inverse Henderson problem: reconstructing an isotropic pair interaction potential $u_{k+1}$5 from particle configurations at a fixed state point $u_{k+1}$6. The method uses two estimators of the radial distribution function $u_{k+1}$7: a distance-histogram estimator $u_{k+1}$8 and a force-based estimator

$u_{k+1}$9

For a trial potential $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$0, pairwise forces generate a trial RDF $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$1, and the potential is updated by the Schommers/IBI rule

$$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$2

or, when $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$3 at small $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$4, by the shifted variant

$$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$5

The procedure is called direct because it reuses a fixed dataset of configurations and does not perform a new simulation at each iteration. The paper emphasizes that the force formula remains valid at high density as long as forces are well-defined, and that the approach naturally applies to any state point, including when the density is large and alternative methods may fail (Coquand et al., 12 Mar 2026).

6. Common structure, misconceptions, and limitations

A recurring misconception is that “direct” means “non-iterative.” The cited literature shows the opposite. InDI in image restoration is explicitly a deterministic multi-step inversion process. Direct Boltzmann inversion is iterative in function space. The non-uniform difference-operator construction is an infinite orbit sum. D-iteration is a sequential residual-diffusion process. Inverse NFFT can be direct only after an iterative or optimization-based precomputation. By contrast, the stochastic-simulation papers reserve “direct inversion” mainly for the absence of per-sample nonlinear root-finding (Delbracio et al., 2023, Coquand et al., 12 Mar 2026, Temple et al., 2018, Hong, 2012, Kircheis et al., 2023, Malham et al., 2024).

A second misconception is that these methods all compute an explicit inverse operator. Several do not. D-iteration computes $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$6 or $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$7 only through the limit of the diffusion process. The squared Bessel and Heston methods compute inverse CDF actions rather than analytic inverses. The Monge–Ampère method recovers an eigenpair by repeated PDE solves, not by an explicit nonlinear inverse formula. The image-restoration method learns a residual flow, not a closed-form inverse map. This suggests a broader interpretation: inversion by direct iteration is often an inverse action, not an inverse object.

The principal limitations are domain-specific and structural. The Monge–Ampère scheme relies on convexity, Aleksandrov solutions, comparison principles, and the special Rayleigh quotient $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$8. The non-uniform difference-operator theory depends on the non-degeneracy of $$\begin{cases} \det D^2 u_{k+1} = R(u_k)(-u_k)^n & \text{in } \Omega,\ u_{k+1}=0 & \text{on } \partial\Omega, \end{cases}$$9, geometric convergence of $$R(u) := \frac{\displaystyle\int_{\Omega} (-u)\, dMu}{\displaystyle\int_{\Omega} (-u)^{n+1}}.$$0 to fixed points, and the solvability condition on orbit sums. The stochastic-simulation methods depend on accurate regime splitting, stable Chebyshev approximations, and manageable acceptance–rejection overhead. D-iteration requires contractive or diffusion-reducing structure, such as strict diagonal dominance by columns. Direct inverse NFFT depends on node geometry, rank conditions, and the quality of the precomputed weights or sparse matrices. InDI requires paired training data and can deteriorate for very large numbers of inference steps. Direct Boltzmann inversion yields a state-point-dependent effective pair potential when many-body effects are folded into a pairwise description. These constraints are not incidental; they are the mechanisms that make direct iteration analytically or computationally viable (Abedin et al., 2020, Temple et al., 2018, Malham et al., 2020, Hong, 2012, Kircheis et al., 2023, Delbracio et al., 2023, Coquand et al., 12 Mar 2026).

Across these settings, a plausible unifying interpretation is that inversion by direct iteration replaces a hard inverse problem by a structured recurrence whose elementary step is either exactly solvable, cheaply evaluable, or learnable from data. The resulting methods differ radically in their mathematics—Aleksandrov theory, orbit dynamics, Chebyshev spectral approximation, graph diffusion, frame-based Fourier inversion, and supervised residual flows—but they share the same operational idea: inverse information is accumulated progressively, with each step constrained by the geometry of the underlying problem.