Pith. sign in

REVIEW 3 major objections 6 minor 15 references

Numerical study of the sharp stratification limit towards bilayer models

T0 review · 3 major / 6 minor · reviewed 2026-07-14 · grok-4.5

Pith's one-line read With shear, Kelvin-Helmholtz growth rates blow up as the pycnocline thins, so continuous stratified Euler cannot fully justify bilayer Euler or bilayer shallow-water models in Sobolev spaces.

desk verdict Solid V=0 justification with a proved rate; the shear-case numerics make a credible case that KH also blocks bilayer SW, but the δ^{-1} scaling remains unproven spectral fidelity. read the letter →

arxiv 2603.15287 v2 pith:6NH7Q6R5 submitted 2026-03-16 physics.flu-dyn math.AP

classification physics.flu-dynmath.AP MSC 35Q3576B7076M2286A05
keywords stratifiedEulerequationssharpstratificationlimitbilayermodelsKelvin-Helmholtzinstabilitydispersionrelationnormalmodespycnoclineshallow-water
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

Ocean models often replace a thin continuous density jump (a pycnocline of thickness δ) by two layers of constant density. Without background shear the paper proves that solutions of the linearized stratified Euler equations converge to the bilayer Euler equations as δ o0, with an explicit error of order δ^{1/2}|log δ|. With a sharp shear profile the same limit is obstructed by Kelvin-Helmholtz instabilities. A modal numerical scheme for the dispersion relation shows that the continuous system itself develops unstable modes whose growth rates scale like 1/δ (and whose wave numbers also scale like 1/δ). Those rates prevent uniform-in-δ energy estimates on any fixed time interval for generic Sobolev data. The same obstruction survives the shallow-water limit, so the bilayer shallow-water equations likewise cannot be fully justified from continuous stratification in Sobolev spaces, even though the bilayer model itself is well-posed.

What carries the argument

Normal-mode decomposition of the linearized stratified Euler equations (Sturm-Liouville eigenfunctions of the Taylor-Goldstein operator together with the resulting finite-dimensional matrices B_k) that converts the continuous dispersion relation into a computable eigenvalue problem whose imaginary parts track the Kelvin-Helmholtz growth rates.

What would settle it

A rigorous spectral analysis (or a high-resolution independent computation) of the Taylor-Goldstein operator for the family of sharp density-and-shear profiles that either confirms or refutes the conjectured scalings k_max,δ≈δ^{-1} and Im(ω*)≈δ^{-1}.

Watch

Extended reading notes

Core claim

In the presence of a shear profile that becomes discontinuous as δ o0, numerical computation of the dispersion relation of the linearized stratified Euler equations reveals Kelvin-Helmholtz modes whose maximal growth rates satisfy Im(ω*_δ)≈δ^{-1}. Consequently neither the bilayer Euler equations nor the bilayer shallow-water equations can be fully justified, in finite-regularity Sobolev spaces on a time interval independent of δ, as reduced models of the continuous system.

Load-bearing premise

The eigenvalues of the truncated modal matrices are assumed to approximate the true continuous spectrum closely enough that the observed 1/δ growth-rate scaling is not an artifact of truncation or continuous-spectrum numerical noise.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

3 major / 6 minor

Summary. The paper studies the sharp-stratification limit of the linearized stratified Euler equations in a strip as the pycnocline thickness δ o0. In the absence of shear, it proves quantitative convergence in Sobolev spaces of solutions toward the linear bilayer Euler equations (Proposition 6.1), using isopycnal coordinates and energy estimates on the difference. In the presence of a sharp shear profile V^δ, it develops a normal-mode discretization, computes the dispersion relation numerically, and reports Kelvin–Helmholtz-type unstable bands whose maximal growth rates appear to scale as Im(ω*_δ)≈δ^{-1} (and k*_δ≈δ^{-1}). From this numerical evidence the authors formulate Conjecture 6.5 and argue that such growth rates obstruct full justification, in finite-regularity Sobolev spaces on a δ-independent time interval, of both the bilayer Euler equations and the bilayer shallow-water equations as reduced models of the continuous system.

