IRSP: Inverse Random Source Problem Overview

• IRSP is a class of inverse problems where the statistical descriptors (mean, variance, and covariance) of random forcing terms in PDEs are recovered from field observations.
• Methodologies include modeling the source as an additive mean-plus-noise term or as a microlocally isotropic Gaussian field, enabling recovery through Fourier analysis, Born approximations, and resolvent techniques.
• Key computational strategies employ regularization, phase retrieval, and multi-frequency data assimilation to mitigate ill-posedness and enhance reconstruction stability.

Inverse random source problem (IRSP) denotes a class of inverse problems in which the forcing term in a PDE is modeled as a random field, and the objective is to recover statistical descriptors of that source from observations of the induced field. In the cited arXiv literature, the unknown may be the mean and variance of an additive random source, the micro-correlation strength in the principal symbol of a covariance operator, a matrix-valued micro-correlation strength, or the principal symbols of both covariance and relation operators. IRSPs have been formulated for time-harmonic Helmholtz and Maxwell equations, elastic and biharmonic wave equations, stochastic wave and diffusion equations, the multi-term time-fractional diffusion-wave equation, the Helium production-diffusion equation, the biharmonic Schrödinger equation, and the fractional Helmholtz equation (Li et al., 2020, Wang et al., 2024, Li et al., 2021, Feng et al., 2019, Feng et al., 2021, Li et al., 2022, Feng et al., 2023, Wang et al., 2024, Li et al., 23 Feb 2026).

1. Canonical formulations and source models

A recurrent feature of IRSP is that the forward field is deterministic only conditionally on a realization of the source; the inverse target is therefore statistical rather than pathwise. In the cited works, one broad formulation writes the source as an additive mean-plus-noise term, for example $f(x)=g(x)+h(x)\,\dot W_x$ in stochastic Helmholtz models, or $f(x)h(t)+g(x)\,\dot B^H(t)$ in stochastic wave and fractional diffusion models. The corresponding inverse problem asks for $g$, $h^2$, $f$, or $|g|$ from expectations, variances, covariances, or phaseless Fourier data (Li et al., 2023, Chen et al., 29 Aug 2025, Feng et al., 2019, Feng et al., 2021, Niu et al., 2018).

A second formulation models the source as a microlocally isotropic generalized Gaussian random field whose covariance operator is a classical pseudo-differential operator. For the stochastic Helmholtz equation in two and three dimensions, the principal symbol is $c(x,\xi)=h(x)|\xi|^{-m}$, where $h\in C^\infty_0(D)$ is the micro-correlation strength (Wang et al., 2024). For the three-dimensional time-harmonic Maxwell system, the source is a centered complex-valued Gaussian vector field whose covariance operator has principal symbol $M(x)|\xi|^{-2s}$, $0$f(x)h(t)+g(x)\,\dot B^H(t)$0 called the micro-correlation strength matrix (Li et al., 2020).

A third formulation, developed for the fractional Helmholtz equation, introduces both covariance and relation operators. In that setting the source is a centered complex-valued microlocally isotropic generalized Gaussian random field of order $f(x)h(t)+g(x)\,\dot B^H(t)$1, with

$f(x)h(t)+g(x)\,\dot B^H(t)$2

so the inverse problem targets both $f(x)h(t)+g(x)\,\dot B^H(t)$3 and $f(x)h(t)+g(x)\,\dot B^H(t)$4 (Li et al., 23 Feb 2026).

Formulation	Random source model	Typical inverse target
Additive mean-plus-noise	$f(x)h(t)+g(x)\,\dot B^H(t)$5 or $f(x)h(t)+g(x)\,\dot B^H(t)$6	Mean, variance, $f(x)h(t)+g(x)\,\dot B^H(t)$7, $f(x)h(t)+g(x)\,\dot B^H(t)$8, phaseless Fourier modes
Microlocally isotropic Gaussian	Covariance principal symbol $f(x)h(t)+g(x)\,\dot B^H(t)$9	Micro-correlation strength $g$0
Complex GMIG / vector Gaussian	Covariance or relation symbols, or matrix symbol $g$1	$g$2, $g$3, or matrix-valued $g$4

This taxonomy clarifies a common misconception: IRSP generally does not attempt to reconstruct a single sample of a random source. In the cited literature, the recoverable object is typically a statistical descriptor encoded in expectation, covariance, variance, or principal-symbol data.

