🔢 mathematics

Approximation of the Lévy-driven stochastic heat equation on the sphere

This paper establishes new regularity results and derives strong and weak convergence rates for a spectral approximation in space combined with Euler–Maruyama schemes in time to solve the Lévy-driven stochastic heat equation on the sphere, with theoretical findings validated by numerical simulations.

Original authors: Annika Lang, Andrea Papini, Verena Schwarz

Published 2026-06-04

📖 4 min read🧠 Deep dive

CC BY 4.0

Original authors: Annika Lang, Andrea Papini, Verena Schwarz

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

Imagine the Earth as a giant, glowing orange. Now, imagine that the surface of this orange isn't just sitting still; it's being heated up and cooled down in a chaotic, unpredictable way. Sometimes the heat spreads smoothly like warm water, but other times, it gets hit by sudden, sharp jolts—like someone throwing ice cubes or hot coals onto it randomly.

This is the problem the authors of this paper are trying to solve: How do we mathematically predict the temperature on a sphere (like Earth) when the heat source is a mix of smooth changes and sudden, random jumps?

Here is a breakdown of their work using simple analogies:

1. The Problem: Smooth vs. Jumpy

For a long time, scientists have been good at modeling "smooth" randomness (like the gentle, rolling waves of the ocean, known as Wiener noise). They use these models for things like weather or finance.

However, real life is often "jumpy." A sudden storm, a stock market crash, or a biological mutation happens all at once. The authors wanted to study a sphere where the randomness comes from Lévy processes. Think of a Lévy process as a "drunkard's walk" that occasionally takes giant, unpredictable leaps instead of just stumbling forward.

2. The Solution: Breaking the Orange into Pieces

To solve the equation for this "jumpy" heat on a sphere, the authors used a method called Spectral Approximation.

The Analogy: Imagine you want to describe the shape of a bumpy orange. Instead of trying to describe every single bump individually, you break the orange down into a set of standard, smooth "bump patterns" (mathematicians call these Spherical Harmonics).
The Method: You start with just the biggest, smoothest patterns. Then, you add smaller and smaller patterns to get more detail. The authors proved that if you keep adding these patterns, your approximation gets closer and closer to the true, chaotic reality. They calculated exactly how fast this happens based on how "rough" or "smooth" the random jumps are.

3. The Time Machine: Forward and Backward Steps

Calculating the heat at every single instant is impossible, so you have to take steps through time. The authors tested two ways to take these steps:

Forward Euler-Maruyama: This is like looking at where you are now and guessing where you'll be in the next second based on your current speed. It's simple, but it can get unstable if your steps are too big (like trying to walk a tightrope without a pole).
Backward Euler-Maruyama: This is like looking at where you want to be in the next second and working backward to figure out the step needed to get there. It's more stable and doesn't care how big your steps are, making it a safer bet for these chaotic equations.

4. The Results: Strong vs. Weak Predictions

The paper distinguishes between two types of accuracy:

Strong Convergence (The "Exact Path"): This asks, "If I run the simulation and the real world, how close are the two specific paths?" The authors found that the accuracy depends heavily on how rough the "jumps" are. If the jumps are very smooth, the math works great. If they are very jagged, you need more patterns to get a good answer.
Weak Convergence (The "Average Outcome"): This asks, "If I run the simulation 1,000 times, does the average result match the real average?"
- The Magic Trick: The authors discovered a "double bonus." While the "exact path" might be off by a certain amount, the "average result" is often twice as accurate. This is a huge win for scientists who care more about general trends than specific, chaotic details.

5. The Proof: Computer Simulations

Finally, the authors didn't just do the math on paper; they ran computer simulations. They created digital versions of the sphere with different types of "noise" (smooth waves, sudden jumps, or a mix of both).

They watched the computer models converge (get closer to the truth) as they added more patterns and took smaller time steps. The results matched their mathematical predictions perfectly, confirming that their formulas work for both smooth and jumpy scenarios.

Summary

In short, this paper provides a new, reliable toolkit for predicting how things change on a sphere (like our planet) when the world is full of sudden, unpredictable shocks. They proved that by breaking the problem into smooth building blocks and using stable time-stepping methods, we can accurately predict both the specific path of the chaos and the average behavior of the system, with the average being surprisingly easy to predict.

Technical Summary: Approximation of the Lévy-Driven Stochastic Heat Equation on the Sphere

Problem Statement
This work addresses the numerical approximation of the stochastic heat equation (SHE) on the unit sphere $S^2$ , driven by an additive, square-integrable Lévy process. While stochastic partial differential equations (SPDEs) on Euclidean spaces driven by Gaussian (Wiener) noise are well-established, real-world phenomena such as turbulent flows and stellar fluids often exhibit abrupt, non-Gaussian behavior. Furthermore, applications on the Earth necessitate extending SPDE theory to spherical domains. This paper bridges the gap between classical Gaussian settings and more general non-Gaussian settings on surfaces, specifically focusing on the SHE:
$dX(t) = \Delta_{S^2}X(t) dt + dL(t)$
where $L$ is an infinite-dimensional Lévy process independent of the initial condition $X_0$ .