Significance. The stable-case result (Proposition 6.1 and Appendix B) is a genuine full justification with an explicit rate O(δ^{1/2}|log δ|(1+t)), carefully written energy estimates, and a clear use of isopycnal coordinates to compare continuous and bilayer unknowns. The modal numerical scheme (Section 5), its semi-discrete convergence proposition (Proposition 5.2), and the open scripts are strengths that make the unstable-case study reproducible. The claim that KH growth may also block justification of the well-posed bilayer shallow-water system is interesting and, if the scalings hold, would be a useful caution for geophysical modeling. The paper correctly separates proved statements from numerical conjectures.

major comments (3)
  1. The obstruction argument in §6.3 (and the bilayer-SW conclusion in particular) rests on Conjecture 6.5, whose δ^{-1} growth and k_max,δ≈δ^{-1} scalings are read from eigenvalues of the truncated matrices B_k (Definition 5.1, (5.8)–(5.9), (5.13)). Appendix C only checks residual consistency of reconstructed (c,w) pairs against the Taylor–Goldstein residual (C.2)–(C.4). As the authors note, the spatial operator T is neither self-adjoint nor skew-adjoint, so a small residual does not imply spectral proximity. Continuous-spectrum pollution for Re(c) in the range of V^δ is visible, and k_max,δ is extracted with an ad-hoc Im(c)>3·10^{-3} threshold (Fig. 17). Without a mode-refinement study that tracks Im(ω*_δ) and k*_δ under simultaneous increase of vertical modes ℓ and Fourier cutoff, or an independent diagnostic (e.g. direct discretization of the TG equation), the claim that growing modes wi
  2. Proposition 5.2 proves convergence of the modal truncation for the evolution problem under V∈W^{2,∞}, but the dispersion computation in §6.2 uses the same truncation for a family of profiles with ||V^δ'||_∞∼1/δ, which is not uniform in δ. The energy estimate of Proposition 3.7 likewise requires V∈W^{1,∞} with constants that blow up as δ o0. The paper should clarify whether the observed unstable band and its δ-scalings remain stable under this non-uniformity, or whether the semi-discrete spectrum could be polluted by the increasingly steep shear layer for the values of δ and ℓ used in Figs. 18–22.
  3. In the stable case, the measured numerical rate for Err versus δ (Fig. 11, §6.1) is roughly δ^{0.56} and the points are not well aligned; the authors note that this does not conclusively confirm the theoretical δ^{1/2}|log δ| rate of Proposition 6.1. Given that the theoretical rate is the main proved contribution, a short discussion of why the numerical rate is inconclusive (resolution of the pycnocline, number of modes, choice of initial data (6.4), or the L^∞-in-time error definition (6.5)) would strengthen confidence that the numerics and analysis are consistent.
minor comments (6)
  1. Several figures (e.g. Figs. 14–16, 23) are hard to read in grayscale; consider distinct markers or line styles in addition to color.
  2. Notation for the number of vertical modes switches between ℓ, N, and ℓ in captions and text (e.g. §6.2.2); unify.
  3. The date line reads “March 17, 2026”; confirm this is intentional.
  4. In (3.14) and Lemma 3.1 the asymptotic n c_n o c is stated without an explicit reference for the constant; a pointer to [AM87] is given later but could be placed at first use.
  5. Typographical: “Saint-Andrew cross” / “Saint Andrew’s cross” appear in both forms; standardize. Occasional missing spaces before parentheses and “Grönwall’s” spelling vary.
  6. Remark 6.6 on analytic spaces is useful; a short pointer to the vortex-sheet literature already cited ([SSBF81], [CO86]) in the introduction would help non-specialist readers.

Circularity Check

1 steps flagged · score 1.0 of 10