2. Direct problems and forward maps

The direct problem in IRSP is the stochastic PDE with a rough source, and much of the literature first establishes well-posedness in mild, distributional, or weighted Sobolev senses. For the stochastic Helmholtz equation, if $g$5, then almost surely there is a unique solution $g$6, $g$7, $g$8, with the pointwise representation

$g$9

and in particular $h^2$0 (Wang et al., 2024).

For Maxwell’s equations, the time-harmonic system

$h^2$1

can be reduced, in the distributional sense, to

$h^2$2

The unique solution is given almost surely by the volume potential

$h^2$3

and the field inherits Sobolev regularity from the rough source $h^2$4 (Li et al., 2020).

For the fractional Helmholtz equation, the direct problem is analyzed through the Lippmann–Schwinger equation

$h^2$5

with $h^2$6 the free resolvent and $h^2$7 the potential term. For sufficiently large wavenumbers, $h^2$8 is invertible and the stochastic fractional Helmholtz problem is well posed in the distributional sense (Li et al., 23 Feb 2026). Closely related resolvent techniques also appear for the biharmonic Schrödinger equation, where meromorphic continuation and resolvent bounds underpin the direct problem and Born series analysis (Wang et al., 2024).

In time-fractional diffusion and wave models, the direct problem is usually formulated in modal or semigroup form. The solutions are represented through Mittag–Leffler kernels, sine propagators, Green’s functions in the frequency domain, or analytic semigroups, depending on the model. These formulations yield existence, uniqueness, and regularity under conditions such as $h^2$9 or $f$0, and they provide the moment formulas needed for inversion (Feng et al., 2019, Feng et al., 2021, Niu et al., 2018, Li et al., 2022, Feng et al., 2023).

The direct theory is not auxiliary. It determines the admissible observation functionals, the regularity class of the data, and the analytic structure that later supports uniqueness, stability, and numerical inversion.

3. Observation models and identifiable quantities

The observation model varies substantially across IRSP formulations, and the identifiable quantity depends on that choice. In stochastic fractional diffusion and stochastic wave problems driven by fractional Brownian motion or white noise, the data are often the expectation and covariance of the final-time modal coefficients. Under positivity conditions on the temporal weight, the collection $f$1 uniquely determines the deterministic source profile and the magnitude of the stochastic amplitude, typically $f$2 and $f$3 (Feng et al., 2019, Feng et al., 2021, Niu et al., 2018, Li et al., 2022).

In the multi-term time-fractional diffusion-wave equation, the measured quantity is boundary data $f$4, or equivalently $f$5. The key identity is

$f$6

with $f$7, so the phaseless Fourier modes are recovered by

$f$8

The remaining phase retrieval step is then handled computationally (Feng et al., 2023).

For Maxwell’s equations with a Gaussian vector source, the first moment vanishes, $f$9, so the inverse problem uses the second moment

$|g|$0

For $|g|$1 in a bounded observation domain disjoint from the source support,

$|g|$2

and hence the matrix-valued high-frequency limit determines $|g|$3. The diagonal entries $|g|$4 are recoverable even from the phaseless amplitudes $|g|$5 (Li et al., 2020).

For the one-dimensional Helmholtz equation with attenuation, the high-frequency energy identity is

$|g|$6

so the micro-correlation strength $|g|$7 is encoded by convolution with $|g|$8 and can be recovered by Fourier inversion (Li et al., 2020).

Multi-frequency correlation data are central in stochastic Helmholtz scattering. Near-field correlations and far-field correlations over a finite frequency band determine the micro-correlation strength $|g|$9, with explicit estimates for $c(x,\xi)=h(x)|\xi|^{-m}$0 from the data and a high-frequency tail (Wang et al., 2024). Statistical phaseless data form a different observation model: by adding reference point sources and using algebraic phase-retrieval formulas, one first reconstructs $c(x,\xi)=h(x)|\xi|^{-m}$1, $c(x,\xi)=h(x)|\xi|^{-m}$2, $c(x,\xi)=h(x)|\xi|^{-m}$3, and $c(x,\xi)=h(x)|\xi|^{-m}$4, and then solves Fredholm integral equations for the source mean $c(x,\xi)=h(x)|\xi|^{-m}$5 and variance $c(x,\xi)=h(x)|\xi|^{-m}$6 (Chen et al., 29 Aug 2025).