Methodology
The authors employ a semi-discrete spectral method in space and fully discrete Euler–Maruyama schemes in time.

Functional Setting: The solution is expanded using real-valued spherical harmonic functions $Y_{\ell,m}$ , which form an orthonormal basis for $L^2(S^2)$ . The driving noise $L$ is similarly expanded into a series of independent (or identically distributed) Lévy processes $L_{\ell,m}$ with coefficients determined by the noise's regularity.
Spatial Approximation: A spectral truncation $X^{(\kappa)}$ is defined by summing the series expansion up to a maximum degree $\kappa$ .
Temporal Approximation: To handle the stochastic convolution which is difficult to simulate for general Lévy noise, the authors introduce forward and backward Euler–Maruyama schemes. The backward scheme is noted for its unconditional stability, whereas the forward scheme requires a coupling between the spatial truncation $\kappa$ and the time step $h$ (specifically $\ell(\ell+1)h \leq C_c$ ) to ensure stability.
Regularity Analysis: New regularity results are derived for the mild solution in Sobolev spaces $H^\rho(S^2)$ and $L^p(\Omega, H^\rho(S^2))$ . These results distinguish between settings where the Lévy noise components are independent versus more general settings involving martingale decompositions (Wiener plus compensated Poisson parts).

Key Contributions and Results

Strong Convergence Rates:
- For the spectral approximation, the strong error $\|X - X^{(\kappa)}\|_{L^2(\Omega, L^2(S^2))}$ converges with rate $\kappa^{-\beta}$ .
- Under general Assumption 2.1 (noise in $H^\eta$ ), $\beta = \eta + 1$ .
- Under the specific isotropic Assumption 2.2 (angular power spectrum decay $\ell^{-\alpha}$ ), $\beta = \alpha/2$ .
- For the time discretization (Euler–Maruyama), the strong error converges with rate $h^{\min\{\beta/2, 1\}}$ .
Weak Convergence Rates:
- The paper establishes weak convergence rates for test functions $\phi$ with polynomial growth.
- In settings with independent noise components (Wiener or independent Lévy series), the weak rate is shown to be up to twice the strong rate ( $\kappa^{-2\beta}$ and $h^{\min\{\beta, 1\}}$ ).
- In the more general setting (Assumption 2.3, non-independent components), the weak rate depends on the integrability parameter $p$ of the noise, yielding a rate of roughly $h^{\min\{(\eta+1+\rho)/2, 1\}}$ where $\rho < \eta + 2/p$ . This highlights a trade-off: a broader class of test functions (larger $p$ ) leads to reduced convergence rates.
Regularity of Solutions:
- New regularity estimates are proven for the solution $X(t)$ in $L^p(\Omega, H^\rho(S^2))$ . These estimates are crucial for deriving the weak convergence rates and are shown to hold for both Gaussian and non-Gaussian Lévy drivers.
- The authors prove that for $p > 2$ , the independence of the noise components is a crucial assumption to obtain maximal spatial regularity rates independent of integrability conditions.

Numerical Verification
The theoretical convergence rates are validated through numerical simulations using both Wiener and Poisson noise (and their sums).

Spectral Method: Simulations confirm the theoretical strong convergence rates for various decay parameters $\alpha$ and verify the behavior of the mean and second moment errors.
Time Discretization: Both forward and backward Euler–Maruyama schemes are tested. The backward scheme demonstrates unconditional stability, while the forward scheme adheres to the predicted rates under the required stability constraints.
Weak Convergence: Monte Carlo simulations with specific test functions confirm that weak rates are indeed up to twice the strong rates in the independent noise setting.

Significance and Claims
The paper claims to extend the numerical analysis of SPDEs from the classical Gaussian setting on spheres to the more general and realistic setting of square-integrable Lévy processes.

It generalizes previous results on spectral approximations on spheres (previously limited to Wiener noise) to include jump processes.
It provides the first weak convergence analysis for SHE on the sphere driven by Lévy noise, showing that the "doubling of rates" phenomenon (weak rate $\approx 2 \times$ strong rate) holds under specific independence assumptions but is sensitive to the noise structure in more general cases.
The work establishes a rigorous framework for handling non-Gaussian noise on manifolds, offering a foundation for modeling complex systems with abrupt changes in spherical geometries.

The authors note that while the current work covers a broad class of Lévy processes (including compound Poisson, variance gamma, and jump-diffusion), processes beyond $L^2(\Omega)$ remain a subject for future work. Additionally, the numerical implementation currently relies on the independence of noise components to satisfy weak convergence assumptions, a limitation acknowledged for future research.