No load-bearing circularity: V=0 convergence is proved from energy identities; KH scalings are read off eigenvalues of matrices built from the PDE, not fitted or defined into existence.

  1. self citation load bearing [§1 Introduction; Prop. B.1 / [Fra24]]
    "the well-posedness of (3.1) together with the boundary conditions (3.2) and suitable initial conditions in Sobolev spaces is studied in [DLS20] and [Fra24]. ... Recall that this is in contrast with well-posedness results on the stratified Euler equations, and in particular with [Fra24], which includes a non-zero shear flow. However the latter result is not uniform in δ."

    The author’s own prior well-posedness paper is cited for the continuous stratified system. This is ordinary background and is not used to force the sharp-limit convergence rate or the numerical KH scalings; those rest on independent energy estimates (App. B) and eigenvalue computations of B_k. Flagged only as minor non-load-bearing self-citation.

full rationale

The paper’s two main strands are independent of circular constructions. For V=0, Proposition 6.1 and Appendix B derive uniform energy estimates and an O(δ^{1/2}|log δ|) difference bound between the continuous isopycnal system (3.6) and the bilayer system (4.10) by direct comparison of the PDEs; nothing is fitted and the result does not rest on a self-citation uniqueness theorem. For V=V^δ the dispersion relation is obtained by truncating the modal system (3.34) to finite vertical modes, forming the matrices B_k from the explicit integral coefficients M,A_j of ρ_δ and V_δ, and reading phase velocities from eigenvalues (Def. 5.1, (5.8)–(5.13)). Those eigenvalues are numerical outputs of the discretized operator, not parameters adjusted to data, and the claimed scalings Im(ω*_δ)≈δ^{-1}, k*_δ≈δ^{-1} are observations (Conjecture 6.5), not forced by normalization. Self-citations (Fra24 well-posedness, Duc22 bilayer formulae) supply background that is independent of the sharp-limit statements; they are not used to forbid alternatives or to smuggle an ansatz that already encodes the conclusion. The skeptic’s concern about fidelity of truncated eigenvalues to the continuous spectrum is a correctness/approximation issue (Appendix C only checks residual consistency of a non-self-adjoint operator), not circularity. Score 1 only for ordinary non-load-bearing self-citation of the author’s prior well-posedness result.

Assumptions & free parameters 3 free parameters · 6 assumptions · 0 invented entities

The paper works entirely within classical inviscid stratified Euler theory. Load-bearing inputs are standard domain assumptions (stable stratification, rigid lid, linearization, 2D strip) plus specific profile families and numerical cutoffs used only for the unstable-case evidence. No new physical entities are postulated.

free parameters (3)
  • Im-part threshold for k_max,δ detection = 3e-3
    Set by hand to 3·10^{-3} to separate continuous-spectrum numerical noise from KH instabilities when measuring k_max,δ vs δ (§6.2.2).
  • Vertical mode count � and Fourier count � = varies (e.g. 40–90 modes, up to 1200 Fourier)
    Truncation parameters of the semi-discrete scheme; results for k_max and growth rates are reported for several choices (e.g. �=80/90, �=1200) but remain discretization-dependent.
  • Shear amplitude v and density jump (ρ−,ρ+) = v=0.25, ρ−=1.5, ρ+=0.75
    Fixed numerical values (v=0.25, ρ−=1.5, ρ+=0.75) used for all unstable dispersion plots; not fitted to data but chosen by hand and affect the bilayer Im(c_max,BL).
assumptions (6)
  • domain assumption Background density is smooth and stably stratified: −ρ′ ≥ c* > 0 (and ρ bounded above and below by positive constants).
    Stated in (3.3) and used throughout for Sturm–Liouville bases, energy estimates, and well-posedness; standard in GFD to avoid Rayleigh–Taylor.
  • domain assumption Rigid-lid, flat-bottom strip domain T_L × [−H,0]; no Coriolis, wind stress, or topography.
    §1 and Fig. 1; simplifies the linear operator and the modal analysis.
  • domain assumption Linearization about a pure shear equilibrium (V(z),0) with continuous or sharp profiles ρ_δ, V_δ.
    Systems (3.6) and (3.32); the whole comparison is linear.
  • standard math Sturm–Liouville eigenfunctions (f_n),(g_n) form orthonormal bases of the weighted L2 spaces and c_n ∼ c/n.
    Lemma 3.1, classical SL theory (AM87, DLS20).
  • domain assumption Bilayer Euler with nonzero shear is ill-posed in Sobolev spaces due to KH (cited).
    Used as background for the unstable-case discussion (ITT97, KL05, Duc22).
  • ad hoc to paper Profile family ρ_δ (and V_δ) of arctan type, or more generally profiles satisfying (B.14)–(B.16).
    Chosen to realize a sharp pycnocline of width δ; Prop. 6.1 extends to the class in Remark B.3, but all numerics use the arctan family.

