Source-linked AI summary
A PDE-theory-guided neural operator and the long-time relaxation of the two-dimensional anisotropic Boussinesq system
Akram Moustafa, Jiahong Wu, Xiaoming Zheng
TL;DR
The paper asks how small perturbations of the two-dimensional anisotropic Boussinesq system behave at long times and what determines their limit. It uses rigorous long-time analysis to design a neural operator that predicts late-time vorticity, then evaluates its accuracy and uses it to search for larger values of the theorem-bounded quantity. The operator’s predictions reproduce that quantity on test runs, and solver-checked search candidates raise its observed lower bound.
Problem
The long-time limit of a small perturbation and the part determined by its initial data are not settled by stability results alone.
Source: §1 · p. 3Method
The analysis guides a neural operator to predict late-time vorticity as the initial horizontal average plus separately scaled changes in the average and oscillation.
Source: §1 · p. 3, §5 · p. 18Results
A single network reproduces the theorem-bounded change in horizontally averaged vorticity with 2.2 percent median error over 2 × 10^3 test runs, and solver-checked search finds a value 4.4 times the test-run maximum.
Source: §5 · pp. 18–19Takeaways & Limitations
The trained operator can search initial data for larger values of the theorem-bounded quantity, with every candidate checked by the solver.
Source: §5 · p. 19The late-time temperature map is severely ill-conditioned, so temperature is not used as a prediction target.
Source: §4.3 · p. 10
Abstract
from arXiv · showhide
We design a neural operator for the two-dimensional Boussinesq system with anisotropic dissipation on the periodic box, guided by a rigorous analysis of its long-time behavior. The operator maps a small perturbation of hydrostatic equilibrium, of amplitude $\varepsilon$, to the vorticity at a late time. We prove that the flow relaxes to a stratified shear flow, with the horizontally averaged vorticity staying within $O(\varepsilon^2)$ of its initial value $\bar\omega_0$ and the oscillation about it of order $\varepsilon$. Accordingly, the network outputs $\bar\omega_0$ plus these two parts, each scaled to its size, which grows as $\varepsilon^2$ and $\varepsilon$, respectively. Without this guidance the training either does not converge or misses the $\varepsilon^2$ growth. With it, one network reproduces, sample by sample on $2\times10^3$ test runs over a factor of $23$ in amplitude, the quantity bounded by the theorem to a median error of $2.2$ percent. Being fast and differentiable with respect to the initial data, the network is also used to search the space of initial data, and with every candidate checked by the solver the search finds initial data for which this quantity is $4.4$ times the largest value among the test runs, raising the lower bound for the constant of the theorem by that factor.
1. Introduction.
The paper asks what small perturbations of the anisotropic Boussinesq system approach at long times, beyond the stability question of whether they remain small. It uses rigorous long-time analysis to shape a neural operator that predicts the vorticity’s stratified limit, separating the initial mean, its small change, and the oscillation.
1. Introduction.
- The paper’s approach uses rigorous solution properties to choose the network’s predictions, output split, and amplitude scaling, then uses the trained operator to search for a lower bound confirmed by the solver.This differs from building network structure from equation-level or formal-limit information alone.Source: §1 · p. 2
- The long-time question is what a globally small perturbation converges to and how much of its limit is determined by its initial data.Earlier stability results do not by themselves settle this question.Source: §1 · p. 3
- With horizontal viscosity and vertical thermal diffusion, the temperature decays to zero while vorticity tends to a stratified limit depending only on x2.The anisotropic dissipation leaves modes undamped, so buoyancy coupling must supply stabilization.Source: §1 · p. 3
- The horizontally averaged vorticity changes by O(ε^2), while the vorticity oscillation is of order ε; the late-time flow is a stratified shear.The mean change arises through a nonlinear Reynolds-stress contribution to its evolution.Source: §1 · p. 3
- The operator predicts the known initial horizontal average plus separately scaled approximations of its change and the oscillation, omitting temperature because it decays to zero.The two learned parts scale with their respective sizes, ε^2 and ε.Source: §1 · p. 3
- A standard whole-field target obscures the small mean change, while different amplitude scaling and ill-conditioned late-time temperature dependence further complicate learning.The paper trains on solver-generated initial-to-final pairs and tests the analytical bound on solver runs and surrogate predictions.Source: §1 · p. 4
2. Governing equations and the mean/oscillation decomposition.
The paper formulates the perturbation equations on a periodic square and decomposes fields into horizontal means and zero-mean oscillations. These components separate the directions reached by the two anisotropic diffusion terms.
2.1. The perturbation system.
- The perturbation system is posed on the periodic box [0, 2π]^2, with horizontal viscosity ν and thermal diffusivity η.The velocity and temperature are perturbations from hydrostatic equilibrium.Source: §2.1 · p. 4
- The vorticity formulation uses ω = ∂1u2 − ∂2u1, alongside the velocity and stream function variables.The perturbation velocity is denoted u, with initial value u0.Source: §2.1 · p. 5, §2.2 · p. 5
2.2. Mean/oscillation decomposition.
- Each periodic field is split into its x1-average and a zero-x1-mean oscillation; a field independent of x1 is called stratified.The mean and oscillation are orthogonal in L2, and their products decompose into mean and oscillatory contributions.Source: §2.1 · p. 5, §2.2 · pp. 5–6
- Horizontal viscosity damps the oscillation but leaves the horizontal mean untouched, whereas vertical diffusion annihilates the x1-average and damps its complement.The horizontal and vertical Poincaré inequalities bound the oscillatory components using derivatives in their respective directions.Source: §2.2 · p. 5
3. Main result: stratification and the limit profile.
For sufficiently small initial data, the system has a global solution that stratifies: the vertical velocity and temperature vanish, while vorticity converges to a horizontal mean profile. The paper also identifies how the mean vorticity evolves and shows that linearized temperature decay has no uniform rate.
3.1. Standing hypotheses on the initial data.
- Under the standing smallness hypotheses, the solution is unique and global, and the initial zero-mean conditions remain preserved.The initial data are assumed to lie in H2 and satisfy the stated divergence and mean conditions.Source: §3.2 · p. 6, §3 · p. 6, §3.1 · p. 6
3.2. Statement of the main result.
- The horizontal mean vorticity evolves through the vertical second derivative of the product of the velocity oscillations, ∂tω̄ = ∂22(eu1eu2).This identity follows from the averaged horizontal momentum equation and the zero-mean conditions.Source: §3.2 · p. 6
- For sufficiently small data, u2 and θ vanish while ω converges to an x2-dependent limiting profile ω̄∞.For fixed u0, admissible initial temperatures yield limits differing by at most 2Cε² in H−2; C depends only on ν.Source: §3.2 · pp. 6–7
3.3. The linear decay rate.
- The linearized temperature decay has no uniform rate because λ+ approaches zero as |k1| grows; the nonlinear theorem likewise provides no decay rate.Source: §3.3 · p. 7
4. A neural-operator surrogate, its training data, and a test of Theorem 3.1(b).
The paper uses a neural operator to map small-amplitude initial vorticity and temperature fields to later vorticity. It specifies an FNO and generates simulation data from independently sampled initial fields across amplitudes.
4. A neural-operator surrogate, its training data, and a test of Theorem 3.1(b).
- The neural operator maps the initial fields (ω0, θ0), across varying amplitudes, to the later vorticity ω(·, T).Source: §4 · p. 8
4.1. Neural operators and the Fourier neural operator.
- The FNO parameterizes translation-invariant kernels in Fourier space, retaining low modes and truncating higher modes for spectral smoothing and reduced parameter count.On an N × N grid, its FFT-based cost is O(N^2 log N).Source: §4.1 · p. 8
4.2. Training data generation.
- Initial vorticity and temperature are independently sampled as zero-mean Gaussian fields, then rescaled to satisfy the theorem’s smallness condition.The vorticity and temperature use distinct covariance parameters; each sample stores initial and final fields on a 128 × 128 grid.Source: §4.2 · pp. 9–10
- The training set contains 10^4 samples across ten amplitudes, with separate training-amplitude and interpolation test sets excluded from the loss.Simulations use Fourier pseudo-spectral spatial discretization and a second-order semi-implicit time integrator.Source: §4.2 · pp. 8–10
4.3. The late-time temperature and why it is not a prediction target.
At late time, the temperature is dominated by horizontally varying modes, while its initial-to-final map is severely ill-conditioned. The network instead predicts structured vorticity corrections, scaled by their distinct amplitude growth.
4.3. The late-time temperature and why it is not a prediction target.
- At T = 100, k2 = 0 modes carry a median 99.9999% of temperature energy, and the predicted spectral average differs by at most 0.6% for 7 ≤ |k1| ≤ 20.Across runs, the temperature norm falls by about 50-fold; the peak wavenumber has median 10, near k*1 = 9.4.Source: §4.3 · p. 10, Figure 1 · p. 11
- The late-time temperature map is severely ill-conditioned: its measured modal multiplier ranges from 4.8 × 10^-11 to 1.00 across the reported modes.A preliminary network trained to predict both final fields did not converge, so the network predicts vorticity alone.Source: §4.3 · pp. 10–11
4.4. Network structure and training.
- The network uses four Fourier layers with 64 hidden channels and is trained for 600 epochs using a relative-error loss that weights the mean correction ten times more heavily.Training takes 10 h 59 min on one NVIDIA Tesla V100S GPU.Source: §4.4 · pp. 11–12
- The output adds learned corrections to the exact input mean vorticity, separating the change in horizontal average from the remaining oscillation.The two output channels provide the mean correction and oscillation components.Source: §4.4 · pp. 11–12
- The output components are scaled by their amplitude-dependent sizes: the mean correction grows as ε^2 and the oscillation as ε.The scales are interpolated between training amplitudes; normalized outputs remain order one.Source: §4.4 · p. 12
4.5. Surrogate accuracy and cost.
The surrogate accurately predicts late-time vorticity and the theorem-bounded mean-vorticity quantity across test amplitudes, while evaluation is far cheaper than the spectral solver. Its component errors differ substantially, and fine-scale errors in the mean shift are less prominent in the H−2 norm.
4.5. Surrogate accuracy and cost.
- The surrogate reproduces final vorticity with median relative L2 errors of 3.8 × 10−5 to 2.4 × 10−4 at interpolation amplitudes and 4.0 × 10−5 to 3.3 × 10−4 at training amplitudes.The interpolation amplitudes were not used in training.Source: §4.5 · p. 12, Table 1 · p. 15
- The median relative L2 error is 4.0 to 4.7 percent for the oscillation and 23 to 27 percent for the mean shift, while both learned parts are small relative to the final field.For the mean shift, 95 to 99 percent of squared L2 error lies at |k2| > 16, compared with 22 to 26 percent for the mean shift itself; its median relative H−2 error is 6 to 10 percent.Source: §4.5 · pp. 12–13
- At batch size 16, surrogate inference takes 1.33 ms per field, about 8 × 10^3 times faster per field than the solver’s 10.2 s simulation.On one NVIDIA Tesla V100S, a batch of 64 takes 73.6 ms, or 1.15 ms per field.Source: §4.5 · p. 13
4.6. The quantity of Theorem 3.1(b): solver and surrogate.
- The mean of the normalized mean-vorticity change stays between 1.33 × 10−6 and 1.55 × 10−6 across amplitudes 0.021 to 0.476, with fitted growth exponent p = 1.99.The 95 percent bootstrap interval for the exponent is [1.95, 2.03], and the maximum stays between 2.45 × 10−6 and 4.24 × 10−6.Source: §4.6 · p. 14
- 2.2 percent is the median relative error for the theorem-bounded mean-vorticity quantity across 2 × 10^3 test runs.At every amplitude, the surrogate’s mean, maximum, and standard deviation agree with the solver’s within 7, 5, and 8 percent, respectively.Source: §4.6 · pp. 15–16, Figure 4 · p. 15
- The normalized full-field quantity R is nearly constant at small amplitudes, then rises to a sample mean of 7.4 × 10−7 at ε = 0.476.The first-order oscillation dominates at small ε; the second-order mean shift becomes comparable at larger ε, and the surrogate’s R matches the solver’s at a median of 0.9 percent.Source: §4.6 · p. 16
4.7. Comparison with Fourier neural operators trained without the analysis.
Fourier neural operators trained without the analysis can fit training fields yet substantially miss the small mean-vorticity change and its second-order growth. The surrogate also guides a search whose solver-checked candidates yield a value 4.4 times the test-set maximum, although later ascent steps overestimate the target.
4.7. Comparison with Fourier neural operators trained without the analysis.
- The plain and scaled networks predict mean-shift errors 30 to 37 times the size of the mean shift, with normalized values 6 to 180 times too large on average.Their predicted growth exponents are ε^0.95 to ε^0.98 for the scaled network and ε^1.17 to ε^1.18 for the plain network, rather than the solver-confirmed exponent 2.00.Source: §4.7 · pp. 16–17
- Both comparison networks fit training fields below 10−3 relative L2 error and reproduce test fields to about 3 percent, but neither adds the initial average exactly or separates and scales the two output parts.The plain network predicts the whole final vorticity; the scaled network predicts the whole field divided by ε.Source: §4.7 · p. 16
4.8. Surrogate investigation of a lower bound for the constant of The-
- The search keeps each amplitude fixed and constrains initial-field norms to ε/2, while modifying phases and the distribution of Fourier-coefficient moduli within annuli.The ascent starts from 40 samples and maximizes log R_NO by Adam for 2000 steps; the constraints keep fields in the class on which the surrogate was trained.Source: §4.8 · pp. 17–18
- Screening 5 × 10^6 new samples took 83 minutes of GPU time, versus about 1.6 years for the solver; solver checks of the top 80 yielded a largest value of 8.32 × 10−6.The surrogate’s largest value among the screened samples was 7.34 × 10−6.Source: §4.8 · p. 17
- The solver finds a mean-vorticity quantity of 1.87 × 10−5 along the ascent, 4.4 times the test-set maximum, before later surrogate overestimation causes the solver value to decrease.The peak is at step 500; recomputation on a 256 × 256 grid or with half the time step changes it by less than 10−5 in relative terms.Source: §4.8 · p. 18
- The reported lower bound C ≥ 1.87 × 10−5 for ν = η = 1 applies only if the tested amplitudes lie below the theorem’s threshold ε*.The candidates’ reported R values are recomputed with the spectral solver.Source: §4.8 · pp. 17–18
5. Discussion and conclusions.
The paper combines long-time analysis with a neural operator that predicts the late-time vorticity by separately scaling its mean change and oscillation. The guided network accurately predicts the theorem-bounded quantity and supports a solver-checked search for larger values.
5. Discussion and conclusions.
- 2.2% median error and at most 12% over 2 × 10^3 test runs for the H−2 norm of the change in horizontally averaged vorticity.One network trained on 10^4 solver runs across ten amplitudes spanning a factor of 23 also captures the quantity’s ε² growth and the vorticity oscillation’s order ε.Source: §5 · p. 18
- The network predicts initial mean vorticity plus separate approximations of its change and the remaining oscillation, scaling the outputs by their respective ε² and ε orders.The loss measures each part’s relative error separately and weights the mean change ten times more than the oscillation; only final-time vorticity at T = 100 is predicted.Source: §5 · p. 18
- About 30× median error on the theorem-bounded quantity and approximately linear growth in ε for the mean change resulted from training with the standard whole-vorticity loss.Including final temperature as a target in a preliminary experiment did not converge.Source: §5 · p. 19
- 4.4 times the largest test-set value: solver-checked gradient ascent found a change in horizontal-average vorticity of 1.87 × 10^-5.The search used the surrogate’s predictions and gradients; Fig. 5 compares solver and surrogate values along ascent paths and marks the test-set maximum.Source: §5 · p. 19, Figure 5 · p. 19
- The approach is proposed for other equations where analysis identifies a nearly conserved quantity and the amplitude order of each solution part.The trained operator can also estimate quantities such as constants in estimates that the analysis does not provide.Source: §5 · p. 19
Appendix A. The uniform stability estimate.
Appendix A records a uniform stability estimate for small initial data on T2 and an auxiliary bound used in the proof of Theorem 3.1. Its proof bounds a remaining term by decomposing a derivative and applying vertical Poincaré.
Appendix A. The uniform stability estimate.
- For ε ≤ ε*, divergence-free H2 initial data satisfying (3.1) yield a unique global solution, with ε* depending on ν and η.Theorem A.1 is attributed to.Source: Appendix A · p. 20
- The thermal diffusivity η enters the stability result only through the smallness threshold ε*, while Cν depends on ν.The supplied estimate gives Cν = 64(2 + ν)/min{1, ν²}.Source: Appendix A · p. 20
- Theorem A.1 and Remark A.2 provide the estimates used in Theorem 3.1’s proof; the part (b) constant depends only on ν.The proof bounds the last term by decomposing ∂1θ into ∂1Q0θ and ∂1θ⊥, applying vertical Poincaré to the zero-mean component, and combining the resulting bound with (A.4).Source: Appendix A · p. 20
Appendix B. Proof of Theorem 3.1.
Appendix B uses Barbalat’s lemma to turn finite dissipation integrals into convergence to zero. The stated criterion is finite integral together with uniform continuity, for example Lipschitz continuity.
Appendix B. Proof of Theorem 3.1.
- Finite dissipation integrals are converted into convergence to zero using Barbalat’s lemma [3, 11].The appendix records the lemma before applying it in the proof.Source: Appendix B · p. 20
- In the proof, finite dissipation integrals yield decay only when the relevant functions also meet the lemma’s uniform-continuity condition.The appendix identifies Lipschitz continuity as an example of that condition.Source: Appendix B · p. 20
- A uniformly continuous function g with finite integral from 0 to infinity satisfies g(t) → 0 as t → ∞.Lipschitz continuity is given as one sufficient regularity condition.Source: Appendix B · p. 20
B.1. Proof of part (a).
The proof of Theorem 3.1(a) uses uniform H2 bounds and Barbalat’s lemma to establish decay of oscillatory velocity and temperature components. It justifies the needed Lipschitz bounds for squared derivative norms through an H−1–H1 pairing.
B.1. Proof of part (a).
- The proof uses a global H2 bound of Cε, with constants depending on ν and η in part (a), while part (b) tracks constant dependence explicitly.The global solution and bound follow from Theorem A.1 and Remark A.2.Source: §B.1 · p. 21
- Squared L2 norms of first derivatives are Lipschitz because their time derivatives can be bounded using ∂tv in L2 and second spatial derivatives in L2.The proof uses an H−1–H1 duality pairing because ∂tDv is only known in H−1, so the ordinary L2 chain-rule pairing is not available.Source: §B.1 · p. 21
- The averaged vertical velocity vanishes, so u2 is entirely oscillatory.The proof uses the zero-mean condition and divergence-free property.Source: §B.1 · p. 21
- Finite dissipation and Barbalat’s lemma give ∥∂1u(t)∥L2 → 0; horizontal Poincaré and interpolation then yield decay of the oscillatory velocity in H1.The relevant derivative norm is shown to be Lipschitz before applying the lemma.Source: §B.1 · p. 21
- The temperature derivatives decay in L2, and the zero mean of θ plus Poincaré give temperature decay; interpolation with uniform H2 bounds yields u2 and θ decay in Hs for s < 2.The derivative limits follow by applying Barbalat’s lemma to the two temperature derivative norms.Source: §B.1 · pp. 21–22
B.2. Proof of part (b).
The vorticity becomes increasingly stratified, and the full field converges to a stratified limit. The horizontally averaged vorticity remains within O(ε^2) of its initial value, with convergence also holding in H^s for every s < 1.
B.2. Proof of part (b).
- The oscillatory vorticity tends to zero in L2, so the vorticity becomes stratified.The proof bounds its L2 norm by the H1 norm of the oscillatory velocity, which tends to zero.Source: §B.2 · p. 22
- The horizontal average is Cauchy in H^-2 and stays within Cε^2 of its initial value for all t ≥ 0.The bound also applies to the difference between the limiting average and the initial average; the constant depends only on ν.Source: §B.2 · p. 22
- The full vorticity converges in H^-2 to the limiting average, and its distance from the initial average tends to a value bounded by Cε^2.This follows by combining decay of the oscillatory part with convergence of the horizontal average.Source: §B.2 · p. 22
- The limiting average belongs to H1, and the vorticity converges to it in H^s for every s < 1.The proof obtains the limit's H1 bound by Fatou's lemma and interpolates convergence from H^-2 and a uniform H1 bound.Source: §B.2 · p. 22
Appendix C. Proof of Proposition 3.2.
The proof rewrites the Fourier system using the Biot–Savart law and treats different Fourier-mode cases separately. It uses energy estimates or the mode equations to establish decay bounds.
Appendix C. Proof of Proposition 3.2.
- The Biot–Savart law converts the Fourier transform of the first system equation into the first equation of (3.6), while the second follows directly.The proof uses ˆu2 = −ik1ˆω/|k|^2 for nonzero k.Source: Appendix C · p. 22
- For modes with k1 = 0 and k2 ≠ 0, ˆu2 vanishes and θ satisfies dˆθ/dt = −ηk2^2ˆθ.This mode equation gives decay of θ for η > 0.Source: Appendix C · p. 23
- For modes with k2 = 0 and k1 ≠ 0, the proof uses the matrix roots of (3.6) to represent ˆθ when the roots are real, distinct, and negative.For modes with both components nonzero, the proof also invokes Grönwall's inequality to obtain a bound.Source: Appendix C · p. 23