{"id":"bad64c8d-c5ee-49bf-93eb-b8f5134e6c18","arxiv_id":"1908.02032","paper_version":4,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"low","formal_verification":"none","parameter_count":2,"one_line_summary":"Explicit Zolotarev-based pole choices for rational Krylov evaluation of Stieltjes matrix functions yield proven convergence rates, including Kronecker-structured cases.","lead":"This paper gives explicit rules for choosing the poles of rational Krylov methods when computing functions like exponentials or negative powers of large positive definite matrices. The rules come with proven convergence rates, and extend to Kronecker-structured problems that arise in fractional diffusion and exponential integrators.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Laplace-Stieltjes bound inherits an unproven saddle-point constant γℓ,κ; if the conjectured γ=1 fails the rate in Cor. 3.16/4.1 still holds but the constant is unvalidated.","rationale":"The reader's weakest_assumption (known spectral bounds) is the most obvious practical limitation, but it is explicitly acknowledged in Section 5 and is shared by essentially all Zolotarev-based analyses. The internally unverified constant in Lemma 3.5 is more load-bearing in the formal sense: it is the one step in the Laplace-Stieltjes chain that is not checked, and the authors themselves note in Remark 3.7 that the bound is not sharp and that γ=1 is plausible. However, the uncertainty is about the size of a constant only; the exponential rates in ρ^{ℓ/2} and ρ^ℓ follow from the rational approximation results (Theorem 3.10 and the Möbius transformations) independently of the constant. A wrong γ would not change the headline convergence-rate claims, and a concrete check of Lemma 3.5 would settle the residual doubt. The Cauchy-Stieltjes bounds are supported by published Zolotarev results and the numerical experiments; no additional formal verification need is apparent beyond the routine constant bookkeeping in Corollaries 3.14 and 4.8. Therefore the conditional acceptance stands, and the only refinement I would suggest is verifying Appendix B before relying on the specific constants in future work. I partially agree with the reader because I see the spectral-interval assumption as less dangerous (it is stated as a hypothesis, and the numerical tests use exact bounds) while the unverified constant is the actual soft spot in the proof chain.","tokens_in":29438,"tokens_out":2109,"duration_ms":20156,"concrete_test":"Independent verification of Lemma 3.5, e.g., enumerate random sets of ℓ poles in [a,b] for several b/a and compute L⁻¹[(1/s)p(s)/p(−s)] numerically by quadrature over s=it; check whether the maximum over t is bounded by the claimed γℓ,κ and whether it exceeds 1. Equivalently, re-derive the bound in Appendix B symbolically and verify each inequality, or implement the proof in a proof assistant to remove the sole unverified step. A second test: repeat the experiment of Figure 3.1 with a deliberately inaccurate spectral bound (e.g., a/2 and 2b) and check whether the observed rate degrades in a way consistent with the theorem applied to the wrongly assumed interval.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central quantitative claim for Laplace-Stieltjes functions (Corollaries 3.16 and 4.1, Theorem 4.10) rests on Theorem 3.6, whose key constant γℓ,κ comes from Lemma 3.5: ‖L⁻¹[(1/s)(p(s)/p(−s))]‖_∞ ≤ 2.23 + (2/π)log(4ℓ√κ/π). The proof in Appendix B is long and elementary, but the reader did not machine-check it, and the paper itself flags in Remark 3.7 that experiments suggest the true sharp constant is γ=1 for all admissible poles, not growing with ℓ. Remark 3.7 is a self-acknowledged gap: the bound is currently a logarithmic function of ℓ that may be loose, and the proof's estimate is not known to be tight. However, the rate ρ^{ℓ/2} is unaffected by this constant; a numerically observed constant larger than the bound would be a contradiction, while a smaller constant only strengthens the bound. The more central issue is therefore not the constant itself but the dependence on the assumption that the spectrum is exactly in [a,b]: the poles Ψ_[a,b] are constructed from a and b, and all bounds are formulated in terms of the same interval used to define the poles. If the spectral interval is overestimated in practice, the rate still applies with the larger interval; if the spectrum is underestimated, the proven theorem does not apply, and the numerical experiments with tridiag(−1,2,−1) use the known exact spectral bounds, so their agreement with the bound does not test the sensitivity to inaccurate bounds. This is a limitation, but it was openly acknowledged in Section 5 and is standard for Zolotarev-based analyses; it does not undermine the mathematical claims. The argument for Cauchy-Stieltjes cases uses the Möbius maps and the inequalities ρ_[â,1] ≤ ρ_[a,4b] (Lemma 3.12) and ρ_[ã,1] ≤ ρ_[a,2b] (Lemma 4.7), which are asymptotically tight for large b/a; for moderate b/a the 4b and 2b factors introduce a small overestimate, not a risk.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript develops pole selection strategies for rational Krylov subspace methods that approximate x = f(M)v when f is a Laplace-Stieltjes or Cauchy-Stieltjes function and M is positive definite. For a Hermitian positive definite A with spectrum contained in [a,b], the authors prove explicit convergence bounds: Laplace-Stieltjes functions converge like O(rho^{ell/2}_{[a,b]}) using Zolotarev poles (Corollary 3.16) and Cauchy-Stieltjes functions converge like O(rho^ell_{[a,4b]}) using a Mobius-transformed Zolotarev set (Corollary 3.18). They also construct nested pole sequences via equidistributed sequences (Section 3.5). For Kronecker-structured matrices M = I otimes A - B^T otimes I with low-rank right-hand sides, they analyze tensorized rational Krylov subspaces and prove analogous rates: O(rho^{ell/2}_{[a,b]}) for Laplace-Stieltjes functions (Corollary 4.1) and O(rho^ell_{[a,2b]}) for Cauchy-Stieltjes functions (Corollary 4.8). Finally, they derive singular-value decay bounds for the solution matrix X, extending known low-rank approximability results for Sylvester equations (Theorems 4.10 and 4.11). Numerical experiments on finite-difference matrices, including fractional diffusion and exponential integrator examples, support the predicted rates.","tokens_in":29819,"tokens_out":14301,"duration_ms":139990,"significance":"If the results are correct, this is a valuable contribution to the rational Krylov literature. The paper provides explicit, non-fitted convergence rates and concrete pole sets for two important classes of matrix functions, and it extends the analysis to Kronecker-structured problems with low-rank right-hand sides, yielding both practical tensorized rational Krylov algorithms and new low-rank approximability bounds. The proofs are built on standard machinery -- Zolotarev theory, skeleton approximation, inverse Laplace transform estimates, and ADI error analysis -- and the constants are explicit rather than fitted to numerical data. The authors are appropriately candid about the limits of the analysis: Remark 3.7 acknowledges that the constant gamma_{ell,kappa} in Lemma 3.5 is likely not sharp, and Section 5 states that the non-normal case is left open. These caveats do not undermine the main rate claims, because gamma_{ell,kappa} only affects constants and not the exponential factor rho^{ell/2}.","major_comments":[],"minor_comments":[{"comment":"The constant gamma_{ell,kappa} is presented as an upper bound, but Remark 3.7 reports numerical evidence suggesting that the true constant is 1 for all admissible poles. Since Corollaries 3.16, 4.1, and Theorem 4.10 all inherit this constant, the authors should state more explicitly that (i) the exponential rate in Table 1.1 is unaffected by gamma, and (ii) the logarithmic factor in gamma is subexponential and therefore does not change the O(rho^{ell/2}) classification when ell is treated as the growing parameter.","section":"Section 3.2, Lemma 3.5 and Remark 3.7"},{"comment":"The phrase 'By setting I_A = I, I_B = -I in the statement of Theorem 4.6' is confusing: I_A and I_B are intervals, not identity matrices. Please rewrite this sentence using the normalized intervals, e.g. I_A = [\\tilde{a},1] and I_B = [-1,-\\tilde{a}], so that the Mobius transformation argument is unambiguous.","section":"Section 4.3, proof of Corollary 4.8"},{"comment":"The numerical experiments on tridiag(-1,2,-1) use the exact spectral interval, so they do not probe how the method behaves when the available bounds on [a,b] are over- or underestimated. Since all theorems are stated in terms of the interval used to construct the poles, adding a sentence in Section 5 on practical estimation of [a,b] (e.g., by Lanczos or Ritz values) would help the reader apply the bounds in realistic settings.","section":"Section 3.6 and Section 5"},{"comment":"There are a few minor typos: the acknowledgment reads 'The author wish' instead of 'The authors wish', and the introduction contains a duplicated 'of of' in the phrase 'the solution of of Sylvester equations' (Section 1.1, item 5). These should be corrected in a final pass.","section":"Throughout"}],"recommendation":"minor_revision","confidential_remarks":"The paper is within the scope of the journal and the central claims appear sound. The main caveats -- the non-sharp constant in Lemma 3.5 and the restriction to Hermitian positive definite matrices -- are explicitly acknowledged by the authors and do not affect the exponential convergence rates. I recommend minor revision rather than rejection."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this is a well-built numerical analysis paper. The main new content is real. For Hermitian positive definite A, rational Krylov with Zolotarev-based poles gives explicit O(rho^{ell/2}) convergence for Laplace-Stieltjes f(A)v and O(rho^ell) for Cauchy-Stieltjes, with the interval-dependent rate. The Kronecker-structured extension for M = I tensor A - B^T tensor I, including low-rank right-hand sides and singular value decay bounds, is also new relative to the cited Benzi-Simoncini and Beckermann-Townsend work. The nested EDS pole sequences are a practical bonus. I checked the logical chain from rational approximation to the corollaries; the proofs use standard machinery (Zolotarev, skeleton approximation, inverse Laplace transform estimates) and the constants are explicit, not fitted to the experiments. That is a real point in the paper's favor.\n\nWhere I would push back: the numerical section does not say how the reference solutions were computed, and no code or data is shipped. That is fixable, but it means the experiments currently read as illustrative rather than independently reproducible. Also the spectral-interval assumption [a,b] is load-bearing: the poles are built from a and b, and if the user underestimates the interval the theory does not apply. The paper says this in Section 5, but the experiments always use exact known bounds, so the sensitivity is untested. Minor. The gamma constant in Lemma 3.5 and Remark 3.7 is acknowledged as likely non-optimal; the stress-test note is right that if gamma is smaller the bound only gets better, and the rate is unaffected. I did not machine-check the long Appendix B estimate, but the structure is elementary and the claimed constant is plausible; I would trust it with one reviewer pass.\n\nOverall: this deserves serious refereeing. The results are useful to people computing f(A)v for completely monotonic functions, especially in the Kronecker-structured and low-rank settings. I would recommend acceptance after the experimental reproducibility issues are addressed and a sentence is added about the unproven sharpness of gamma.","headline":"Solid rational Krylov pole-selection paper with explicit bounds for Laplace- and Cauchy-Stieltjes functions; the math is credible, and the main caveats are experimental documentation and an unproven sharp constant that does not affect the rates.","tokens_in":30483,"tokens_out":1782,"would_cite":true,"duration_ms":19263,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65F60","65F30","41A20"],"pacs":[],"model":"deepseek-v4-flash","headline":"Rational Krylov poles chosen from the Zolotarev problem make the action of Laplace-Stieltjes and Cauchy-Stieltjes matrix functions converge exponentially, with rates expressed purely in terms of the spectral interval and with analogous…","keywords":["rational Krylov methods","Stieltjes matrix functions","Laplace-Stieltjes functions","Cauchy-Stieltjes functions","Zolotarev poles","pole selection","Kronecker structure","low-rank approximation"],"falsifier":"Pick the 1D Laplacian $A$ of size $10^5$ with known spectrum in $[a,b]\\approx[1/n,4]$, set $f(z)=e^{-z}$, and run rational Krylov with the paper's Zolotarev poles for $\\ell=5,10,15,20$. The log-error should fall with slope at least $-\\pi^2/(2\\log(4b/a))$; a markedly flatter slope, or failure of the error to decay once past the constant regime, would falsify the bound.","tokens_in":29242,"feed_emoji":"📉","tokens_out":9264,"duration_ms":80311,"temperature":0.7,"pith_summary":"This paper proves that rational Krylov methods with carefully chosen Zolotarev poles evaluate Stieltjes matrix functions on vectors with exponential convergence, where the rate is determined only by the matrix spectrum. For Laplace-Stieltjes functions, such as matrix exponentials and completely monotonic functions, the error after $\\ell$ steps is $O(\\rho_{[a,b]}^{\\ell/2})$; for Cauchy-Stieltjes functions, such as inverse fractional powers, it is $O(\\rho_{[a,4b]}^{\\ell})$, with $\\rho_{[\\alpha,\\beta]}=\\exp(-\\pi^2/\\log(4\\beta/\\alpha))$. The same analysis extends to Kronecker-sum matrices $M=I\\otimes A-B^T\\otimes I$ acting on vectorized low-rank matrices, where tensorized rational Krylov subspaces preserve structure and the error rates become $O(\\rho_{[a,b]}^{\\ell/2})$ and $O(\\rho_{[a,2b]}^{\\ell})$. This matters because fractional diffusion solvers, exponential integrators, and Sylvester-type matrix equations are precisely the settings where these functions and structures arise.","feed_headline":"Stieltjes matrix functions converge exponentially with the right poles","feed_subtitle":"The per-step error rate depends only on the eigenvalue interval, e.g. exp(−π²/log(4b/a)) for Laplace-Stieltjes functions.","key_machinery":"The engine of the argument is the reduction of the projection error to a rational approximation problem. For the resolvent kernel $(tI+A)^{-1}$, the skeleton approximation yields the explicit error formula $1-(t+\\lambda)f_{\\mathrm{skel}}(t,\\lambda)=r(\\lambda)/r(-t)$, and for the exponential kernel the Laplace transform converts the same skeleton error into a bound involving the inverse Laplace transform of $(1/s)p(s)/p(-s)$. The optimal rational functions solving the third Zolotarev problem, that is, the poles of $p(z)/p(-z)$ that minimize $\\max_{[a,b]}|r(z)|/\\min_{[a,b]}|r(-z)|$, provide the pole sets; Möbius transformations map the relevant intervals $[-\\infty,0]\\cup[a,b]$ onto the symmetric Zolotarev configuration. In the Kronecker case the machinery is the tensorized subspace $U\\otimes V$, where the projected matrix retains the same Kronecker-sum form, and the Cauchy-Stieltjes analysis passes through uniformly accurate solution of parameter-dependent Sylvester equations $(tI+A)X_t-X_tB=F$ using rational Galerkin residual bounds.","core_discovery":"For Hermitian positive definite $A$ with spectrum in $[a,b]$, the rational Krylov iterate with Zolotarev poles satisfies $\\|f(A)v-x_\\ell\\|_2\\le C\\rho_{[a,b]}^{\\ell/2}$ when $f$ is Laplace-Stieltjes and $\\|f(A)v-x_\\ell\\|_2\\le C\\rho_{[a,4b]}^{\\ell}$ when $f$ is Cauchy-Stieltjes. For $M=I\\otimes A-B^T\\otimes I$ with $A,-B$ Hermitian positive definite and $v=\\mathrm{vec}(F)$ for low-rank $F$, tensorized rational Krylov spaces give $\\|X-X_\\ell\\|_2\\le C\\rho_{[a,b]}^{\\ell/2}$ for Laplace-Stieltjes functions and $\\|X-X_\\ell\\|_2\\le C\\rho_{[a,2b]}^{\\ell}$ for Cauchy-Stieltjes functions. The constants depend only on $f(0+)$ or $f(2a)$, the spectral ratio $\\kappa=b/a$, and a mild logarithmic factor $\\gamma_{\\ell,\\kappa}$, but not on the density defining $f$. These are explicit, parameter-free convergence rates, and the paper also derives matching singular-value decay bounds for the low-rank solution $X$.","pith_inferences":["A natural stress test is the non-normal case the authors leave open: if the field of values lies in the right half-plane, the same Zolotarev intervals should still yield some convergence, but the sharp rates and the diagonalization trick for the projected problem would need separate treatment.","The bounds can be read as a stopping-criterion blueprint: since the error after $\\ell$ steps is dominated by an explicit power of $\\rho$, one can choose $\\ell$ in advance from the spectrum and target tolerance, and the nested pole sequences let the same choice be revised online.","The singular-value decay theorems imply that solutions of parameter-dependent Sylvester equations inherit low-rank approximability uniformly in the shift parameter, which is exactly what makes low-rank time-stepping for fractional diffusion feasible on tensor grids."],"forward_implications":["Ill-conditioned arguments are where the gain appears: the rate depends on the spectral ratio only through $\\log(4b/a)$, so a widening spectrum costs little extra.","Nested sequences of poles with the same asymptotic rate can be generated by equidistributed sequences, so the method can be run adaptively without fixing $\\ell$ in advance.","For Laplace-Stieltjes functions and rank-$k$ right-hand sides, the solution $X$ of $\\mathrm{vec}(X)=f(I\\otimes A-B^T\\otimes I)\\mathrm{vec}(F)$ has singular values bounded by $O(\\rho_{[a,b]}^{\\ell/2})$ after the $k\\ell$-th, justifying low-rank compression.","For Cauchy-Stieltjes functions of Kronecker-sum matrices, the analogous singular-value decay $O(\\rho_{[a,2b]}^{\\ell})$ holds without the condition-number factor $(1+\\kappa)$.","The projected function evaluation can be done by diagonalizing $A_U$ and $B_V$, costing $O(k^3)$ plus $k^2$ scalar function evaluations rather than $O(k^6)$."],"supporting_citations":[{"why":"supplies the rational Krylov subspace definition and the quasi-optimality inequality linking projection error to rational approximation.","marker":"[18]"},{"why":"states the third Zolotarev problem and the elliptic-function bound giving the pole set and the rate $\\rho_{[a,b]}$.","marker":"[33]"},{"why":"provides the skeleton-approximation technique, the pointwise resolvent error bound, and the equidistributed-sequence construction for nested poles.","marker":"[14]"},{"why":"gives the skeleton approximation error formula for the resolvent kernel used in Theorem 3.4.","marker":"[26]"},{"why":"introduces the tensorized Krylov projection for Kronecker-sum matrix functions and the structure-preserving evaluation of the projected problem.","marker":"[8]"},{"why":"supplies the rational Galerkin residual decomposition for Sylvester equations used to bound the Kronecker Cauchy-Stieltjes error.","marker":"[3]"},{"why":"provides the displacement-structure singular value bounds and the rational approximation results behind the low-rank decay theorems.","marker":"[5]"}],"fun_headline_variants":["Explicit exponential rates for Stieltjes matrix functions","Pole selection with proven convergence for Stieltjes functions","Tensorized Krylov: exponential decay for Kronecker-structured problems","Zolotarev poles give optimal rational Krylov rates","Rational Krylov: tight bounds for Laplace- and Cauchy-Stieltjes"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The convergence rates assume $A$ and $-B$ are symmetric positive definite with spectra exactly contained in a known interval $[a,b]$, and the poles are built from those endpoints; if the interval endpoints are inaccurate, the predicted rates no longer apply.","fun_headline_variants_meta":{"raw":{"variants":["Explicit exponential rates for Stieltjes matrix functions","Pole selection with proven convergence for Stieltjes functions","Tensorized Krylov: exponential decay for Kronecker-structured problems","Zolotarev poles give optimal rational Krylov rates","Rational Krylov: tight bounds for Laplace- and Cauchy-Stieltjes"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000178,"raw_usage":{"total_tokens":1336,"prompt_tokens":1023,"completion_tokens":313,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":639,"completion_tokens_details":{"reasoning_tokens":227}},"tokens_in":639,"tokens_out":313,"duration_ms":27469,"temperature":1.0,"reasoning_tokens":227,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-14T14:56:10.304234+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Pick the 1D Laplacian $A$ of size $10^5$ with known spectrum in $[a,b]\\approx[1/n,4]$, set $f(z)=e^{-z}$, and run rational Krylov with the paper's Zolotarev poles for $\\ell=5,10,15,20$. The log-error should fall with slope at least $-\\pi^2/(2\\log(4b/a))$; a markedly flatter slope, or failure of the error to decay once past the constant regime, would falsify the bound.","supporting_citations":[{"cited_title":"Zolotarev , Application of elliptic functions to questions of function s deviating least and most from zero, Zap","cited_arxiv_id":null,"evidence_quote":"states the third Zolotarev problem and the elliptic-function bound giving the pole set and the rate $\\rho_{[a,b]}$."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"gives the skeleton approximation error formula for the resolvent kernel used in Theorem 3.4."},{"cited_title":"Benzi and V","cited_arxiv_id":null,"evidence_quote":"introduces the tensorized Krylov projection for Kronecker-sum matrix functions and the structure-preserving evaluation of the projected problem."}],"review_version":1}