how reviews work

0 comments
Cite this review

Pith. "Pith review of Numerical study of the sharp stratification limit towards bilayer models." pith.science (2026). https://pith.science/paper/6NH7Q6R5

@misc{pith2026260315287,
  author       = {Pith},
  title        = {Pith review of: Numerical study of the sharp stratification limit towards bilayer models},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/6NH7Q6R5}},
  note         = {Machine review of arXiv:2603.15287}
}
read the original abstract

In the study of oceanic flows at the geophysical scale, the phenomenon of density stratification plays a central role in the dynamics of the system. Two categories of mathematical models are commonly used to describe the role played by the density stratification: on the one hand, continuously stratified models - such as the stratified Euler equations in a strip, considered in the present article - offer an accurate description of vertical effects, but come with a high level of complexity, both at the theoretical and numerical levels. On the other hand, bilayer models approximate the stratification by a piecewise constant profile. In the latter case, the main point is to study the evolution of the free interface between both layers, which leads to a substantially simplified model. In the present article, we compare both approaches in the framework of the linearized stratified Euler equations around density profiles that are close to piecewise constant profiles, and prove the convergence towards the bilayer Euler equations. However, in the presence of a shear flow, bilayer models have a range of validity limited by the presence of Kelvin-Helmholtz instabilities. In this case, we use a suitable normal modes decomposition to compute numerically the dispersion relation of this linearized model, and provide numerical evidence that the Kelvin-Helmholtz instabilities limit the applicability of two widely used bilayer models, namely the bilayer Euler equations and the bilayer shallow-water equations.

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

