{"id":"d0dcca24-7496-4271-9263-79c158a9f1a5","arxiv_id":"2505.22480","paper_version":2,"verdict":"CONDITIONAL","confidence":"HIGH","novelty_score":4.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":3,"one_line_summary":"An ESDIRK-HHO scheme for incompressible Navier-Stokes is validated on manufactured solutions and applied to the Taylor-Green vortex at Re=1600, with the turbulence run made using an older code variant lacking the proposed L2 projections.","lead":"This paper tests a hybrid high-order numerical method for incompressible fluid flow that combines a hybrid velocity-pressure discretization with high-order implicit Runge-Kutta time stepping. The authors validate the scheme on benchmark problems and apply it to the Taylor-Green vortex at Reynolds 1600, claiming it can reproduce key turbulence statistics on coarse meshes.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The TGV turbulence evidence was produced by an older code lacking the L2-projected convective term, so the paper's central DNS-like-precision claim is not yet demonstrated for the proposed formulation.","rationale":"The reader's weakest assumption identifies exactly the load-bearing concern: the TGV demonstration was generated with an older code that omits a component of the proposed formulation, and the paper asserts only that the impact is marginal. My read agrees with this assessment. The rest of the paper—convergence tests, pressure-robustness, inviscid double shear layer—is internally consistent and supports the mathematical properties of the ESDIRK-HHO discretization, so the issue is not with the method per se but with the central turbulence claim being attributed to the wrong code. Because the reader's verdict is already CONDITIONAL and this concrete test is a natural condition for acceptance, I recommend no change to the verdict. If the rerun shows negligible differences, the concern is resolved and the paper's headline claim would be supported; if it shows material differences, the TGV section would need to be revised or labeled as a preliminary study with the older variant.","tokens_in":32686,"tokens_out":4228,"duration_ms":53600,"concrete_test":"Rerun the TGV case on the coarse 24×4^3 grid at k=7 and k=9 (and, if feasible, on the fine grid at k=8) with the updated code that includes the L2 projections in the convective term, keeping all other settings (ESDIRK5, tolerance tola, mesh) identical. Compare relative enstrophy E/E0, kinetic energy K/K0, and numerical dissipation (Eh−E)/E0 against the old-code results in Figures 6–7 and DNS. If the peak enstrophy or the kinetic-energy curve shifts by more than about 2–3% relative to DNS, the 'marginal impact' assertion is refuted and the TGV results cannot be cited for the proposed formulation.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Section 5 states that all Taylor-Green vortex computations were run with an older code that omits the L2 projections of test functions in the convective term (Eq. 2 in §3.4), and asserts without demonstration that the impact is 'marginal at the higher polynomials degrees'. This is load-bearing because the strongest claim, DNS-like precision on coarse meshes, rests entirely on these TGV results, and the omitted L2 projections are the very mechanism that gives the proposed formulation its Reynolds-semi-robust convergence in the convection-dominated regime (§3.4). The paper itself notes that omitting the projections costs a full order of convergence in the convection-dominated regime, while §5 characterizes the loss as only 'half an order'; this internal inconsistency underscores that the old code is not the analyzed scheme. Since TGV at Re=1600 is convection-dominated and the reported diagnostic includes numerical dissipation, the omitted projections could materially change the enstrophy and kinetic-energy evolution shown in Figures 6 and 7. The claim that high k makes the omission negligible is plausible but entirely unquantified.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper proposes an ESDIRK-in-time Hybrid High-Order (HHO) discretization of the incompressible Navier–Stokes equations, based on hybrid velocity and pressure spaces, L2 projections of test functions in the convective term, upwind stabilization, static condensation, p-multilevel preconditioning, and adaptive time stepping. The method is validated on manufactured and analytical solutions (travelling waves, Ethier–Steinmann, LLMS pressure-robustness), on the inviscid double shear layer, and then applied to the Taylor–Green vortex at Re=1600 on tetrahedral meshes with polynomial degrees up to k=9. The main claim is that the formulation can reach DNS-like precision in relevant flow features (kinetic energy and enstrophy) on coarse meshes. The manuscript explicitly states that all TGV computations were performed with an older version of the code that omits the L2 projections in the convective term, and argues that the impact is marginal at high polynomial degree.","tokens_in":32839,"tokens_out":3051,"duration_ms":40345,"significance":"If the claims are substantiated, the paper would make a useful contribution: it combines several attractive properties (pressure-robustness, cellwise mass conservation, Reynolds-semi-robust convergence, high-order implicit time stepping, static condensation) into one HHO framework and provides extensive convergence and robustness evidence. The manufactured-solution tests in Section 4 are broad and mostly match the theoretical rates of [8,9]; the pressure-robustness study and the numerical-dissipation quantification for the double shear layer are valuable. The central limitation is that the turbulence-capability claim rests on TGV runs made with an older code lacking the proposed convective L2 projections, so the DNS-like-precision conclusion is not yet demonstrated for the formulation as presented. No reproducibility package or data archive is mentioned.","major_comments":[{"comment":"The TGV results at Re=1600, which carry the paper's central turbulence claim, were produced with an older code that omits the L2 projections of test functions in the convective term (Eq. (2)), while Section 3.4 states that omitting these projections costs a full order of convergence in the convection-dominated regime; Section 5 instead describes the loss as only 'half an order'. This internal inconsistency, combined with the absence of any quantitative estimate of the projection omission's effect on enstrophy or kinetic energy at high k, makes the assertion that the impact is 'marginal' unsupported. Since TGV at Re=1600 is convection-dominated and the reported numerical dissipation (E_h - E) is exactly the quantity most likely to be affected by the modified convective term, the curves in Figures 6 and 7 cannot be taken as evidence for the proposed formulation until either the key high-k runs are repeated with the updated code or a quantitative sensitivity study (e.g., comparing both variants on a coarse grid at moderate k) is provided.","section":"Section 5, paragraph after Fig. 6; Section 3.4, Eq. (2)"},{"comment":"The evidence presented does not fully support the conclusion that the method 'reaches DNS-like precision in the relevant flow features' on the two tetrahedral grids. Only two grid resolutions are shown, with no error bars or repeated runs; on the coarser grid even the k=9 total-dissipation curve is shifted left relative to DNS, and on the finer grid the kinetic energy still deviates in the enstrophy-decreasing phase. The text itself notes that the maximum relative numerical dissipation below unity implies a 3% discrepancy in kinetic energy relative to a fully resolved computation. The manuscript should state a quantitative accuracy target for 'DNS-like precision' and either demonstrate that this target is met or soften the claim to agreement in trend rather than DNS-like precision.","section":"Section 5, Figs. 6 and 7; Section 6"}],"minor_comments":[{"comment":"The phrase 'would endow losing a full order of convergence' should be corrected to 'would lead to losing' or similar.","section":"Section 3.4, Eq. (2)"},{"comment":"The captions use 'entrophy' instead of 'enstrophy' in several places; please correct the spelling.","section":"Section 5, Figs. 6 and 7 captions and text"},{"comment":"The pressure error in L2 stagnates around 1e-9 for the finest time steps; a brief explanation (e.g., round-off or spatial-error floor) would avoid the appearance of an unexplained saturation.","section":"Section 4.1, Table 1"},{"comment":"The safety factor mu is introduced and set to sqrt(10), but its influence on the tolerance criterion is not discussed; a sentence clarifying why this value is chosen would help.","section":"Section 3.7, Eq. (15)"},{"comment":"For the smallest viscosities and some k values, the observed convergence rates in the convection-dominated regime are noticeably below the claimed k+3/2 (e.g., several k=0 rows in Table 2); a comment on whether this reflects pre-asymptotic behavior or a limitation of the manufactured solution would strengthen the validation.","section":"Section 4.1, Tables 2-6"}],"recommendation":"major_revision","confidential_remarks":"The paper is technically rich and the validation program is impressive, but the central turbulence claim depends on results obtained with a code variant that omits the key convective L2 projections. Given that the authors themselves flag this, the appropriate path is a major revision requiring either reruns with the current formulation for the decisive TGV cases or a quantitative justification of the omission's negligible effect. The scope of the journal is appropriate; no concern about novelty overlap beyond the authors' own prior work."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this paper does what it promises for the first four sections, but the turbulence headline in Section 5 is not yet supported by the evidence. The TGV results were produced by an older code that omits the L2 projections in the convective term, which is the very ingredient that the paper says gives the scheme its Reynolds-semi-robust convergence. The authors disclose this, but then assert without quantification that the impact is \"margical at the higher polynomial degrees\" — and the paper is internally inconsistent, saying in Section 3.4 that omitting the projections costs a full order of convergence, while Section 5 claims only half an order is lost. That is the difference between the analyzed method and the simulated method, and it is load-bearing because the central claim rests entirely on those TGV runs.\n\nWhat is genuinely good: the convergence study is extensive and convincing. The manufactured-solution tests in 2D and 3D cover a wide range of Reynolds numbers, and the observed rates match the theory from the authors' prior work. The pressure-robustness tests, including the non-simplicial-mesh remark, are a real strength. The double shear layer results are also useful, and the static condensation plus p-multigrid strategy is sensible. The paper is honest about many limitations, which matters.\n\nThe soft spots, in proportion: the TGV evidence has only two coarse grids and no error bars, which is fine for a first look but not enough for the strong conclusion. The enstrophy decay phase is visibly off even at k=8, and the kinetic energy discrepancy is about 3% at the finest runs. The conclusion's phrase \"DNS-like precision\" overstates what Figures 6 and 7 actually show. The paper also does not compare against other modern implicit LES solvers, so the performance claim is not contextualized. No code or data are provided, which makes it harder to verify the old-code caveat.\n\nWho this is for: researchers working on high-order methods for incompressible flows, especially HHO and HDG communities. They will find the validation useful and the turbulence section a promising but unfinished story. It deserves a serious referee, but the referee should insist on either re-running TGV with the current code or providing a quantified estimate of the projection's effect on the diagnostics, and should ask for at least one comparison with another implicit LES approach before publication.","headline":"Solid numerical validation of an ESDIRK-HHO scheme, but the Taylor-Green turbulence claim is not yet demonstrated because the runs used an older code without the proposed convective L2 projections.","tokens_in":33395,"tokens_out":2186,"would_cite":false,"duration_ms":26782,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65M60","65N30","76M10","76D05"],"pacs":["47.11.Fg","47.27.E-"],"model":"deepseek-v4-flash","headline":"The paper claims that a hybrid velocity-pressure HHO discretization combined with ESDIRK time stepping reproduces DNS-level kinetic energy and enstrophy for the Taylor-Green vortex at Re=1600 on coarse tetrahedral meshes.","keywords":["Hybrid High-Order methods","incompressible Navier-Stokes","ESDIRK time integration","pressure-robustness","static condensation","Taylor-Green vortex","turbulence modelling","high-order methods"],"falsifier":"Rerun the Taylor-Green vortex at Reynolds 1600 on the finer $24\\times 8^3$ tetrahedral mesh at polynomial degree 8 using the current code that includes the $L^2$ projections in the convective term, and compare the relative enstrophy $E/E_0$ and kinetic energy $K/K_0$ against the DNS reference; if the discrepancies grow beyond the reported 3% kinetic-energy error at peak dissipation, the marginal-impact assumption fails.","tokens_in":32433,"feed_emoji":"🌀","tokens_out":11675,"duration_ms":113632,"temperature":0.7,"pith_summary":"The paper proposes a Hybrid High-Order (HHO) discretization of the incompressible Navier-Stokes equations in which both velocity and pressure are hybrid unknowns living on cells and faces, coupled with ESDIRK implicit time marching up to fifth order. It aims to establish that this combination is pressure-robust, conserves mass cell-by-cell to machine precision, remains stable in the inviscid limit, and can be solved efficiently at high polynomial degree thanks to static condensation and p-multilevel preconditioners. The central numerical evidence is the Taylor-Green vortex at Reynolds 1600, where on coarse tetrahedral meshes with polynomial degree up to 9 the method reproduces DNS reference values for kinetic energy and enstrophy closely enough that the authors describe the precision as DNS-like. A sympathetic reader would care because this points to a practical implicit high-order route toward under-resolved turbulence simulation without relying on explicit turbulence models.","feed_headline":"High-order hybrid scheme reaches DNS-like turbulence on coarse meshes","feed_subtitle":"Taylor-Green vortex at Reynolds 1600 matched by degree-9 ESDIRK-HHO with a fraction of DNS degrees of freedom.","key_machinery":"The central mechanism is full hybridization: the velocity space $V_h^k$ uses $P^{k+1}$ polynomials in cells and $P^k$ on faces, while the pressure space $Q_h^{k+1}$ uses $P^k$ in cells and $P^{k+1}$ on faces, with gradient reconstruction operators $G_T^k$ and $g_T^{k+1}$ tying the two. This structure allows static condensation of both cell velocity and cell pressure, so the global algebraic system contains only face unknowns, and p-multilevel preconditioners operate on this condensed system. The convective term uses upwind fluxes with $L^2$ projections of test functions to retain Reynolds semi-robustness, a time-derivative stabilization term handles the vanishing-viscosity limit, and embedded ESDIRK error estimators drive local time-step adaptation.","core_discovery":"The paper's central claim is that its ESDIRK-HHO formulation—hybrid velocity and pressure spaces, an upwind-stabilized convective term with $L^2$ projections of test functions, and high-order ESDIRK time stepping with local time-step adaptation—delivers pressure-robust, pointwise divergence-free solutions with controlled numerical dissipation up to the inviscid limit. In the Taylor-Green vortex at Reynolds 1600 this combination reproduces the DNS reference evolution of kinetic energy and enstrophy on coarse tetrahedral meshes using polynomial degrees 1 through 9 (up to 8 on the finer mesh), with numerical dissipation decreasing systematically as the polynomial degree and mesh resolution increase. The authors note that the TGV computations were carried out with an older code version that omits the $L^2$ projections in the convective term, and they assert that the impact of this omission is marginal at high polynomial degree.","pith_inferences":["If the marginal-impact assumption about the older TGV code holds, then the updated scheme with $L^2$ projections should gain half an order of convergence in the convection-dominated regime, so the reported TGV results are likely conservative rather than optimistic.","The same hybrid framework could be carried to wall-bounded and non-periodic turbulent flows; the planar-face requirement for static condensation favors tetrahedral meshes, suggesting skeleton-based adaptive mesh refinement as a natural extension.","The systematic reduction of numerical dissipation with polynomial degree suggests that p-refinement alone, without mesh refinement, may be sufficient to resolve the dissipation peak for canonical flows, a trend that runs beyond degree 9 on finer meshes could test directly.","Hybridizing the pressure opens a route to skeleton-based a posteriori error estimation and goal-oriented time-step control that the paper does not explore."],"forward_implications":["The method offers an implicit, fully discrete high-order approach to under-resolved turbulence simulation on coarse meshes, with memory cost dominated by face unknowns rather than cell unknowns.","Cell-by-cell mass conservation to machine precision removes the pressure-velocity coupling errors that typically limit long-time incompressible simulations.","Pressure-robustness means large irrotational body forces affect only the pressure field and do not contaminate the velocity approximation, which matters for buoyancy and rotating flows.","Robustness in the inviscid limit together with measurable, decreasing numerical dissipation makes the scheme usable at very high Reynolds numbers without added artificial viscosity.","Local time-step adaptation with embedded error estimators automatically controls temporal accuracy even when the required step size varies by two orders of magnitude over a simulation."],"supporting_citations":[{"why":"Supplies the base HHO scheme for the incompressible Navier-Stokes and Euler equations that this work extends with ESDIRK time stepping.","marker":"[3]"},{"why":"Provides the stability, convergence, and pressure-robustness theory for hybrid velocity and pressure spaces on which the formulation relies.","marker":"[8]"},{"why":"Introduces the Reynolds-semi-robust variant with $L^2$ projections of test functions in the convective term that the updated scheme follows.","marker":"[9]"},{"why":"Supplies the DNS reference data for the Taylor-Green vortex at Reynolds 1600 against which the method's turbulence results are compared.","marker":"[39]"},{"why":"Defines the Taylor-Green vortex problem and the energy-cascade behavior that makes it the central turbulence test case.","marker":"[21]"},{"why":"Provides the exact three-dimensional solution used to measure the scheme's convergence rates in 3D.","marker":"[19]"},{"why":"Supplies the double shear layer test case used to quantify numerical dissipation and inviscid-limit robustness.","marker":"[36]"},{"why":"Supplies the manufactured Stokes solution used to demonstrate pressure-robustness across viscosities.","marker":"[20]"},{"why":"Provides the ESDIRK coefficients for the high-order implicit time integrators used in the scheme.","marker":"[24]"}],"fun_headline_variants":["HHO-ESDIRK matches DNS Taylor-Green on coarse meshes","Pressure-robust HHO reproduces TGV at Re 1600","Hybrid high-order scheme reproduces TGV at Re 1600","Degree-9 HHO nails Taylor-Green without DNS cost"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The Taylor-Green vortex results were produced with an older version of the code that omits the $L^2$ projections of test functions in the convective term, and the paper assumes, without demonstrating, that this omission has only a marginal effect at high polynomial degree.","fun_headline_variants_meta":{"raw":{"variants":["HHO-ESDIRK matches DNS Taylor-Green on coarse meshes","Pressure-robust HHO reproduces TGV at Re 1600","Hybrid high-order scheme reproduces TGV at Re 1600","Degree-9 HHO nails Taylor-Green without DNS cost"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000483,"raw_usage":{"total_tokens":2374,"prompt_tokens":921,"completion_tokens":1453,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":537,"completion_tokens_details":{"reasoning_tokens":1374}},"tokens_in":537,"tokens_out":1453,"duration_ms":11720,"temperature":1.0,"reasoning_tokens":1374,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T13:06:05.687394+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Rerun the Taylor-Green vortex at Reynolds 1600 on the finer $24\\times 8^3$ tetrahedral mesh at polynomial degree 8 using the current code that includes the $L^2$ projections in the convective term, and compare the relative enstrophy $E/E_0$ and kinetic energy $K/K_0$ against the DNS reference; if the discrepancies grow beyond the reported 3% kinetic-energy error at peak dissipation, the marginal-impact assumption fails.","supporting_citations":[{"cited_title":"Botti, F","cited_arxiv_id":null,"evidence_quote":"Supplies the base HHO scheme for the incompressible Navier-Stokes and Euler equations that this work extends with ESDIRK time stepping."},{"cited_title":"Botti, M","cited_arxiv_id":null,"evidence_quote":"Provides the stability, convergence, and pressure-robustness theory for hybrid velocity and pressure spaces on which the formulation relies."},{"cited_title":"Beir˜ ao da Veiga, D","cited_arxiv_id":null,"evidence_quote":"Introduces the Reynolds-semi-robust variant with $L^2$ projections of test functions in the convective term that the updated scheme follows."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the DNS reference data for the Taylor-Green vortex at Reynolds 1600 against which the method's turbulence results are compared."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Provides the exact three-dimensional solution used to measure the scheme's convergence rates in 3D."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the double shear layer test case used to quantify numerical dissipation and inviscid-limit robustness."}],"review_version":1}