Single-realization recovery is also possible in some high-frequency settings. For the biharmonic wave equation and for random elastic scattering, the principal symbol of the covariance operator is uniquely determined, with probability one, from a single realization of the magnitude of the wave field averaged over a frequency band (Li et al., 2021, Li et al., 2018).

4. Analytical mechanisms of uniqueness

The main uniqueness mechanisms in IRSP are microlocal asymptotics, harmonic analysis, and stochastic averaging. In microlocally isotropic models, high-frequency asymptotics reduce correlation data to a deterministic potential of the unknown strength. For Maxwell, stationary phase and microlocal analysis yield the Newtonian-type integral equation for $c(x,\xi)=h(x)|\xi|^{-m}$7, while for the biharmonic wave equation the frequency-averaged intensity converges to

$c(x,\xi)=h(x)|\xi|^{-m}$8

with

$c(x,\xi)=h(x)|\xi|^{-m}$9

Standard deconvolution then recovers $h\in C^\infty_0(D)$0 (Li et al., 2020, Li et al., 2021).

Born linearization and far-field asymptotics play a similar role in elastic, biharmonic Schrödinger, and fractional Helmholtz models. In the elastic case, the Born approximation isolates the leading source term, and averaging $h\in C^\infty_0(D)$1 over frequency extracts a single-layer potential of the principal symbol $h\in C^\infty_0(D)$2 (Li et al., 2018). For the biharmonic Schrödinger equation, the zeroth-order far-field correlation is the leading term and higher Born terms are shown to be lower order in $h\in C^\infty_0(D)$3, enabling stability from finite-band correlation far-field data (Wang et al., 2024). For the fractional Helmholtz equation, the far-field pattern admits a decomposition $h\in C^\infty_0(D)$4; correlating $h\in C^\infty_0(D)$5 and $h\in C^\infty_0(D)$6 over $h\in C^\infty_0(D)$7 recovers $h\in C^\infty_0(D)$8, while a related formula using $h\in C^\infty_0(D)$9 recovers $M(x)|\xi|^{-2s}$0 (Li et al., 23 Feb 2026).

A distinct mechanism is analytic continuation in the frequency variable. For the multi-frequency stochastic Helmholtz problem, the relevant boundary-integral quantities extend holomorphically to a sector in $M(x)|\xi|^{-2s}$1, and a Phragmén–Lindelöf-type lemma controls the high-frequency tail from data on a finite band. This is the basis of the near-field and far-field increasing-stability estimates (Wang et al., 2024). The biharmonic Schrödinger analysis uses the same general strategy after establishing a meromorphic continuation and bounds for the resolvent (Wang et al., 2024).

Fixed-frequency correlation inversion requires different tools. The 2026 acoustic work establishes conditional Hölder-logarithmic stability through complex geometrical optics solutions and derives a variational source condition for spectral regularization methods (Mickan et al., 24 Feb 2026). This shifts the IRSP analysis from purely high-frequency asymptotics to a fixed-frequency operator-theoretic framework.

Finally, some single-realization reconstructions depend on ergodicity in frequency. In the Maxwell problem, the frequency average

$M(x)|\xi|^{-2s}$2

converges almost surely to the ensemble quantity $M(x)|\xi|^{-2s}$3 under mild continuity assumptions and covariance decay in $M(x)|\xi|^{-2s}$4 (Li et al., 2020). An analogous frequency-domain ergodicity argument appears in the biharmonic wave problem (Li et al., 2021). This resolves another common misunderstanding: a “single realization” result in IRSP is usually a statement about an averaged frequency functional, not about inversion from one field sample at one frequency.

5. Stability, ill-posedness, and regularization theory

Uniqueness in IRSP rarely implies robust reconstruction. Many models are explicitly ill posed. In the multi-term time-fractional diffusion-wave equation, the inversion

$M(x)|\xi|^{-2s}$5

is ill posed of polynomial type because $M(x)|\xi|^{-2s}$6 decays at high frequency (Feng et al., 2023). In the time-fractional diffusion equation driven by a fractional Brownian motion, the inverse map is described as severely ill posed: the kernels $M(x)|\xi|^{-2s}$7 and $M(x)|\xi|^{-2s}$8 decay for large eigenvalues, so small data errors are strongly amplified (Feng et al., 2019). Closely related instability mechanisms appear in the stochastic wave equation, the Helium production-diffusion equation, and the stochastic fractional diffusion equation, where recovering high modes entails amplification by factors depending on $M(x)|\xi|^{-2s}$9 (Feng et al., 2021, Li et al., 2022, Niu et al., 2018).

