{"id":"ac6b1eda-d3d0-469e-b93d-25966c34f0aa","arxiv_id":"2602.08639","paper_version":4,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":6,"one_line_summary":"With proper limiting, DG and FEM schemes can conserve mass, preserve bounds, and dissipate energy for the Cahn-Hilliard-Navier-Stokes system with degenerate mobility, with SWIP-L the most robust in the tests.","lead":"This paper compares numerical schemes for simulating two-phase fluids with diffuse interfaces, focusing on schemes that preserve mass, keep the phase-field inside physical bounds, and dissipate energy. It introduces refined discontinuous-Galerkin variants and a limiting procedure, then tests them on benchmarks including a rising bubble.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Theorem 3.3's energy-dissipation proof ignores the DG advection term; the claimed identity (68) cannot hold in the advective CHNS setting, so the energy-dissipation column of Table 5 lacks theoretical support.","rationale":"The paper's central claim is that the limited DG schemes, especially SWIP-L, are structure-preserving in the CHNS setting. Structure preservation for energy rests on Theorem 3.3. That theorem is stated for the full DG formulation including advection, but the proof ignores c(u,ψ,υ). The resulting identity (68) cannot be derived from the weak form; the advective term enters the energy balance. This is a concrete mathematical error, not merely an unproven assumption. The cell-average concern raised by the reader is real but is explicitly acknowledged in Remark 3.12 and is an empirical condition; the energy theorem is presented as a proven result yet is not valid as stated. A single numerical check—computing the discrete energy residual in an advective CH test—would settle whether the theorem's omission matters. If the residual matches c(u,ψ,υ), the theorem must be revised, for example by proving cancellation with the kinetic energy in the coupled system or by restricting the claim to the u=0 case. The reader's conditional verdict remains appropriate; no change in verdict is needed, but the revision requirements should include correcting Theorem 3.3.","tokens_in":33088,"tokens_out":17553,"duration_ms":180807,"concrete_test":"Implement the SIPG/SWIP scheme (Eqs. (64)-(65)) for the advective Cahn-Hilliard equation on a periodic domain with constant divergence-free velocity u and non-symmetric initial data (e.g., a translated bubble). At each time step, compute the discrete phase-field energy E[ψ_h] from Eq. (67) and the right-hand side of Eq. (68). If their difference is nonzero and tracks c(u,ψ_h,υ_h) from Eq. (63), Theorem 3.3 is false. Repeat with the limiter active to test whether the limited scheme can ever increase E.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Theorem 3.3 ('Recovery of energy dissipation') is the paper's theoretical basis for the 'energy dissipative' entries in Table 5. The proof asserts ∂t E[ψ_h] = -(CnWePe)^{-1} b(M(ψ_h),υ_h,υ_h) by coercivity, but the DG formulation (64)-(65) contains the discrete advection term c(u,ψ_h,φ) from Eq. (63). Setting φ=υ_h in the phase-field equation gives ⟨∂tψ_h,υ_h⟩ = c(u,ψ_h,υ_h) - Pe^{-1} b(...). Consequently the energy identity for E[ψ_h] in Eq. (67) must include c(u,ψ_h,υ_h) unless this term vanishes for divergence-free u. The upwind DG advection operator is not skew-symmetric for the nonlinear energy; c(u,ψ_h,υ_h) is generally nonzero. Hence Eq. (68) is not a valid equality in the advective (CHNS) case. This is not a missing detail: phase-field energy alone is not monotonically decreasing in the coupled system; only the total energy is. The theorem therefore cannot justify the energy-dissipation claims, and the limiter used in SWIP-L further modifies ψ_h after each step, moving the solution away from any discrete energy law.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper presents a systematic numerical comparison of structure-preserving discretizations for the Cahn-Hilliard-Navier-Stokes (CHNS) system with degenerate mobility, in the setting of conforming adaptive mesh refinement. The schemes considered are continuous FEM (with and without a cut-off post-processing, with and without a scaling limiter), DG-SIPG and DG-SWIP (both with and without a scaling limiter), and the upwind ASU scheme. The authors provide a DG formulation with a weighted-SWIP variant and a scaling limiter, and they report convergence rates, mass conservation, phase-field boundedness, and energy dissipation for pure CH and coupled CHNS test cases, including a rising-bubble benchmark. The central recommendation is that only the limited schemes (FEM-L, SIPG-L, SWIP-L) and ASU preserve both mass and bounds, while the standard unlimited schemes violate the bounds and the cut-off FEM loses mass at O(10^-4). SWIP-L is highlighted as the most robust among the DG variants.","tokens_in":33468,"tokens_out":9864,"duration_ms":107857,"significance":"If the comparison and the associated software are reliable, this is a useful practical contribution: it gives a head-to-head assessment of recently proposed structure-preserving schemes in a UFL/Dune-Fem framework, with reproducible code, adaptive-mesh benchmarks, and a clear recommendation for users who need boundedness and mass conservation in multiphase flow simulations. The numerical study is carefully structured, with convergence-rate tables, mass/energy/bounds tables, and adaptive-grid experiments. The paper also explicitly acknowledges some open theoretical points (e.g., ASU energy dissipation in the decoupled CHNS setting, and the cell-average assumption for the limiter), which is a strength. However, the theoretical support for the energy-dissipation claims in the advective case is not correct as stated, and the boundedness proof for the limited schemes is explicitly conditional on an unproven cell-average property. These issues affect the interpretation of the headline claims, although the numerical evidence itself appears solid.","major_comments":[{"comment":"The claimed energy identity is not valid for the advective DG formulation. The phase-field equation in Eqs. (64)–(65) contains the discrete advection term −c(u,ψ_h,φ) from Eq. (63). Setting φ=υ_h in that equation gives ⟨∂tψ_h,υ_h⟩ = c(u,ψ_h,υ_h) − Pe^{-1} b(M(ψ_h),υ_h,υ_h). The energy derivative ∂tE[ψ_h] therefore contains the c-term, which does not vanish in general for the upwind DG discretization with a solenoidal u. The proof of Theorem 3.3 simply asserts the dissipation without deriving this identity, and Eq. (68) is not a consequence of the discrete scheme. In the continuous CHNS setting only the total energy (24), including the kinetic energy, is dissipated; the phase-field energy alone is not monotone under advection. Consequently, Theorem 3.3 cannot justify the 'energy dissipative' entries for the coupled runs in Table 5. The numerical energy plots in Sections 4.4–4.5 may well b","section":"§3.3, Theorem 3.3 and Eq. (68)"},{"comment":"The boundedness of FEM-L, SIPG-L, and SWIP-L is conditional on the cell averages of the phase-field remaining in [-1,1] after each time step. This is stated explicitly in Corollary 3.6.1 ('provided that the FEM phase-field for each element K satisfies 1/|K|∫_K ψ̃_h dx ∈ [-1,1]') and in Remark 3.12, where the authors say this condition is not proven and only 'no violation ... has been encountered' in the experiments. The scaling limiter in Eq. (83) scales around the mean value, so if the mean itself lies outside the interval the limiter cannot restore the bounds while preserving mass. Thus Table 5's 'boundedness' column for these schemes is not a proven structure-preservation property, only an observed one. Since boundedness is one of the three advertised structure-preserving properties and a central part of the paper's recommendation, this conditional should be prominently qualified in t","section":"§3.5, Corollary 3.6.1 and Remark 3.12"},{"comment":"The coercivity theorem requires the mobility-dependent penalty Λ_e(M(ψ_h)) to be at least of order λ⋆, the global contrast defined in Eq. (47), cf. Eq. (52). With the regularized mobility M_δ (the text uses δ=10^{-20} in §3.2 and δ=10^{-14} in the UFL code), the contrast λ⋆ can be very large near ψ=±1, since the denominator min_K M(ψ_h) is of order δ. However, Remark 3.3 fixes Λ_e=5 and η=k(k+d−1)(θ+1)^2/4 for all computations. These values do not in general satisfy the sufficient inequality (52) unless λ⋆ is bounded by a small mesh-independent constant, which is not established. The authors note the dependence on smoothness of ψ_h but do not verify it for their examples. This is not a fatal flaw for a numerical comparison, but the paper presents Theorem 3.2 as a theoretical foundation for the SWIP-L/SIPG-L schemes; the connection between the theorem and the actual computational paramete","section":"§3.2, Theorem 3.2 and Remark 3.3"}],"minor_comments":[{"comment":"Please check the sign convention of the advection form c(u,ψ_h,φ). Eq. (64) writes ⟨∂tψ_h,v⟩ − c(u,ψ_h,φ) + ..., while the continuous equation (17) contains +∇·(uψ). If c is defined by Eq. (63) as ∫_Th u·∇φ ψ_h dx − Σ(...)[φ]ds, then −c corresponds to the standard weak form only after an integration by parts; the sign conventions should be stated explicitly and made consistent in all formulae and UFL listings.","section":"Eq. (63)–(65)"},{"comment":"The DG energy lifting in Eq. (67) writes the last term as '−(1+θ){∇ψ_h·n+}[ψ]' without the surface measure ds. Please add the missing measure and check the boundary terms. Also clarify how this discrete energy relates to the continuous total energy (24) in the coupled CHNS runs shown in Figs. 12, 15, and 18.","section":"Eq. (67)"},{"comment":"In Table 5 the footnote for FEM-L is not referenced or explained in the text. Also, in Section 4.5.2 the statement 'The simulations using Ex. 4.4 did not converge using the ASU scheme' is a notable negative result that is not discussed. A brief explanation (e.g., conditioning, coupling, or time-step restrictions) would help the reader interpret the comparison.","section":"Table 5 and Section 4.5.2"},{"comment":"Refs. [13] and [39] appear to be the same article (one published, one preprint), with the same DOI in the bibliography. Please consolidate the citations and avoid duplicate entries.","section":"References [13] and [39]"}],"recommendation":"major_revision","confidential_remarks":"The numerical study is well executed and the software availability is a strong point. The main problem is that one of the three advertised structure-preserving properties—energy dissipation—is supported by a theorem (Thm. 3.3) that is plainly incorrect in the advective case, and the boundedness property is explicitly conditional on an unproven cell-average assumption. These are fixable within the scope of the paper by either proving the correct total-energy statement, clearly restricting the theory to the non-advective case, and qualifying the boundedness claims. I recommend major revision rather than rejection because the central benchmark conclusions are likely unaffected, but the theoretical framing must be corrected."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Bottom line: this is a useful numerical paper whose central recommendation is probably right, but one of its theoretical supports is broken. The comparison is honest and broad: FEM, FEM-L, SIPG(-L), SWIP(-L), and ASU across pure CH, rotating bubbles, and the Hysing-Turek rising bubble, with adaptive meshes. The new pieces — the SWIP weighted interior-penalty form in Eq. (59), the FEM-L mass-lumped limiter pipeline, and the careful mass-vs-solver-tolerance table — are real contributions, and the conclusion that limited schemes beat FEM-C is backed by the experiments.\n\nThe soft spot is Theorem 3.3. The proof of energy dissipation does not account for the DG advection term in Eq. (63), and in the coupled CHNS problem the phase-field energy alone is not monotone; only the total energy is. So Eq. (68) cannot hold as stated for the advective case, and the \"energy dissipative\" entries in Table 5 for the coupled runs are not supported by that theorem. The stress-test note is right. This needs either a proof that the discrete advection term is harmless for the chosen energy, a switch to the total discrete energy, or a clear statement that energy dissipation in the coupled runs is observed numerically, not proven. It is not a reason to reject the numerical comparison, but it is a reason not to call the schemes structure-preserving in that strong sense.\n\nSecond soft spot: bound preservation for the limited schemes is conditional on the cell average lying in [-1,1] (Cor. 3.6.1 and Remark 3.12). The authors say they never saw a violation; that is fine, but it should be reported as an empirical property for these tests, not as a theorem. The ASU scheme's non-convergence in the rising bubble case is reported but not discussed; a few sentences would help since ASU is a baseline competitor. Minor: the proofs of Theorems 2.1 and 2.2 are sketches, but those are standard and not load-bearing.\n\nThe math and citation pattern look okay; the DG coercivity proof is more detailed than the rest, and the UFL formulations are a nice reproducibility gesture, though no code archive is linked.\n\nWho this is for: numerical PDE people choosing a phase-field scheme, and DG practitioners. It deserves serious refereeing. I would send it out, with a request to fix Theorem 3.3 and to soften the structure-preserving claims to match what is actually proven.","headline":"Useful benchmark comparison and two new scheme variants, but the energy-dissipation theorem ignores the advection term and should be fixed before the paper is treated as a theoretical reference.","tokens_in":33947,"tokens_out":3859,"would_cite":true,"duration_ms":47016,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65M60","65M50","76T99","35Q35"],"pacs":[],"model":"deepseek-v4-flash","headline":"This paper establishes that, for the Cahn–Hilliard–Navier–Stokes system with degenerate mobility, an element-wise scaling limiter applied to discontinuous Galerkin or finite element discretizations yields schemes that simultaneously conserv","keywords":["Cahn-Hilliard-Navier-Stokes","structure-preserving schemes","degenerate mobility","bound preservation","mass conservation","energy dissipation","discontinuous Galerkin","adaptive mesh refinement"],"falsifier":"Run the SWIP-L scheme on a coarse adaptive mesh with a strong advection field, starting from initial data with cell averages inside [−1, 1], and check whether any pre-limiter cell average exits [−1, 1] after the first IMEX step. If it does, the scaling limiter in Eq. (83) cannot restore the bound and the boundedness claim fails. Alternatively, test whether the discrete energy in Eq. (67) is monotonically non-increasing for a range of penalty parameters; a violation would contradict Theorem 3.3.","tokens_in":32928,"feed_emoji":"💧","tokens_out":8342,"duration_ms":77549,"temperature":0.7,"pith_summary":"This paper asks which numerical schemes can simulate the Cahn–Hilliard–Navier–Stokes (CHNS) equations for two-phase flow while preserving three physical properties: the phase-field mass stays constant, the phase-field variable stays within [−1, 1], and the total free energy never increases. The authors compare existing schemes and propose modifications: they apply an element-wise scaling limiter to standard discontinuous Galerkin (SIPG, SWIP) and continuous Galerkin (FEM) discretizations, and they introduce a weighted DG treatment, SWIP, that reduces mobility contrast at element faces. They prove discrete energy dissipation for the DG schemes and give a conditional bound-preservation result for the limited schemes. In all tested cases, the limited schemes (FEM-L, SIPG-L, SWIP-L) and the existing ASU scheme conserve mass, stay in bounds, and dissipate energy, while unlimited FEM/SIPG/SWIP violate the bounds and the cut-off FEM-C loses mass at the 10⁻⁴ level. The authors conclude that SWIP-L is the most reliable of the schemes, especially with adaptive mesh refinement.","feed_headline":"Limiter makes Cahn-Hilliard schemes conserve mass, bounds, energy","feed_subtitle":"An element-wise scaling limiter turns DG and FEM solvers into structure-preserving schemes for two-phase flow.","key_machinery":"The element-wise scaling limiter: on each element K, the phase-field polynomial ψ_h is scaled around its cell mean ψ̄_K by a factor α_K chosen so that all quadrature points lie within [ψ_min, ψ_max]=[−1, 1]. Because the cell mean is unchanged, the limiter is mass-conservative by construction, and it enforces the bound provided the cell averages are already inside the interval. Supporting this, the paper proves coercivity of a weighted interior-penalty bilinear form; the harmonic-average weighting (SWIP) is the key modification that lowers the mobility contrast and improves conditioning of the linear systems.","core_discovery":"The central claim is that structure preservation for the CHNS system with degenerate mobility M(ψ)=1−ψ² is achievable without sacrificing accuracy or mass conservation. Applying an element-wise scaling limiter to either DG (SIPG/SWIP) or standard FEM discretizations yields schemes that conserve phase-field mass, keep the phase field in [−1, 1], and dissipate the discrete free energy in every test performed. The paper further claims that the proposed SWIP variant, which weights the numerical fluxes using harmonic averages of the degenerate mobility, reduces the mobility contrast over element faces, leading to a better-conditioned and more reliable scheme than standard SIPG while preserving th","pith_inferences":["The discrete bound-preservation theorem is conditional: it holds only if the cell average of the un-limited phase-field lies in [−1, 1]. A rigorous proof for the regularized degenerate-mobility scheme would close the gap the paper leaves open.","The harmonic-average weighting idea could transfer to other degenerate-mobility or strongly heterogeneous diffusion problems, where face-averaged coefficients create large contrast and degrade conditioning.","The mass drift observed for FEM-L in the coupled CHNS tests may stem from the interaction between the limiter and the non-divergence-free Taylor–Hood velocity field; a divergence-free velocity discretization might remove this drift.","A systematic parameter sweep in the penalty parameter η and mobility regularization δ could map the actual coercivity margin of SWIP and reveal when the scheme ceases to be energy dissipative."],"forward_implications":["The limited schemes FEM-L, SIPG-L, and SWIP-L, together with the ASU scheme, achieve all three structure-preserving properties (mass conservation, bound preservation, energy dissipation) in every numerical test; unlimited FEM, SIPG, and SWIP violate the bounds.","The scaling limiter avoids the O(10⁻⁴) mass loss of the cut-off method (FEM-C), so applications that need in-bounds phase fields can use the limited schemes instead of artificial post-processing.","The DG schemes converge at O(h²) in L2 for piecewise-linear elements, while the ASU scheme converges at O(h); the DG schemes achieve comparable accuracy to FEM at coarser meshes, reducing computational time.","The SWIP-L scheme is better conditioned than SIPG-L, making it the recommended choice among the DG variants, especially in combination with adaptive mesh refinement.","The structure-preserving properties survive conforming adaptive mesh refinement with the indicator function (1−ψ_h²)⁴, which focuses refinement at the interface."],"fun_headline_variants":["Limiter enforces mass, bounds, energy in two-phase flow solvers","Element-wise scaling limiter makes CHNS schemes structure-preserving","New limiter keeps phase field in bounds while conserving mass and energy","SWIP variant beats SIPG for degenerate-mobility CHNS with adaptive mesh","Harmonic mobility weighting improves structure-preserving CHNS scheme"],"cache_read_input_tokens":2304,"weakest_assumption_plain":"The scaling limiter can only enforce ψ∈[−1, 1] if the cell average of the computed phase-field already lies in [−1, 1]; the paper assumes, based on numerical observation but without a discrete proof, that the degenerate-mobility scheme keeps cell averages in bounds.","fun_headline_variants_meta":{"raw":{"variants":["Limiter enforces mass, bounds, energy in two-phase flow solvers","Element-wise scaling limiter makes CHNS schemes structure-preserving","New limiter keeps phase field in bounds while conserving mass and energy","SWIP variant beats SIPG for degenerate-mobility CHNS with adaptive mesh","Harmonic mobility weighting improves structure-preserving CHNS scheme"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000719,"raw_usage":{"total_tokens":3035,"prompt_tokens":684,"completion_tokens":2351,"prompt_tokens_details":{"cached_tokens":256},"prompt_cache_hit_tokens":256,"prompt_cache_miss_tokens":428,"completion_tokens_details":{"reasoning_tokens":2259}},"tokens_in":428,"tokens_out":2351,"duration_ms":16860,"temperature":1.0,"reasoning_tokens":2259,"cache_read_input_tokens":256,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-03T03:12:23.077557+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Run the SWIP-L scheme on a coarse adaptive mesh with a strong advection field, starting from initial data with cell averages inside [−1, 1], and check whether any pre-limiter cell average exits [−1, 1] after the first IMEX step. If it does, the scaling limiter in Eq. (83) cannot restore the bound and the boundedness claim fails. Alternatively, test whether the discrete energy in Eq. (67) is monotonically non-increasing for a range of penalty parameters; a violation would contradict Theorem 3.3.","supporting_citations":[],"review_version":1}