{"id":"f37e7325-24b0-43b6-814a-f0328a2875ef","arxiv_id":"2412.16307","paper_version":1,"verdict":"CONDITIONAL","confidence":"HIGH","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":5,"one_line_summary":"A numerical study of a sulphation PDE-ODE model with a Pearson-diffusion boundary condition, including a positivity-preserving discretization and simulations of noise propagation.","lead":"This paper builds a numerical scheme for a stochastic model of marble sulphation, where the sulphur dioxide level at the stone surface follows a random Pearson diffusion. The authors show the scheme preserves positivity and stability, and they explore how boundary noise spreads into the material in slow and fast reaction regimes.","discovery_kind":"new_application","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Proposition 4.3's positivity proof silently assumes the discrete boundary sequence ~ψ^n lies in [0,~η); the algorithm (41) never specifies how the continuous-time integral is discretized, so the theorem is not yet connected to the implemented scheme.","rationale":"The reader's weakest assumption identifies a real presentation gap: Proposition 4.3 requires bounded discrete boundary values, but the algorithm for computing those values via (41) is underspecified. I do not think this gap is fatal, because the desired bound ~ψ^n ∈ [0,~η] follows from the elementary observation that ψ^n = η sin^2(y_n/2) ∈ [0,η] and that any nonnegative quadrature of a nonnegative integrand keeps the exponential factor in (0,1], so the denominator in (41) cannot drop below φ(c0). Thus the concern is about completeness and reproducibility rather than correctness of the underlying mathematics. The second concern raised by the reader, Theorem 3.2's assumption ~Ψ(0)=0, is not load-bearing for the numerical experiments because the simulations appear to take Ψ_0=0. Proposition 4.3 itself is a valid discrete maximum principle once the boundary bound is supplied. For these reasons, the paper's central claim is plausible and the identified issues are addressable, so the existing CONDITIONAL verdict remains appropriate; my stress-test does not move the verdict. I would nonetheless require the authors to add the missing boundary-bound lemma and to state explicitly the quadrature rule used for the integral in (41), as a condition of acceptance.","tokens_in":26524,"tokens_out":21957,"duration_ms":189018,"concrete_test":"Add a short lemma: for the LSST boundary values ψ^n = η sin^2(y_n/2) and any nonnegative quadrature Q_n approximating ∫_0^{t_n} Ψ_u du, one has c0 e^{-λQ_n} ≤ c0, hence φ(c0 e^{-λQ_n}) ≥ φ(c0), and therefore ~ψ^n ≤ η/φ(c0) = ~η. Then instrument the code to report max_n ~ψ^n and max_{m,n} s_m^n over the N=500 samples of Section 6, and verify numerically that these maxima stay below ~η and \\bar c0 respectively for both the slow (λ=1) and fast (λ=100) regimes.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central claim, Proposition 4.3, is a discrete maximum principle for the scheme (54)-(55). Its proof requires, for every n, the boundary value s_0^n = ~ψ^n to belong to [0,~η). The continuous boundary function is bounded in (36) by Proposition 3.1, but the actual simulation uses (41), where ~ψ^n is defined through the continuous-time integral exp(-λ∫_0^{t_n} Ψ_u du) of the SDE path. The paper never specifies the quadrature rule for this integral, nor does it prove that the computed discrete boundary sequence satisfies the needed bound. This is a genuine gap between the theorem's hypotheses and the algorithm as implemented. The gap is likely easy to close: since ψ^n = η sin^2(y_n/2) always lies in [0,η], and any nonnegative quadrature of a nonnegative integrand gives an exponential factor in (0,1], the denominator φ(c0 exp(-λQ)) is at least φ(c0), so ~ψ^n ≤ η/φ(c0) = ~η. However, this argument is absent, so the theorem as stated does not formally cover the numerical method actually run. A second, smaller gap is Theorem 3.2's unflagged assumption ~Ψ(0)=0, which is inconsistent with (35) for general Ψ_0>0; the numerical experiments appear to use Ψ_0=0, so this is a presentation issue rather than a fatal one.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"This paper presents a numerical study of a half-line reaction-diffusion/ODE system modelling marble sulphation, with a stochastic Dirichlet boundary condition given by a Pearson diffusion. The authors recall well-posedness and boundedness of the boundary SDE, simulate it through a Lamperti transform combined with a sloping smooth truncation (LSST), and then discretize the PDE by a splitting strategy: an FTCS heat equation carrying the stochastic boundary data, coupled with an explicit finite-difference scheme for the nonlinear reaction part with deterministic boundary condition. The main theoretical contributions are the pathwise stability of the heat component (Proposition 4.2) and the positivity/boundedness result for the full discrete system (Proposition 4.3) under explicit conditions (47), (56), and (57). The paper closes with extensive numerical experiments on noise propagation, porosity effects, slow and fast reaction regimes, moving fronts, and statistical summaries over 500 sample paths.","tokens_in":26876,"tokens_out":9237,"duration_ms":85256,"significance":"If the stated results are fully connected to the implemented algorithm, the paper provides a useful, positivity-preserving and stable numerical method for a stochastic boundary-value problem that is relatively unexplored numerically. The strengths are the explicit and checkable coefficient conditions in Proposition 4.3, the pathwise character of the stability analysis, the use of a Lamperti-type transform to preserve the boundary domain, and the careful qualitative experiments that identify the moving-front regime. The honest deferral of full convergence analysis to the companion paper [5] is a clear limitation but is not by itself disqualifying for a numerical study. However, two load-bearing gaps currently separate the theorems from the code: the inverse Lamperti transform in the displayed schemes misses the factor eta, and the boundary sequence used in the PDE scheme is defined through an unspecified continuous-time integral with no proof that its discretization satisfies the hypothesis of Proposition 4.3.","major_comments":[{"comment":"The inverse Lamperti transform is correctly stated in Eq. (17) as Psi = eta sin^2(Y/2), but the LSST sampling scheme (25) and the time-continuous version (26) define psi^n = sin^2(y^n/2) and psi_t = sin^2(y_t/2) without the factor eta. Since all later bounds and experiments use eta = 1.5, the boundary process generated by the scheme as written lies in [0,1) rather than in [0,eta), which would either scale the boundary data by 1/eta relative to the model or indicate a typo in the displayed equations. This factor must be restored, because the boundary data enter every PDE simulation and the bound (36) depends on eta.","section":"Section 2.3, Eqs. (25)-(26)"},{"comment":"The discrete boundary sequence feeding the scheme is defined in (41) through the continuous-time integral exp(-lambda integral_0^{t_n} Psi_u du), but no quadrature rule for that integral is specified, and no proof is given that the computed values tilde-psi^n remain in [0,tilde-eta). Proposition 4.3's proof assumes exactly that interval property for every n, and its statement covers only m = 1,...,M, although the boundary value m = 0 is needed in the heat scheme (42). The gap is easy to close: with a nonnegative quadrature of a nonnegative integrand, exp(-lambda Q) in (0,1], so phi(c0 exp(-lambda Q)) >= phi(c0) and tilde-psi^n <= eta/phi(c0) = tilde-eta. This argument should be stated and the quadrature chosen explicitly so that the theorem formally covers the implemented algorithm.","section":"Section 4.1, Eq. (41) and Proposition 4.3"},{"comment":"The statement of Theorem 3.2 assumes tilde-Psi(0)=0 as part of the hypotheses, but this is not a consequence of condition (10) and is inconsistent with (35) when Psi_0>0. In that case the splitting u(0,x)=0, v(t,0)=0 creates a boundary discontinuity unless s_0=0. The theorem should either state the compatibility condition explicitly and restrict the admissible data, or the splitting should be modified to absorb a nonzero initial boundary value. The numerical experiments appear to use Psi_0=0, so this may be a presentation issue in practice, but the theorem as written is formally incomplete for general data.","section":"Theorem 3.2"}],"minor_comments":[{"comment":"The summation in the variation-of-constants formula appears off by one: for the recurrence U^{n+1}=A U^n + bar-Delta tilde-U_0^n, the solution should sum j=0,...,n-1 with A^{n-1-j}, not j=0,...,n with A^{n-j}. As written, U^n depends on tilde-U^n, which is not yet available at time t_n.","section":"Section 4.1, Eq. (45)"},{"comment":"Some captions and inline references use inconsistent parameter names: Figure 2 refers to sigma2 = 1 where Table 3 and the surrounding text use sigma3 = 1, and Figure 9(a) is labelled lambda1 = 0.1 while Table 3 sets lambda1 = 1. These should be harmonized.","section":"Figures 2 and 9"},{"comment":"The text says the accelerated regime considers lambda_i, i=2,3, but the discussion also refers to lambda1; please clarify whether the slow case is lambda1=1 or the lambda1=0.1 shown in Figure 9(a).","section":"Section 5.2"},{"comment":"The sentence introducing the numerical spatial accuracy estimate should state more explicitly that only the spatial error is measured for a fixed time step, and that the full convergence analysis is deferred to [5]; this would avoid the impression that the scheme's complete error is being estimated.","section":"Section 4.2.2"}],"recommendation":"major_revision","confidential_remarks":"The paper is a reasonable numerical study with fixable gaps. The missing factor eta in the displayed LSST scheme and the missing quadrature/boundary-invariance argument are load-bearing for the claim that the implemented scheme satisfies Proposition 4.3, but both can be repaired without changing the overall approach. I do not see grounds for rejection."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this is a solid, honest applied paper. It assembles a Lamperti–LSST boundary sampler and an FTCS splitting scheme for a stochastic sulphation model, proves explicit positivity/boundedness/stability conditions for the PDE scheme, and documents noise propagation in slow and fast regimes. The numerical exploration is genuinely useful for the cultural-heritage modeling community.\n\nWhat is new: the model and splitting come from earlier work by the same group ([35] and prior), and the LSST sampler is from Chen & Gan [15]. The genuine contributions are Proposition 4.3's discrete invariant-rectangle conditions (56)-(57), the assembly of the scheme, and the experiments showing how boundary noise spreads through the domain, with percentiles, standard deviations, and RMSD against the deterministic case. The proof of Proposition 4.3 is coherent and the conditions are checkable; the experiments support the qualitative claims. Credit where due: the paper is transparent about deferring convergence analysis to the companion paper [5], and it does not oversell.\n\nThe soft spots are real but not fatal. First, Proposition 4.3 assumes the discrete boundary sequence ~ψ^n stays in [0,~η). In the actual algorithm, ~ψ^n is computed via (41), which involves a continuous-time integral of the SDE path; no quadrature is specified, and the theorem does not cover the discretized version of that integral. This is a genuine gap between theorem and implementation. It is probably easy to close—any nonnegative quadrature of the nonnegative integrand keeps the exponential factor in (0,1], so ~ψ^n stays below ~η—but that argument is absent. Second, Theorem 3.2 silently assumes ~Ψ(0)=0, which contradicts (35) when Ψ_0>0; the experiments start from Ψ_0=0, so this is a presentation bug, not a numerical one. Third, no code or data are provided, and the spatial accuracy is only estimated, not proven (that is explicitly deferred). These are all addressable in revision.\n\nWho is it for: applied mathematicians working on stochastic boundary conditions, positivity-preserving finite differences, or cultural-heritage degradation models. It deserves a serious referee, not a desk reject. The referee should ask for the quadrature to be specified and analyzed, the Theorem 3.2 assumption to be fixed or flagged, and ideally for code or data to support the experiments.","headline":"Useful applied numerical study with a checkable stability result, but the main theorem's hypotheses are not yet connected to the implemented scheme; worth refereeing.","tokens_in":27382,"tokens_out":2519,"would_cite":false,"duration_ms":22742,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65M06","65C30","35K57","60H15","65M12"],"pacs":[],"model":"deepseek-v4-flash","headline":"The paper constructs a fully discrete, positivity-preserving scheme for the stochastic sulphation model and proves pathwise boundedness under explicit mesh and initial-data conditions.","keywords":["sulphation","stochastic dynamical boundary condition","Pearson diffusion","Lamperti transform","FTCS finite difference scheme","positivity preservation","moving front","cultural heritage degradation"],"falsifier":"Run the implemented scheme with boundary data built from a concrete quadrature for (41), starting with $\\Psi_0$ close to $\\eta$ and a coarse $\\Delta t$; if any computed $\\tilde{\\psi}^n$ leaves $[0,\\tilde{\\eta})$ or any $s$-value turns negative, Proposition 4.3's conclusion fails for that quadrature and the proof's hypothesis is violated.","tokens_in":26319,"feed_emoji":"🪨","tokens_out":11174,"duration_ms":94417,"temperature":0.7,"pith_summary":"The paper sets out to show that the randomly forced sulphation system—sulphur dioxide eating into calcium carbonate with a stochastic boundary concentration—can be simulated by a fully discrete scheme that respects the physics of the model. The construction splits the solution into a heat component that carries the random boundary and a nonlinear component with a deterministic boundary, then samples the boundary with a Lamperti-transformed Pearson diffusion and advances both parts with explicit finite differences. The authors prove pathwise positivity, boundedness, and stability under explicit conditions on the spatial grid, the time step, and the initial calcite density. If these claims hold, the scheme is a working laboratory for qualitative questions: how boundary noise spreads into the material, how porosity changes penetration, and how fast reactions sharpen the degradation front. The wider value is a template for random-boundary reaction-diffusion problems where concentrations must never leave their physical range.","feed_headline":"A positivity-preserving scheme tracks noise in stone sulphation","feed_subtitle":"A split isolates the random boundary, and explicit stability conditions keep simulated concentrations physical.","key_machinery":"The load-bearing mechanism is the splitting $s=u+v$, in which $u$ solves the heat equation with the stochastic boundary condition and $v$ solves a nonlinear, nonlocal equation with zero boundary data; diffusion smooths the irregular boundary noise before it reaches the reaction part. Boundary sampling is carried by a Lamperti transform $Y=2\\arcsin(\\sqrt{\\Psi/\\eta})$ that converts the non-Lipschitz Pearson diffusion into additive noise, followed by a sloping smooth truncation of the drift so that the explicit Euler\\,--\\,Maruyama update is stable and monotone. The PDE discretization is a forward-time, centred-space (FTCS) scheme, and the proof of Proposition 4.3 runs on coefficient positivity: under (47), (56), and (57), every coefficient multiplying the old values in the $s$-update is nonnegative, so the interval $[0,\\tilde{\\eta})\\times[0,\\bar c_0]$ is invariant pathwise.","core_discovery":"On the paper's own terms, the central discovery is that the splitting $s=u+v$ used for the well-posedness of the continuous problem can be transplanted to a numerical scheme with a rigorous positivity statement. The discrete heat component $u$ inherits the stochastic boundary values $\\tilde{\\psi}^n$, produced by the Lamperti\\,--\\,sloping-smooth-truncation sampler, while the nonlinear component $v$ and the calcite density $c$ are advanced by an explicit scheme. Proposition 4.3 proves by induction that every pathwise solution of the $s$-update (54) and $c$-update (55) stays in $[0,\\tilde{\\eta})\\times[0,\\bar c_0]$ whenever $\\bar\\Delta\\le 1/2$, the initial data satisfy (56), and the time step satisfies (57); the induction rests on checking that each coefficient of the update is nonnegative under these conditions. Corollary 4.4 then bounds $v$ between $-\\tilde{\\eta}$ and $\\tilde{\\eta}$, and Proposition 4.5 converts those bounds into $L^2$ and maximum-norm stability, pathwise and in mean. The numerical sections use the scheme to describe slow and fast reaction regimes, including the formation of a moving front with variance concentrated around it.","pith_inferences":["Editorial extension: no quadrature rule for the integral in (41) is specified, so the implemented boundary sequence could leave the interval $[0,\\tilde{\\eta})$; proving a discrete analogue of Proposition 3.1 for that quadrature is the missing step that would make Proposition 4.3 airtight.","Editorial extension: the same Lamperti\\,--\\,LSST-plus-FTCS construction should transfer to other bounded Pearson boundary noises, such as Wright\\,--\\,Fisher processes with $\\eta=1$, where positivity of concentrations or frequencies is a hard modelling constraint.","Editorial extension: the concentration of sample variance near the moving front suggests a stochastic front-speed law; one could test whether the front position scales like $\\sqrt{\\lambda t}$ with noise-induced corrections by estimating it from many realizations.","Editorial extension: if the scheme is correct, it gives a direct way to quantify uncertainty in cultural heritage risk assessments, since it returns full space-time distributions of pollutant concentration and material loss for given environmental noise parameters."],"forward_implications":["If the stability bounds hold, the scheme offers a practical recipe for choosing time steps when simulating the stochastic sulphation model: condition (57) ties $\\Delta t$ to the spatial mesh, the reaction rate $\\lambda$, and the initial calcite bound.","The simulations imply that boundary noise does not stay at the boundary: for slow reactions the variance of both $\\rho$ and $c$ spreads gradually over the whole domain, while for fast reactions it concentrates in a thin region around the moving front.","In the fast-reaction regime the discrete solutions display the same qualitative front formation known from the deterministic fast-reaction limit, now with stochastic boundary data; this supports using the scheme to study black-crust formation under random pollution histories.","The numerically estimated spatial accuracy order is about one for both sulphur dioxide and calcite, which is consistent with the explicit first-order finite-difference structure and gives an error scale for the reported statistics."],"supporting_citations":[{"why":"supplies the splitting $s=u+v$ and the well-posedness theory of the continuous system that the discrete scheme is built to mimic.","marker":"[35]"},{"why":"provides the sloping smooth truncation construction and the strong convergence rates used for the boundary sampler in Proposition 2.8.","marker":"[15]"},{"why":"defines the deterministic sulphation model and its numerical and asymptotic setting that the stochastic system extends.","marker":"[6]"},{"why":"furnishes the deterministic boundary-regularity setting and the fast-reaction scaling used in Sections 3 and 5.","marker":"[30]"},{"why":"introduces the Pearson diffusion as the stochastic boundary condition and the parameter values drawn from pollutant time series.","marker":"[3]"},{"why":"presents the macroscale stochastic model and the splitting strategy on which the discrete construction is based.","marker":"[4]"},{"why":"gives the inverse-moment bounds for Pearson diffusions used to prove bounded inverse moments in Proposition 2.4.","marker":"[44]"},{"why":"supplies the taming and explicit-scheme framework for non-globally Lipschitz SDEs that motivates the sloping smooth truncation of the drift.","marker":"[31]"}],"fun_headline_variants":["Stochastic boundary noise tamed in stone sulfation model","Positivity-preserving solver for noisy stone decay","Splitting scheme yields stable stochastic sulfation simulation","Random boundary dynamics captured in stone decay numerics"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the discrete boundary values fed into the scheme stay inside the allowed interval as the proof assumes, and the paper does not verify this for the quadrature it uses to compute them; a second fragile point is that the well-posedness theorem silently assumes the boundary starts at zero, which is not true for general positive initial data.","fun_headline_variants_meta":{"raw":{"variants":["Stochastic boundary noise tamed in stone sulfation model","Positivity-preserving solver for noisy stone decay","Splitting scheme yields stable stochastic sulfation simulation","Random boundary dynamics captured in stone decay numerics"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000194,"raw_usage":{"total_tokens":1350,"prompt_tokens":939,"completion_tokens":411,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":555,"completion_tokens_details":{"reasoning_tokens":349}},"tokens_in":555,"tokens_out":411,"duration_ms":4202,"temperature":1.0,"reasoning_tokens":349,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-11T10:43:27.330012+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Run the implemented scheme with boundary data built from a concrete quadrature for (41), starting with $\\Psi_0$ close to $\\eta$ and a coarse $\\Delta t$; if any computed $\\tilde{\\psi}^n$ leaves $[0,\\tilde{\\eta})$ or any $s$-value turns negative, Proposition 4.3's conclusion fails for that quadrature and the proof's hypothesis is violated.","supporting_citations":[{"cited_title":"Maurelli, D.Morale, S.Ugolini, Well-posedness of a reaction-diffusion model with stochas- tic dynamical boundary conditions arXiv preprint, arXiv:2 308.06847 (2023)","cited_arxiv_id":null,"evidence_quote":"supplies the splitting $s=u+v$ and the well-posedness theory of the continuous system that the discrete scheme is built to mimic."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"provides the sloping smooth truncation construction and the strong convergence rates used for the boundary sampler in Proposition 2.8."},{"cited_title":"Aregbda-Driollet, F","cited_arxiv_id":null,"evidence_quote":"defines the deterministic sulphation model and its numerical and asymptotic setting that the stochastic system extends."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"furnishes the deterministic boundary-regularity setting and the fast-reaction scaling used in Sections 3 and 5."},{"cited_title":"Arceci, L.M","cited_arxiv_id":null,"evidence_quote":"introduces the Pearson diffusion as the stochastic boundary condition and the parameter values drawn from pollutant time series."},{"cited_title":"A react ion diffusion model with a stochas- tic boundary condition","cited_arxiv_id":null,"evidence_quote":"presents the macroscale stochastic model and the splitting strategy on which the discrete construction is based."},{"cited_title":"Sorensen, J.L","cited_arxiv_id":null,"evidence_quote":"gives the inverse-moment bounds for Pearson diffusions used to prove bounded inverse moments in Proposition 2.4."},{"cited_title":"Hutzenthaler, A","cited_arxiv_id":null,"evidence_quote":"supplies the taming and explicit-scheme framework for non-globally Lipschitz SDEs that motivates the sloping smooth truncation of the drift."}],"review_version":1}