Multi-frequency data can partially offset this ill-posedness. For the one-dimensional stochastic Helmholtz equation driven by white noise, the stability estimates contain a data discrepancy term plus a tail term of order $0Li et al., 2016). The stochastic Helmholtz stability theory in two and three dimensions makes this more explicit: the estimate consists of a Lipschitz type data discrepancy and a logarithmic stability term, and the latter decreases as the upper bound of the frequency increases, exhibiting increasing stability (Wang et al., 2024).

Some problems remain logarithmically unstable even with more structured measurements. For stochastic Maxwell equations in an inhomogeneous medium, a logarithmic stability estimate is established for the source strength $0Wang et al., 13 Mar 2025).

The 2026 acoustic correlation paper develops a more refined stability theory. Under Sobolev smoothness assumptions, it proves conditional Hölder-logarithmic stability and, via a variational source condition, Hölder-logarithmic convergence rates for spectral regularization methods. It also states that, at fixed frequency, the exponents in the logarithmic stability and convergence estimates grow unboundedly as the Sobolev regularity of the source increases (Mickan et al., 24 Feb 2026).

Taken together, these results show that IRSP stability theory is now stratified by data modality: final-time moment inversion is typically severely ill posed; finite-band multi-frequency scattering may exhibit increasing stability; and fixed-frequency correlation inversion often leads to logarithmic or Hölder-logarithmic estimates.

6. Computational methodologies and open directions

The computational literature mirrors the analytical diversity of IRSP. For Fredholm integral equations of the first kind arising from mean and variance data, regularized Kaczmarz methods are a standard choice. They are used for inverse random source scattering in elastic waves and also as the first stage of a data-assisted two-stage method for the Helmholtz equation (Bao et al., 2016, Li et al., 2023). Tikhonov regularization, spectral cut-off, and truncation-based regularization appear repeatedly in stochastic fractional diffusion, stochastic wave, and Helium production-diffusion models (Niu et al., 2018, Feng et al., 2021, Li et al., 2022).

Phaseless inverse problems require an additional phase-retrieval layer. For stochastic time-fractional diffusion and for the multi-term time-fractional diffusion-wave equation, PhaseLift with random masks is used after reconstructing Fourier magnitudes from variance formulas. The convex lifted formulation

$0

provides a computational route from phaseless Fourier data to the unknown temporal source profile (Gong et al., 2020, Feng et al., 2023).

More recent work incorporates statistical learning. The two-stage Helmholtz method first generates initial reconstructions by regularized Kaczmarz and then treats enhancement as an image-to-image translation problem, comparing PCA plus linear regression, dynamic mode decomposition, U-Net, and pix2pix (Li et al., 2023). Bayesian inversion has also entered the field: for statistical phaseless data in the two-dimensional Helmholtz equation, Gaussian-prior Bayesian formulations are used to reconstruct $0Chen et al., 29 Aug 2025).

A notable computational trend is the reduction of frequency requirements. A single-frequency, non-iterative method for elastic waves reconstructs the variance from boundary correlation data by solving Fourier-domain $0Gu et al., 1 Nov 2025). This stands in contrast to the earlier emphasis on multi-frequency iteration.

Reported open problems remain substantial. In the Maxwell setting, explicit questions include removal of the assumption $0Li et al., 2020). In the multi-term time-fractional diffusion-wave equation, higher-dimensional domains, spatially distributed fractional noise, inverse random potential problems, and refined regularization strategies are listed as extensions (Feng et al., 2023). Inhomogeneous-medium problems suggest further work on multi-frequency data, additional polarization states, and higher statistical moments (Wang et al., 13 Mar 2025). The data-assisted literature adds joint reconstruction of the medium and source statistics, operator-learning alternatives, and uncertainty quantification for learned reconstructions (Li et al., 2023).

These directions indicate that IRSP has evolved from uniqueness theory for idealized data toward a broader program that combines microlocal analysis, stochastic PDE theory, stability estimates, regularization theory, and computation under finite-band, phaseless, or single-realization measurements.