Numerical Analysis of differential equations on weighted Sobolev spaces: beyond classical orthogonal polynomials

Technical Summary: Numerical Analysis of Differential Equations on Weighted Sobolev Spaces

Problem Statement
The paper addresses the numerical resolution of differential equations of the form $Lu = f(u)$, where $L = V'(x)\partial_x - \partial_{xx}$, posed on weighted Sobolev spaces $H^1(\nu)$ associated with a Gibbs measure $\nu(dx) = e^{-V(x)}dx/Z$. While the functional analysis of these spaces is well-established for studying Schrödinger-type equations and diffusion processes, their numerical treatment has historically been limited to classical orthogonal polynomial bases (e.g., Hermite, Laguerre, Jacobi). These classical bases are preferred because the derivative of a basis function remains within the same family, allowing for a simple (often diagonal) representation of the differentiation operator. However, for non-classical potentials $V$ (specifically those growing faster than quadratically, such as quartic or sextic potentials), the associated orthogonal polynomials do not possess this property, rendering standard spectral methods ineffective or difficult to implement.

Methodology
The authors propose a novel framework to bridge this gap by introducing a new class of Sobolev orthogonal polynomials tailored to the specific structure of the weighted Sobolev space $H^1(\nu)$.

1. New Inner Product: Instead of the standard Sobolev inner product (which involves a parameter $\Lambda$ and is difficult to construct), the authors define a natural inner product on $H^1(\nu)$:
   $$(u, v) := \int_{\mathbb{R}} u'v' d\nu + \left(\int_{\mathbb{R}} u d\nu\right)\left(\int_{\mathbb{R}} v d\nu\right)$$
   This product is well-defined provided a Poincaré inequality holds for the potential $V$. The associated orthonormal basis $Q = \{q_n\}$ is constructed simply by integrating the standard $L^2(\nu)$-orthonormal polynomials $P = \{p_n\}$ and normalizing, such that $q_n' = p_{n-1}$ and $\int q_n d\nu = 0$.

2. Operator Factorization: The core theoretical insight is the factorization of the operator $L$ with respect to these bases. Let $D_{>0}$ be the upper-triangular matrix representation of the differentiation operator $\partial_x$ restricted to mean-zero functions in the $L^2(\nu)$ basis. The authors demonstrate that:

   • In the $L^2(\nu)$ basis $P$, $L$ admits a Cholesky-like factorization: $[L]_{P^* \to P^*} = D_{>0}^T D_{>0}$.
   • In the Sobolev basis $Q$, $L$ admits a reverse Cholesky factorization: $[L]_{Q^* \to Q^*} = D_{>0} D_{>0}^T$.
   • Crucially, for polynomial potentials of degree $2k$, the operator $D_{>0}$ is banded (specifically, it has $k$ non-zero superdiagonals for a quartic potential), rather than diagonal.

3. Quantitative Compactness Estimates: To justify the numerical convergence, the authors derive explicit bounds on the compactness of the embedding $H^1(\nu) \hookrightarrow L^2(\nu)$. This relies on analyzing the growth of the recurrence coefficients $(a_n)$ of the polynomials $P$. For a quartic potential, these coefficients satisfy a discrete Painlevé I equation. The authors use computer-assisted proofs to rigorously bound the growth of these coefficients, establishing that the embedding decays at a rate of $O(n^{-3/4})$ for quartic potentials and $O(n^{-(2k-1)/2k})$ for general even polynomial potentials of degree $2k$.

4. Computer-Assisted Proofs (CAPs): The theoretical estimates are integrated into a rigorous numerical framework. The authors construct approximate solutions and use the derived compactness constants to prove the existence of true solutions in a certified neighborhood, effectively enclosing the error bounds.

Key Results
The paper provides rigorous mathematical foundations and numerical results for two specific applications:

• Gross–Pitaevskii Equation with Sextic Potential: The authors solve the stationary Gross–Pitaevskii equation with a sextic potential ($W(x) = x^6/4 - \kappa x^4/2 + \dots$). They compute approximate solutions and prove the existence of true solutions with extremely tight error bounds (e.g., $\|\phi^\star - \bar{\phi}\|_{H^1} \leq 10^{-100}$). This demonstrates that the proposed Sobolev basis outperforms classical bases (like Hermite polynomials) for non-quadratic potentials.
• Stochastic Resonance in the Duffing Oscillator: The authors rigorously demonstrate the phenomenon of stochastic resonance for a Brownian particle in a double-well potential subject to periodic forcing and noise. By solving the time-periodic Fokker–Planck equation on the unbounded real line, they compute the stationary periodic density with explicit error bounds. They prove the existence of a specific noise strength $\sigma_{SR}$ that maximizes the resonance indicator, providing the first rigorous mathematical validation of this phenomenon for the stochastic Duffing oscillator on an unbounded domain.

Significance and Claims
The paper claims to lay the mathematical foundations for the numerical analysis of differential equations on weighted Sobolev spaces beyond classical orthogonal polynomials. Its primary contributions are:

1. A Tractable Basis: It introduces a "natural" notion of Sobolev orthogonal polynomials that are simpler to construct and more tractable than traditional Sobolev orthogonal polynomials, which often depend non-trivially on arbitrary parameters.
2. Factorization Strategy: It shows that even when $L$ is not diagonal, the specific factorization $L = D^T D$ (or $D D^T$) allows for efficient inversion and numerical solution.
3. Explicit Compactness: It provides the first derivation of fully explicit compactness estimates for non-classical Sobolev orthogonal polynomial bases, linking the rate of decay to the growth of recurrence coefficients (solutions to discrete Painlevé equations).
4. Rigorous Validation: It successfully applies this framework to complex, non-linear problems on unbounded domains, achieving error bounds orders of magnitude smaller than Monte Carlo simulations and providing the first rigorous proof of stochastic resonance for the Duffing oscillator in this setting.

The authors emphasize that while the specific quantitative bounds were derived for quartic and sextic potentials, the underlying methodology (factorization via $D_{>0}$ and the analysis of recurrence coefficients) generalizes to any even polynomial potential. The work validates the use of non-classical bases for high-precision numerical analysis where classical bases fail.