{"id":"5cfca797-b59b-4564-ab1e-f5df3b534502","arxiv_id":"2502.10069","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":0,"one_line_summary":"A new bound-preserving PAMPA scheme for hyperbolic conservation laws on polygonal meshes is presented, with a proof outline and an eigenvalue-based coefficient formula for the Euler equations.","lead":"This paper describes a new numerical method for simulating wave and fluid dynamics on polygonal grids, and argues that it provably keeps solutions inside physical limits. Readers in simulation-heavy fields might care because such methods could make high-accuracy flow calculations on complex shapes more reliable.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Euler positivity depends on an unproved and apparently misprinted spectral-radius formula in §6; until it is verified, the 'provably bound preserving for Euler' claim is unsupported.","rationale":"The reader's weakest_assumption identifies exactly the unproved closed-form eigenvalue formula in §6 as the load-bearing point for the Euler positivity claim. Reading the paper in good faith, the reduction to a Rayleigh quotient and spectral radius is plausible for a symmetric positive definite A, but the formula itself is presented without derivation and the printed matrices do not obviously correspond to the quadratic forms they are meant to represent. This is the single most load-bearing concern because the entire bound-preserving property for the Euler system rests on the computed θ values. I considered whether the lack of polygonal-mesh numerical tests or the summary nature of the paper is more important, but those are weaknesses in presentation and empirical support, not in the correctness of the central proof. The unverified spectral-radius formula is a genuine correctness risk: if it fails, the scheme may produce negative density or energy. Because the reader already flagged this exact gap and issued a conditional verdict, my stress-test does not change the verdict; it sharpens the concern with a concrete algebraic inconsistency and proposes a direct numerical check that can settle it. If the test fails, the appropriate verdict would be REJECT or at least a demand for a corrected derivation before any acceptance.","tokens_in":11315,"tokens_out":5147,"duration_ms":47021,"concrete_test":"Write a script that, for d=2 and at least 10^3 randomly sampled admissible intermediate states u^⋆∈D and flux differences ΔΦ, forms A and B from the actual quadratic-form coefficients (off-diagonal −a/2, diagonal α/2 and β), computes ρ(A^{-1/2}BA^{-1/2}) by dense eigendecomposition, and compares it with the paper's closed form max(|β0|/α0, |λ+|, |λ−|). Any mismatch—or any failure of the printed A to be positive definite—demonstrates the formula is incorrect or misprinted, and the Euler positivity claim must be revised.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The Euler bound-preserving argument in §6 reduces the blending coefficients to the spectral radius ρ(A^{-1/2}BA^{-1/2}) and then asserts a closed form, max(|β0|/α0, |λ+|, |λ−|), 'after some easy calculations'. This formula is never derived. More concerning, the matrices A and B as printed do not match the quadratic forms they are supposed to represent. The denominator shown immediately above is α0/2 ||w||^2 − a0^T w w_{d+1} + β0 w_{d+1}^2, whose symmetric matrix has off-diagonal block −a0/2, not −2a as printed; similarly for B. Unless a and b are implicitly rescaled—which the text does not state—the claimed formula computes the spectral radius of a different matrix. The positive definiteness of A is also asserted to follow from u^⋆∈D, but this holds for the correctly formed quadratic form, not necessarily for the printed A. Because the Euler positivity proof consists entirely of choosing θ from this formula, any error or misprint here invalidates the central claim that the scheme is provably bound preserving for the Euler equations. The numerical results in §7 do not test this: they show shock cases on triangular meshes but do not isolate the validity of the closed form.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript presents the PAMPA discretization for hyperbolic conservation laws on polygonal meshes, a globally continuous active-flux-type scheme whose degrees of freedom are cell averages and boundary point values of a virtual-element-type quadratic space. The average DoFs are updated via the conservative form with a high-order quadrature flux blended against a low-order local Lax-Friedrichs flux, and the point values via the non-conservative form with blended residuals; the blending coefficients are chosen so that each Euler-forward update is a convex combination of admissible states, which under SSP-RK3 time stepping yields an invariant-domain-preserving scheme. A bound-preserving proof is sketched for scalar problems and, for the Euler equations, via the GQL characterization of Wu and Shu together with a claimed closed-form upper bound on the spectral radius of a symmetric-definite pencil. Numerical results are shown for the Zalesak rotation test and the four-shock (KT) interaction, both on triangular meshes, and the paper states that it is a summary of the companion article [1].","tokens_in":11602,"tokens_out":39598,"duration_ms":323290,"significance":"If correct, the scheme would be a valuable object: a third-order, globally continuous method on general polygons with few boundary degrees of freedom and a provable invariant-domain property. The convex-blending formulation is elegant, and the reduction of the whole analysis to per-intermediate-state admissibility is a clean and usable framework; the scalar analysis follows a standard and sound template, and the reported KT result is visually cleaner than the earlier MOOD variant. However, the Euler half of the central claim is not established: Section 6's closed-form spectral-radius formula is not only unproved but demonstrably incorrect in the generic two-dimensional case, and the printed matrices do not match the quadratic forms they are supposed to represent. Since the Euler blending coefficients are chosen from that formula, the abstract's claim that the scheme is 'provably bound preserving ... for the Euler equations' is unsupported by the manuscript as it stands.","major_comments":[{"comment":"The formula rho(A^{-1/2}BA^{-1/2}) = max(|beta0|/alpha0, |lambda+|, |lambda-|) is not merely asserted without derivation; it is wrong for the generic two-dimensional case. The expression corresponds to the generalized eigenproblem of the printed pencil restricted to the two-dimensional subspace spanned by (a,0) and (0,1), which is exact only when b is collinear with a. For the Euler equations, the vector a in the denominator is determined by the momentum of the intermediate state, while b is determined by momentum-flux differences (Delta(rho v tensor v + pI)); these are not proportional in general. A direct calculation in d = 2 with alpha0 = 3, a = (1,2), alpha_{d+1} = 4, beta0 = 2, b = (1,0), beta_{d+1} = 3 gives, for the printed matrices, det(B - lambda A) = (2 - 3lambda)(20lambda^2 - 26lambda + 8), so rho(A^{-1/2}BA^{-1/2}) = 0.8, whereas the printed formula yields lambda+ = (13 + sqrt(137))/(-4), hence max(|beta0|/alpha0, |lambda+|, |lambda-|) = (13 + sqrt(137))/4 = 6.18, an overestimate by a factor of about 7.7. Because the Euler positivity proof consists entirely of choosing theta from this value, the bound-preserving property for the two-dimensional Euler equations does not follow from the manuscript as written. The authors should replace the closed form by a correct evaluation of rho(A^{-1/2}BA^{-1/2}) (for instance the iterative symmetric eigensolver that the text itself mentions) or derive a correct formula with its exact validity conditions stated.","section":"Section 6, closed-form spectral-radius formula"},{"comment":"The matrices A and B printed in Section 6 do not represent the quadratic forms displayed immediately above them. The denominator alpha0/2 ||w||^2 - a0^T w w_{d+1} + beta0 w_{d+1}^2 equals (1/2) z^T A~ z with A~ = [[alpha0 I_d, -a0],[-a0^T, 2 beta0]], whereas the printed A has off-diagonal block -2a and lower-right block 2 alpha_{d+1}; the analogous discrepancy holds for B, whose entries should be built from the numerator coefficients (alpha1, a1, beta1) of the flux difference, not from (beta0, b, beta_{d+1}). The symbols a, alpha_{d+1}, b, and beta_{d+1} are never defined in terms of the intermediate state u^{fl,*}, the residual Delta Phi, or the flux difference Delta f-hat. Consequently, the assertion that positive definiteness of A is equivalent to the displayed denominator being positive for all nonzero z cannot be checked, and the reduction to the Rayleigh quotient is broken as printed. The section must be rewritten with a consistent notation that ties every entry of A and B to the coefficients of u^{fl,*} psi(w)^T and (Delta f-hat^{fl})^T psi(w).","section":"Section 6, definition of A and B"},{"comment":"The scalar bound-preserving argument is standard in outline, but as printed it cannot be checked. In Eq. (18) the added term is (alpha_E/|E|_sigma)(u^n_sigma - u^n_sigma), which is identically zero, so the displayed identity leading to the convex combination is algebraically false; the intended difference, presumably between u^n_sigma and the relevant edge or cell state, must be stated. The low-order residual in Eq. (14) defines u^E_{i-1,i} = (u_sigma_i + u_sigma_{i-1} + u_E)/6 but u^E_{i,i+1} with weight 1/3, has a missing closing parenthesis, and writes u^P_{i-1,i} in the first line versus u^E_{i,i+1} in the second; the later derivation of u^{E,*}_sigma as a convex combination of thirds is not consistent with these definitions as printed. Since the claim that u^{E,*}_sigma lies in D, and hence the whole scalar bound-preserving result, depends on these displays, the derivation needs to be written out cleanly and self-contained rather than reconstructed from [16].","section":"Section 4, Eqs. (14), (15), (18)"},{"comment":"Several load-bearing definitional and stability statements are deferred rather than proven or even stated precisely. The high-order residual in Eq. (17) uses the projector pi^perp, which is never defined in this paper; the non-degeneracy of N_sigma K^+_{sigma,E} is asserted with a reference to [16]; boundary treatment is skipped with a pointer to [1]; and the formal third-order accuracy and linear stability are cited from [12,13], which concern related but not identical schemes. For a summary paper this would be acceptable for secondary points, but for the abstract's central claim the manuscript should state explicitly which components the bound-preserving proof relies on and which are imported from the companion papers. In particular, 'provably bound preserving ... for the Euler equations' cannot be verified from the present text because of the issues in Section 6 described above.","section":"Sections 3-5, deferred components of the central claim"}],"minor_comments":[{"comment":"The symbols |E|_sigma, |E_sigma|, |mathcal E|_sigma, and |C_sigma| are used for what appears to be the same area associated with the point value sigma; the notation should be unified.","section":"Section 4, Eqs. (18)-(22)"},{"comment":"The sentence 'the problem reduces to study the Rayleigh quotient max_{z in R^{d+1}} z^T B z / |z^T B z|' has identical numerator and denominator up to the absolute value and is therefore identically plus or minus one; it should read max_z |z^T B z| / (z^T A z) or the equivalent.","section":"Section 6"},{"comment":"The text refers to 'a function of V^2_k(P)', but only V_2(P) has been defined; this is presumably a typo.","section":"Section 2.1"},{"comment":"The projector pi^perp appears without any definition in this paper; the high-order residual is only understandable by consulting [16].","section":"Eq. (17)"},{"comment":"Despite the abstract and title, all numerical examples are on triangular meshes; the conclusion itself acknowledges this. The paper also reports no quantitative confirmation of the bound-preserving property, such as the min and max of the scalar solution for Zalesak or the minimum density and internal energy for the KT case, which would be the most direct numerical evidence for the central claim.","section":"Section 7"},{"comment":"The abstract contains a grammar error ('Several cases demonstrates'), and the conclusion's phrase 'rigorously bound preserving' is stronger than what Section 6 as printed supports; the wording should be aligned with the actual content.","section":"Abstract and Section 8"}],"recommendation":"major_revision","confidential_remarks":"The manuscript leans heavily on the authors' own in-preparation companion paper [1] and on prior work [16,20,25] for definitions, stability, and proof details; this is acceptable for a summary paper only if the deferred claims are correct, and the counterexample in Section 6 shows that the Euler claim is not. The numerical tests do not exercise a genuinely two-dimensional state with non-collinear vectors a and b, so the failure of the formula is not detected numerically. If the authors replace the closed form with a correct spectral-radius evaluation and fix the notation, the framework is salvageable; otherwise the Euler part of the claim should be withdrawn or restricted to the one-dimensional case. I would also note that reference [1] is 'in preparation', so the published record would depend on an unreviewed source for the central proofs."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: the paper extends the one-dimensional bound-preserving PAMPA scheme to polygonal meshes and adds an eigenvalue-based recipe for choosing the convex blending coefficients for Euler. That is a real step beyond the authors' earlier work, and the scalar part of the argument is coherent. But the Euler positivity proof, which is the load-bearing novelty, is not actually present in this text: the key spectral-radius formula in §6 is asserted 'after some easy calculations', the matrices A and B as printed do not match the quadratic forms they are supposed to represent, and the derivation is not supplied. Until that formula is verified and the typos fixed, the claim 'provably bound preserving for Euler' is unsupported.\n\nWhat is actually new and good: the convex-limiting framework is standard IDP machinery, but the formulation in terms of point-value and cell-average DoFs on polygons is new, and the reduction of the Euler case to a Rayleigh quotient is a genuinely useful idea if it checks out. The numerical results on the KT case show a visible improvement over their earlier MOOD-stabilized PAMPA, and the Zalesak test confirms the solution stays in bounds. So the paper is not empty.\n\nSoft spots, in order of severity. The Euler section is the main problem. The denominator and numerator quadratic forms are written with coefficients that are inconsistent with the printed A and B: the off-diagonal block should be scaled differently. Either the scalings are wrong or there are typos; either way, the claimed closed-form λ± cannot be checked from the text. Positive definiteness of A is asserted, not shown. Equation (18) also has an obvious typo (u^n_σ minus itself), which undermines the convex rewriting for point values even though the intent is clear. The paper is explicitly a summary of an unpublished companion, so large chunks of the analysis are deferred rather than given. The numerical section only shows triangular meshes despite the polygonal claim, and no convergence-order test is included here.\n\nThe citation pattern is heavy on the authors' own prior work, which is natural for a continuation, but the companion paper is doing a lot of the work. The stress-test concern about the spectral-radius formula is correct; I looked at the same spot and the mismatch is real.\n\nBottom line: this deserves a serious referee, but the Euler bound-preserving claim should not be accepted on the strength of this text. The scalar part and the framework are publishable if the companion delivers the missing derivation. I would hand it to a referee who knows IDP and active flux, and ask for a corrected, complete Euler argument before final acceptance.","headline":"Useful extension of bound-preserving PAMPA to polygons, but the Euler proof is missing the key spectral-radius derivation and the printed matrices don't match the forms they claim.","tokens_in":12095,"tokens_out":3015,"would_cite":false,"duration_ms":27818,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65M08","65M60","65M12"],"pacs":[],"model":"deepseek-v4-flash","headline":"The paper shows that the PAMPA discretization on polygonal meshes is bound preserving: scalar solutions stay in the initial range, and Euler solutions keep positive density and internal energy.","keywords":["hyperbolic conservation laws","bound preserving","invariant domain","polygonal meshes","active flux","PAMPA","Euler equations","convex limiting"],"falsifier":"Take any admissible intermediate state in $\\mathbb{R}^4$, form the Section 6 matrices $A$ and $B$, compute the largest eigenvalue of $A^{-1/2}BA^{-1/2}$ with a standard symmetric eigensolver, and compare it with $\\max(|\\beta_0|/\\alpha_0, |\\lambda_+|, |\\lambda_-|)$; one mismatch would disprove the Euler bound-preserving proof. A complementary test is to run the KT two-dimensional Riemann problem on a Voronoi polygonal mesh at the stated CFL and monitor $\\min \\rho$ and $\\min e$: any negative value would falsify the method's claimed invariant-domain preservation.","tokens_in":11120,"feed_emoji":"📐","tokens_out":10567,"duration_ms":98367,"temperature":0.7,"pith_summary":"The paper proposes a high-order discretization for hyperbolic conservation laws, named PAMPA, that works on general polygonal meshes and is provably bound preserving for scalar problems and the Euler equations. The central move is to update cell averages and edge point values with a convex blend of a high-order flux or residual and a low-order local Lax-Friedrichs one, so that each Euler step can be rewritten as a convex combination of previous values and intermediate states. If those intermediate states lie in the invariant domain and a time-step restriction holds, the updated solution automatically stays in the domain. A sympathetic reader would care because this gives high-order accuracy and shock robustness without an expensive stabilization loop, and it works on the flexible polygonal meshes used for local refinement.","feed_headline":"Polygonal-mesh scheme keeps high-order solutions in physical bounds","feed_subtitle":"The PAMPA method blends high- and low-order fluxes so density and internal energy stay positive on general polygon meshes.","key_machinery":"The load-bearing construction is the convex blending of high-order and low-order updates. For a point value $\\sigma$ in polygon $P$, the residual is $\\Phi^P_\\sigma = \\Phi^{LO}_\\sigma + \\theta^P_\\sigma \\Delta \\Phi^P_\\sigma$, and for a face $f$ the numerical flux is $\\hat{f} = \\hat{f}^{LO} + \\theta_f \\Delta \\hat{f}$; the $\\theta$'s are limited coefficients. Because the low-order terms are chosen with a local Lax-Friedrichs structure and a wave-speed bound $\\alpha$, the update $u^{n+1}$ can be rewritten as a convex combination of $u^n$ and intermediate states $u^\\star$; if $u^\\star \\in D$ and the CFL condition $\\Delta t \\sum |f| \\alpha_f / |E| \\le 1$ holds, then $u^{n+1} \\in D$ follows by convexity. For the Euler equations, the invariant domain is represented through geometric quasilinearization as $\\{ u : u^T \\psi(w) > 0 \\}$ with $\\psi(w) = (\\|w\\|^2/2, -w, 1)$, and the optimal $\\theta$ is obtained by minimizing a ratio $u^{\\star T}\\psi(w) / |\\Delta \\Phi^T \\psi(w)|$, which reduces to the generalized eigenvalue problem $\\rho(A^{-1/2}BA^{-1/2})$ for two symmetric matrices. The paper states a closed-form expression for this spectral radius in terms of $\\lambda_\\pm$, which is what turns the bound-preserving condition into a computable formula.","core_discovery":"The paper's central claim is that the PAMPA scheme, a globally continuous approximation that stores cell averages and point values at Gauss-Lobatto nodes on element boundaries, can be made invariant-domain preserving on polygonal meshes without any polynomial reconstruction. For scalar conservation laws, the preserved set is the range of the initial data; for the Euler equations, it is the set of states with positive density and positive internal energy. The proof rewrites the forward-Euler update of every degree of freedom as a convex combination of the old value and intermediate states built from low-order fluxes, and it shows that each intermediate state is in the invariant domain under a standard CFL condition. The blending coefficients are then chosen, through convex limiting, so that the convex combination remains in the domain; for the Euler equations, the coefficients come from a generalized eigenvalue problem derived from the geometric quasilinearization of the invariant domain. Numerical tests on the Zalesak and KT cases show the expected third-order accuracy on smooth problems and cleaner shock resolution than the earlier MOOD-based version.","pith_inferences":["A natural next experiment is to verify the Section 6 closed-form eigenvalue formula against a dense eigensolver for random intermediate states; this would isolate the one unproved algebraic step in the Euler proof.","The convex-blending-plus-geometric-quasilinearization construction should extend to other hyperbolic systems, such as shallow water or ideal MHD, whenever their invariant domain is described by finitely many quadratic inequalities; the main new task is the corresponding generalized eigenvalue problem.","Because the limiting uses only the stored degrees of freedom, the bound-preserving mechanism is likely to carry over to higher-order VEM-like extensions and to hanging-node refinements without changing the convexity argument."],"forward_implications":["A high-order, globally continuous scheme on polygonal meshes can be made invariant-domain preserving by convex blending alone, with no MOOD-type fallback for the tested shock cases.","The same analysis converts the bound-preserving question into checking intermediate states and computing local blending coefficients, so adding new physics mainly means re-solving a small eigenvalue problem.","Combined with SSP Runge-Kutta time stepping, the forward-Euler bound-preservation property carries over to the full time step.","The polygonal setting supports local refinement more naturally, since the approximation space is defined through virtual-element degrees of freedom on polygons rather than mapped quadrilaterals."],"supporting_citations":[{"why":"Full version of this summary; carries the boundary treatment and complete proofs of the bound-preserving analysis.","marker":"[1]"},{"why":"Previous PAMPA algorithm on polygons; supplies the virtual-element approximation space, sub-triangulation, and high-order residuals this paper modifies.","marker":"[16]"},{"why":"Geometric quasilinearization framework; gives the conic description of the Euler invariant domain used to derive the eigenvalue-based blending coefficients.","marker":"[19]"},{"why":"One-dimensional bound-preserving PAMPA; source of the convex-blending and coefficient-limiting strategy extended here to multidimensional polygonal meshes.","marker":"[25]"},{"why":"Virtual element method; defines the degrees of freedom and projection operators that the polygonal approximation space relies on.","marker":"[26]"},{"why":"Fast upper-bound estimate for the maximum wave speed in the Riemann problem; supplies the dissipation coefficients alpha used by the low-order fluxes.","marker":"[27]"},{"why":"One-dimensional PAMPA bound-preserving variant whose data-only stabilization is the multidimensional extension described in this paper.","marker":"[20]"}],"fun_headline_variants":["PAMPA: provably bound-preserving scheme for polygonal meshes","Polygon-mesh scheme keeps positivity without reconstruction","Bound preserving meets Active Flux: PAMPA on polygons","High-order convex-limiter scheme stays in invariant domain","Euler and scalar bounds hold on arbitrary polygon meshes"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"For the Euler equations, the proof leans on a stated-but-unproved closed-form formula for the largest eigenvalue of $A^{-1/2}BA^{-1/2}$ and on the assertion that $A$ is positive definite whenever the intermediate state is admissible; if either is wrong, the chosen blending coefficients may fail to keep density and internal energy positive.","fun_headline_variants_meta":{"raw":{"variants":["PAMPA: provably bound-preserving scheme for polygonal meshes","Polygon-mesh scheme keeps positivity without reconstruction","Bound preserving meets Active Flux: PAMPA on polygons","High-order convex-limiter scheme stays in invariant domain","Euler and scalar bounds hold on arbitrary polygon meshes"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000527,"raw_usage":{"total_tokens":2481,"prompt_tokens":820,"completion_tokens":1661,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":436,"completion_tokens_details":{"reasoning_tokens":1579}},"tokens_in":436,"tokens_out":1661,"duration_ms":13650,"temperature":1.0,"reasoning_tokens":1579,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T19:30:55.650455+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Take any admissible intermediate state in $\\mathbb{R}^4$, form the Section 6 matrices $A$ and $B$, compute the largest eigenvalue of $A^{-1/2}BA^{-1/2}$ with a standard symmetric eigensolver, and compare it with $\\max(|\\beta_0|/\\alpha_0, |\\lambda_+|, |\\lambda_-|)$; one mismatch would disprove the Euler bound-preserving proof. A complementary test is to run the KT two-dimensional Riemann problem on a Voronoi polygonal mesh at the stated CFL and monitor $\\min \\rho$ and $\\min e$: any negative value would falsify the method's claimed invariant-domain preservation.","supporting_citations":[{"cited_title":"Bound preserving virtual finite element and hyperbolic problems: the PAMPA algorithm","cited_arxiv_id":null,"evidence_quote":"Full version of this summary; carries the boundary treatment and complete proofs of the bound-preserving analysis."},{"cited_title":"Geometric quasilinearization framework for analysis and design of bound-preserving schemes","cited_arxiv_id":null,"evidence_quote":"Geometric quasilinearization framework; gives the conic description of the Euler invariant domain used to derive the eigenvalue-based blending coefficients."},{"cited_title":"Beir˜ ao da Veiga, Franco Brezzi, L","cited_arxiv_id":null,"evidence_quote":"Virtual element method; defines the degrees of freedom and projection operators that the polygonal approximation space relies on."},{"cited_title":"Fast estimation from above of the maximum wave speed in the Riemann problem for the Euler equations","cited_arxiv_id":null,"evidence_quote":"Fast upper-bound estimate for the maximum wave speed in the Riemann problem; supplies the dissipation coefficients alpha used by the low-order fluxes."},{"cited_title":"Bound preserving point-average-moment polynomial-interpreted (PAMPA) scheme: one-dimensional case, 2024","cited_arxiv_id":null,"evidence_quote":"One-dimensional PAMPA bound-preserving variant whose data-only stabilization is the multidimensional extension described in this paper."}],"review_version":1}