15 extracted references · 2 linked inside Pith

  1. [1]

    [ABD24]Mahieddine Adim, Roberta Bianchini, and Vincent Duchêne,Relaxing the sharp density stratification and columnar motion assumptions in layered shallow water systems, C. R. Math. Acad. Sci. Paris362(2024), 1597–1626. MR 4834569 [AFM24]Helmut Abels, Julian Fischer, and Maximilian Moser,Approximation of classical two-phase flows of viscous incompressibl...

  2. [2]

    MR 4793824 [AGP24]Helmut Abels, Harald Garcke, and Andrea Poiatti,Mathematical analysis of a diffuse interface model for multi-phase flows of incompressible viscous fluids with different densities, J. Math. Fluid Mech.26(2024), no. 2, Paper No. 29,

  3. [3]

    MR 4726194 [AM87]F . V . Atkinson and A. B. Mingarelli,Asymptotics of the number of zeros and of the eigenvalues of general weighted Sturm-Liouville problems, J. Reine Angew. Math.375�376(1987), 380–393. MR 882305 [Ban05]Aravind Banerjee,Taylor-goldstein equation and stability, arXiv preprint physics/0510114 (2005). [BD24]Roberta Bianchini and Vincent Duc...

  4. [4]

    [BGC18]Sylvie Benzoni-Gavage and David Chiron,Long wave asymptotics for the Euler-Korteweg system, Rev. Mat. Iberoam.34(2018), no. 1, 245–304. MR 3763346 [BGDDJ05]Sylvie Benzoni-Gavage, Raphaël Danchin, Stéphane Descombes, and Didier Jamet,Structure of korteweg models and stability of diffuse interfaces, Interfaces and free boundaries7(2005), no. 4, 371–4...

  5. [5]

    MR 2676605 41 [BV20]Ricardo Barros and José Felipe Voloch,Effect of variation in density on the stability of bilinear shear currents with a free surface, Physics of Fluids32(2020), no

  6. [6]

    6, 807–838

    [CO86]Russel E Caflisch and Oscar F Orellana,Long time existence for a slightly perturbed vortex sheet, Communications on Pure and Applied Mathematics39(1986), no. 6, 807–838. [CZN25]Michele Coti Zelati and Marc Nualart,Limiting absorption principles and linear inviscid damping in the Euler-Boussinesq system in the periodic channel, Comm. Math. Phys.406(2...

  7. [7]

    1, 105–143

    [DGKR12]Wolfgang Dreyer, Jan Giesselmann, Christiane Kraus, and Christian Rohde,Asymptotic analysis for Korteweg models, Interfaces Free Bound.14(2012), no. 1, 105–143. MR 2929127 [DLS20]Benoît Desjardins, David Lannes, and Jean-Claude Saut,Normal mode decomposition and dispersive and nonlinear mixing in stratified fluids, Water Waves3(2020), no. 1, 153–1...

  8. [8]

    Zimmerman,An introduction to internal waves, Lecture Notes, Royal NIOZ, Texel207(2008),

    [GZ08]Theo Gerkema and Joseph T .F . Zimmerman,An introduction to internal waves, Lecture Notes, Royal NIOZ, Texel207(2008),

Show all 15 references
  1. [9]

    17, Springer-Verlag, New York-Berlin, 1982, Graduate Texts in Mathematics,

    [Hal82]Paul Richard Halmos,A Hilbert space problem book, second ed., Encyclopedia of Mathematics and its Applications, vol. 17, Springer-Verlag, New York-Berlin, 1982, Graduate Texts in Mathematics,

  2. [10]

    2, 199–223

    MR 675952 [ITT97]Tatsuo Iguchi, Naoto Tanaka, and Atusi Tani,On the two-phase free boundary problem for two-dimensional water waves, Math- ematische Annalen309(1997), no. 2, 199–223. [Jam01]Guillaume James,Internal travelling waves in the limit of a discontinuously stratified ...

  3. [11]

    MR 4901758 [Lan13a]David Lannes,A stability criterion for two-fluid interfaces and applications, Arch. Ration. Mech. Anal.208(2013), no. 2, 481–567. [Lan13b] ,The water waves problem, American Mathematical Society,

  4. [12]

    LaCasce and Sjoerd Groeskamp,Baroclinic modes over rough bathymetry and the surface deformation radius, Journal of Physical Oceanography50(2020), no

    [LG20]Joseph H. LaCasce and Sjoerd Groeskamp,Baroclinic modes over rough bathymetry and the surface deformation radius, Journal of Physical Oceanography50(2020), no. 10, 2835 –

  5. [13]

    shallow water

    [Lin03]Zhiwu Lin,Instability of some ideal plane flows, SIAM J. Math. Anal.35(2003), no. 2, 318–356. MR 2001104 [Ovs79]Lev Vasilevich Ovsjannikov,Models of two-layered “shallow water”, Zh. Prikl. Mekh. i Tekhn. Fiz.,(2) (1979), 3–14. [Sch14]Stefan Schaubeck,Sharp interface lim...

  6. [14]

    Sutherland and Paul F

    [SL99]Bruce R. Sutherland and Paul F . Linden,An experimental�numerical study of internal wave transmission across an evanescent level, Inst. Math. Appl. Conf. Ser. New Ser., vol. 68, Oxford University Press, 1999, pp. 251–262. [SNK13]David N. Sibley, Andreas Nold, and Serafim...

  7. [15]

    [VC25]Philipp P Vieweg and Colm-cille P Caulfield,Anisotropy of emergent large-scale dynamics in forced stratified shear flows, arXiv preprint arXiv:2507.14991 (2025). UNIV. BORDEAUX, CNRS, BORDEAUXINP , IMB, UMR 5251, F-33400 TALENCE, FRANCE UNIVRENNES, CNRS, IRMAR - UMR 6625...

Pith tools

Reviewed July 14, 2026 · model on record in the stance chip above.