{"id":"9df01f2f-01f7-40c9-ab3d-2c968f50f3ef","arxiv_id":"2411.09485","paper_version":2,"verdict":"ACCEPT","confidence":"MODERATE","novelty_score":7.0,"correctness_risk":"low","formal_verification":"none","parameter_count":0,"one_line_summary":"An exact quadrature for rational polynomial integrals on 2D triangles enables exact MATLAB implementations of singular Zienkiewicz and Guzmán-Neilan elements.","lead":"The paper derives an exact recursive formula for integrating a class of rational functions on triangles, and uses it to implement two rational finite elements in MATLAB. This removes a long-standing practical barrier to using structure-preserving elements that enforce exact divergence-free or H2 conformity.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Algorithm 2's exact integration is sound; the m=0 integrability threshold (3.10) is correct, though Lemma 3.2 is overbroad for m≥1 and its sharpness proof is not fully demonstrated.","rationale":"The paper's core contribution is a recursive exact quadrature for rational functions of barycentric coordinates. The recursion depends on algebraic identities (Lemma 3.8), base formulas (Lemma 3.4 and Algorithm 1), and an exact finiteness test (3.10). I checked the recursion algebra and the termination measure; both hold. The only fragile-looking spot is Lemma 3.2's sharpness. Replacing its asserted cancellation argument by a direct vertex-asymptotic computation gives (3.10) exactly for m=0, which is the only case Algorithm 2 uses. The lemma's extension to m≥1 is indeed false (λ1 counterexample), so the paper contains an overstatement; but this is outside the integration algorithm's dependency chain. Numerical experiments are illustrative rather than central; they do not affect the acceptance. Therefore the reader's ACCEPT verdict stands.","tokens_in":34340,"tokens_out":26133,"duration_ms":233321,"concrete_test":"Perform a direct asymptotic expansion of R_α^β near each vertex for all α,β with entries 0..3, verifying I(α,β)<∞ iff min_i(|α|−α_i−β_i) > −2; this independently confirms (3.10) without relying on Lemma 3.2's cancellation claim. If the check passes, the reader's flagged concern is resolved and Algorithm 2's acceptance criterion is sound.","verdict_should_be":"UNCHANGED","load_bearing_attack":"No significant objection to the central claim. Algorithm 2's recursion identities and the finiteness checks used in Proposition 3.9 are internally consistent. The reader's flagged assumption, the sharp threshold (3.10), can be verified directly: near vertex v_i the integrand behaves like r^{e_i} with e_i = |α| − α_i − β_i, so I(α,β) < ∞ iff e_i > −2 for every i, which is exactly ∥α+β∥∞ ≤ |α|+1. This argument does not need the 'linearly independent singularities' assertion in Lemma 3.2. That assertion is in fact too strong: Lemma 3.2 as stated for m ≥ 1 is false, e.g. R = λ_1 with α=(1,0,0), β=0 satisfies λ_1 ∈ W^{1,2} although |α|−∥α+β∥∞ = 0 < 1. However, Algorithm 2 and Proposition 3.9 only invoke the m=0 case (3.10), so this overstatement does not touch the exact integration claim. The numerical experiments support the qualitative claims about inexact quadrature; no load-bearing flaw found.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper develops an exact quadrature rule for rational functions of the form λ^α/(1−λ)^β on triangles. The central algorithmic contribution is a recursion (Algorithm 2) that computes the integral mean I(α,β) for all multi-indices with finite integral, with termination guaranteed by a sharp integrability criterion. The method is applied to the rational bubble functions appearing in the singular Zienkiewicz and Guzmán–Neilan finite elements, and MATLAB implementations are provided. Numerical experiments demonstrate that inexact quadrature causes stagnation in biharmonic eigenvalue computations and destroys pressure robustness in the Guzmán–Neilan Stokes discretization.","tokens_in":34558,"tokens_out":15438,"duration_ms":126192,"significance":"If the central claim holds, the paper fills a genuine gap in the literature: no exact quadrature for multivariate non-tensor-product rational functions on triangles was previously available. The derivation is parameter-free and the recursion identities are proven in detail, including a termination argument. The code is publicly available and the numerical experiments are reproducible and directly support the qualitative claims about the effects of inexact quadrature. The main theorem (Algorithm 2) is sound; the only mathematical defect is an overstatement in Lemma 3.2 that does not affect the main result but must be corrected.","major_comments":[{"comment":"The claimed characterization of W^{m,p}-regularity is false for m ≥ 1. A counterexample is R^(1,0,0)_(0,0,0) = λ_1 on the reference triangle: for m = 1 and p = 2, one has |α| − ∥α+β∥∞ = 0, which violates the stated condition |α| − ∥α+β∥∞ > m − 2/p, yet λ_1 ∈ W^{1,2}(T). The sharpness argument invoking 'linearly independent singularities' is not valid in general, because cancellations under differentiation can occur. Since Algorithm 2 and Proposition 3.9 only use the m = 0 case, namely condition (3.10), this error does not affect the exact integration claim; nevertheless, the lemma should be corrected (for example, restricted to m = 0) or replaced by a valid statement for m ≥ 1.","section":"§3.1, Lemma 3.2"}],"minor_comments":[{"comment":"The sorting instruction reads 'Sort ((αj,βj))j=1,2 such that β0 ≤ β1 ≤ β2', but the indices should range over j = 0, 1, 2; please clarify the notation.","section":"§3.2.3, Algorithm 2"},{"comment":"The notation (DF)^{−⊤} is inconsistent with the code in Figure 4, line 6, which uses (DF)^{-1}, and with the preceding derivation Dλref = Dλ DF; the formula should read (DF)^{-1}.","section":"§4, Eq. (4.2)"},{"comment":"The termination proof states that in the final recursion 'α0 or |β| is reduced', but the earlier branches using Proposition 3.9(c) also reduce α0. For completeness, the proof should specify a global decreasing measure, such as α0 + |β|, for all branches of Algorithm 2.","section":"§3.2.3, Theorem 3.10"}],"recommendation":"minor_revision","confidential_remarks":"The false statement in Lemma 3.2 is a real mathematical error, but it is localized and does not touch the main algorithmic result. I am confident the authors can correct it (e.g., by restricting the lemma to the m = 0 case or by stating a correct version for m ≥ 1). The manuscript is otherwise a solid contribution with reproducible code and convincing numerical experiments."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Here's my read on Diening–Storn–Tscherpel. The core contribution is real: Algorithm 2 recursively integrates rational functions R_α^β on triangles exactly, terminating via a clean integrability criterion (3.10), and the authors back it with MATLAB code and numerical experiments. I'm fairly convinced it is the first exact quadrature for this class; the Duffy-transform route in Schneier's paper is inexact, and the tensor-product NURBS approaches do not cover these non-tensor rationals. The recursion identities in Lemma 3.8 are simple but effective, the termination argument in Theorem 3.10 is sound, and the derivation is from first principles—no fitting parameters anywhere.\n\nThe one substantive gap is Lemma 3.2. The sharpness proof for m ≥ 1 leans on a claim about 'linearly independent singularities' that is not demonstrated, and as stated the lemma is actually false for m ≥ 1: take α = (1,0,0), β = 0; then R = λ_1 is in W^{1,2} even though |α| − ∥α+β∥∞ = 0 < 1. The stress-test note is right that this overstatement does not damage the paper, because Algorithm 2 and Proposition 3.9 only invoke the m = 0 case (3.10), which is correct and can be verified directly by checking the exponent near each vertex. The authors should tighten Lemma 3.2 to m = 0 or add the direct proof; it is a minor patch, not a load-bearing flaw.\n\nThe numerical experiments are honest and useful: they show stagnation of eigenvalue error with inexact quadrature and destruction of pressure robustness in the Guzmán–Neilan Stokes discretization. This is a practical finding that justifies the exact rule. The code is available via Code Ocean, though I would have liked a commit hash; that is a nit. Citation patterns look appropriate—Schneier as the inexact baseline, Guzmán–Neilan and Ciarlet for element properties, standard quadrature literature.\n\nThis paper will be most valuable to anyone implementing rational finite elements or working on structure-preserving discretizations; it removes a real practical barrier. I would send it to a serious referee. The central contribution is formally grounded, reproducible, and the numerical claims are independently checkable. My recommendation: accept after a minor revision that fixes the scope of Lemma 3.2.","headline":"Solid, useful paper: exact recursive integration for rational finite elements, with working code; the central claim holds, and the flagged Lemma 3.2 overstatement does not touch the algorithm.","tokens_in":35060,"tokens_out":1722,"would_cite":true,"duration_ms":16503,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65D32","65N30","76D07","76M10"],"pacs":[],"model":"deepseek-v4-flash","headline":"A recursive integration rule computes exact integrals of rational functions on triangles, enabling quadrature-free implementations of rational finite elements.","keywords":["exact quadrature","rational functions","singular Zienkiewicz element","Guzmán-Neilan element","biharmonic equation","Stokes equations","pressure robustness","MATLAB implementation"],"falsifier":"Implement Algorithm 2 and compare its output on all multi-indices with entries up to, say, 3 to a high-precision Duffy-transform quadrature (e.g., 50-point Gauss rules per direction) on the reference triangle; a single mismatch in either the finite/infinite classification or the numerical value beyond rounding error would falsify the recursion or the integrability threshold. More directly, proving that the two terms in identity (3.7) cancel for some $(\\alpha,\\beta)$ would break the sharpness argument of Lemma 3.2.","tokens_in":34149,"feed_emoji":"🧮","tokens_out":8976,"duration_ms":77257,"temperature":0.7,"pith_summary":"The paper establishes the first exact integration rule for the class of rational functions $R_\\alpha^\\beta = \\lambda^\\alpha/(1-\\lambda)^\\beta$ on triangles in 2D, the functions that appear in rational finite element bases. It gives a recursive algorithm that returns the exact integral mean $I(\\alpha,\\beta)$ for every multi-index pair with finite integral, and proves termination. This matters because such rational elements—the singular Zienkiewicz element for biharmonic problems and the lowest-order Guzmán–Neilan element for Stokes flow—achieve $H^2$-conformity or exact divergence-freeness with only 12 degrees of freedom per triangle, but had previously relied on inexact Duffy-transform quadrature. The paper shows numerically that inexact quadrature can stall eigenvalue convergence and destroy pressure robustness, so the exact rule is what makes these elements fully trustworthy.","feed_headline":"Recursion computes exact integrals for rational FEM bases","feed_subtitle":"For singular Zienkiewicz and Guzmán–Neilan elements, this removes quadrature errors that break structure preservation.","key_machinery":"The central object is the rational function $R_\\alpha^\\beta(\\lambda)=\\lambda^\\alpha/(1-\\lambda)^\\beta$ on the 2D simplex, where $\\lambda$ are the barycentric coordinates and $\\alpha,\\beta\\in\\mathbb{N}_0^3$. The argument is carried by three recurrence identities (Lemma 3.8) that express $R_\\alpha^\\beta$ as a linear combination of rational functions with lower total denominator exponent $|\\beta|$ or lower numerator exponent component $\\alpha_0$, together with two closed-form base cases: the polynomial integral formula (Lemma 3.3) and the single-pole formula for $\\beta=(0,0,\\beta_2)$ (Lemma 3.4). Algorithm 2 uses the integrability threshold of Lemma 3.2 as its first check, then routes the recursion through the identities and, when $\\alpha_0=0$, through the one-dimensional closed-form evaluation of Algorithm 1.","core_discovery":"The central claim is that the integral mean of $R_\\alpha^\\beta$ over any triangle can be evaluated exactly by a finite recursion. The recursion uses three algebraic identities that rewrite the integrand into sums of rational functions with smaller parameters, and it is organized by an exact regularity condition: $I(\\alpha,\\beta)<\\infty$ if and only if $\\lVert\\alpha+\\beta\\rVert_\\infty \\le |\\alpha|+1$ (Lemma 3.2). Theorem 3.10 proves Algorithm 2, which implements this recursion, terminates and returns the exact value in all finite cases. As a consequence, the paper provides simple MATLAB implementations of the singular Zienkiewicz and the lowest-order Guzmán–Neilan elements whose system matrices are computed exactly, and its numerical experiments demonstrate that inexact quadrature causes stagnation of the biharmonic eigenvalue error and a loss of pressure robustness in the Stokes discretization, both of which are remedied by the exact integration.","pith_inferences":["The same recursion identities may extend to rational functions in higher dimensions or with denominator factors that are products of general affine functions, since the algebraic decompositions in Lemma 3.8 are dimension-agnostic in form.","If the numerical stagnation observed with inexact quadrature is generic, then any rational finite element method intended for eigenvalue or long-time structure-preserving simulations should use exact integration from the outset.","The sharp integrability threshold is a reusable analytical tool: it provides a simple criterion for designing rational bubbles with prescribed Sobolev regularity on simplices.","Memoizing intermediate values of $I(\\alpha,\\beta)$ and exploiting symmetries, as the performance remark suggests, should make the exact rule fast enough for production-scale use."],"forward_implications":["The singular Zienkiewicz element can be assembled with exactly computed biharmonic stiffness and mass matrices during an offline phase, eliminating quadrature error from the discretization.","The Guzmán–Neilan Stokes element preserves its exact divergence-free constraint when system matrices are formed with the exact rule, so discrete solutions remain exactly divergence-free.","Inexact quadrature is shown to cause stagnation of eigenvalue error on uniform meshes and a reduced convergence order on graded meshes; exact integration removes this barrier.","The provided MATLAB template demonstrates a general pipeline for implementing rational finite elements, so other rational elements can be studied without developing new quadrature."],"supporting_citations":[{"why":"Defines the singular Zienkiewicz element, its local space, and C^1-conformity properties that the paper builds on.","marker":"[Cia02]"},{"why":"Introduces the Guzmán–Neilan element and its exact de Rham complex, the structure preservation that the new quadrature must maintain.","marker":"[GN14b]"},{"why":"The prior inexact Duffy-transform quadrature for rational bubble functions that the paper replaces.","marker":"[Sch15]"},{"why":"Establishes the pressure-robustness test problem and explains why exact divergence-freeness matters for Stokes discretizations.","marker":"[JLMNR17]"},{"why":"Provides the optimal convergence rate for the biharmonic eigenvalue scheme used to interpret the numerical results.","marker":"[CP23]"},{"why":"Supplies the closed-form polynomial integration formula that serves as a base case in the recursion.","marker":"[GM78]"}],"fun_headline_variants":["Exact integrals for rational finite elements via recursion","Recursive exact integration for Zienkiewicz and Guzman-Neilan FEM","Exact rational FEM integrals: recursion instead of quadrature","Exact integration for singular FEM bases, no quadrature error","Recursive exact integrals for singular Zienkiewicz and Guzman-Neilan"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The paper's recursion relies on the exactness of its integrability test: it assumes that no exceptional rational function with $\\lVert\\alpha+\\beta\\rVert_\\infty = |\\alpha|+1$ is actually integrable because its two singular pieces cancel. If such a cancellation existed, the algorithm would wrongly return infinity for a finite integral.","fun_headline_variants_meta":{"raw":{"variants":["Exact integrals for rational finite elements via recursion","Recursive exact integration for Zienkiewicz and Guzman-Neilan FEM","Exact rational FEM integrals: recursion instead of quadrature","Exact integration for singular FEM bases, no quadrature error","Recursive exact integrals for singular Zienkiewicz and Guzman-Neilan"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000674,"raw_usage":{"total_tokens":2984,"prompt_tokens":780,"completion_tokens":2204,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":396,"completion_tokens_details":{"reasoning_tokens":2114}},"tokens_in":396,"tokens_out":2204,"duration_ms":15608,"temperature":1.0,"reasoning_tokens":2114,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-12T20:36:35.881677+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Implement Algorithm 2 and compare its output on all multi-indices with entries up to, say, 3 to a high-precision Duffy-transform quadrature (e.g., 50-point Gauss rules per direction) on the reference triangle; a single mismatch in either the finite/infinite classification or the numerical value beyond rounding error would falsify the recursion or the integrability threshold. More directly, proving that the two terms in identity (3.7) cancel for some $(\\alpha,\\beta)$ would break the sharpness argument of Lemma 3.2.","supporting_citations":[],"review_version":1}