Pith. sign in

REVIEW 5 minor 42 references

The mmatrix toolbox: componentwise accurate algorithms for M-matrices with triplet representation

T0 review · 0 major / 5 minor · reviewed 2026-07-12 · grok-4.5

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read A Matlab toolbox uses subtraction-free GTH algorithms on triplet M-matrices to keep every entry accurate even when the matrix is severely ill-conditioned.

desk verdict Solid Lapack-style software packaging of known GTH accuracy results; the blocked/recursive kernels and OO Matlab class are the real additions, and the experiments back the claims. read the letter →

arxiv 2607.05437 v1 pith:YKJUBLMH submitted 2026-07-03 cs.MS

Bruno Iannazzo , Mehdi Najafi Kalyani , Federico Poloni This is my paper · ORCID
classification cs.MS MSC 65F3065Y1565F0515B48
keywords M-matrixmmatrix-toolboxcomponentwiseerrorGTHalgorithmtripletrepresentationLapack-stylefactorizationaccurateSVDmatrixsquareroot
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

The paper introduces the mmatrix toolbox: a Matlab package whose Fortran core implements unblocked, recursive and blocked variants of the GTH algorithm for M-matrices that are given by a left or right triplet. A triplet stores the off-diagonal entries together with two nonnegative vectors that reconstruct the diagonal without cancellation; the same representation lets every subsequent elimination step stay subtraction-free. Consequently the toolbox returns componentwise-accurate LU factors, linear-system solutions (nonnegative right-hand sides), inverses, Schur complements, singular values, square roots and solutions of associated Riccati equations, even for matrices whose ordinary condition numbers are huge. The Matlab layer wraps these kernels in an object that overloads ordinary linear-algebra operations, so users can keep working with familiar syntax while the library silently preserves the triplet. Numerical experiments confirm that the new routines stay near machine precision in every entry while standard Lapack-style solvers lose all correct digits in the small ones.

What carries the argument

The GTH algorithm: a subtraction-free Gaussian elimination that, at every pivot step, recovers the diagonal entry from the current triplet by a formula that only adds and multiplies nonnegative quantities, then updates the triplet of the Schur complement so that the invariant is preserved.

What would settle it

Run the toolbox and ordinary Lapack solvers on a family of triplet M-matrices whose condition number grows without bound (for example the finite-element matrix of Example 3 with successively smaller sigma) and check whether the componentwise relative error of every entry of the inverse or of the singular-value vector stays near machine epsilon while the classical solvers lose all digits.

Watch

Extended reading notes

Core claim

When an M-matrix is supplied through a left or right triplet, the GTH algorithm (and its blocked and recursive Lapack-style variants) can compute LU factors, inverses, Schur complements and several matrix functions without any subtractive cancellation, thereby guaranteeing high componentwise accuracy even for severely ill-conditioned instances.

Load-bearing premise

A usable left or right triplet must already be known analytically; recovering one from the matrix entries alone is generally an ill-conditioned operation.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

0 major / 5 minor

Summary. The paper presents the mmatrix toolbox, a Matlab package with a Lapack-style Fortran core for componentwise-accurate computations on M-matrices given by left or right triplet representations. It implements unblocked, recursive (cache-oblivious), and blocked GTH algorithms for LU factorization, linear systems with nonnegative right-hand sides, inversion, Schur complements, singular values (via accurate LDU + Jacobi SVD), matrix square roots (Cyclic Reduction with a safety factor), MARE solutions (SDA/ADDA with auxiliary triplets), and the smallest eigenvalue (shifted inverse iteration). Componentwise perturbation theory is collected and mildly extended (Theorem 3 and Appendix Theorem 5), and an object-oriented mmatrix class overloads standard Matlab operations while preserving triplets when possible. Extensive numerical experiments (Examples 3–8) show near-machine componentwise accuracy on ill-conditioned problems where standard Matlab/Lapack routines lose digits, with timings competitive for the blocked/recursive kernels.

Significance. If the claims hold, the work supplies a practical, reusable software artifact that makes classical subtraction-free GTH technology and its extensions usable for a range of M-matrix problems that arise in Markov chains, networks, and discretizations. Strengths include the Lapack-style Fortran kernels (left and right triplets, multiple blocking strategies), the systematic collection of componentwise perturbation results with a mild generalization allowing u-perturbations, the safety-factor argument that keeps non-subtraction-free steps (square root, SDA) under control, and reproducible numerical evidence that componentwise errors stay near machine precision while standard solvers do not. The open Matlab interface lowers the barrier for users who already possess analytic triplets. The main practical limitation—the need for a usable a-priori triplet—is already stated by the authors and does not undermine the algorithmic contribution.

minor comments (5)
  1. In §2.6.4 the SVD path uses diagonal pivoting on the recovered diagonal entries rather than complete pivoting when u is not the all-ones vector; the text correctly notes that κ(L) and κ(Û) may then be larger. A short remark or numerical check quantifying the degradation for the ill-conditioned-u case of Example 6 would help users decide when the method remains reliable.
  2. Table 1 and §3.7 list work-buffer sizes; it would be useful to state explicitly whether the recursive right-triplet kernel (dmmtrf2) re-uses the same n-length buffer across recursion levels or allocates additional temporary space, for readers who wish to call the Fortran routines directly.
  3. The default safety factor γ = 1.5 is used for square-root and SDA; a one-sentence sensitivity note (or a pointer to the earlier literature) on how larger γ trades accuracy for slower quadratic convergence would make the parameter choice more transparent.
  4. Minor typographical inconsistencies appear (e.g., “Themmatrixtoolbox” / “mmatrix-toolbox” / “mmatrix toolbox” in the title and abstract; “themmatrixclass” in §4). Uniform naming would improve readability.
  5. In Algorithm 1 the update v ← u ⊙ (p − min(p)) is subtraction-free only after the min is taken; a brief parenthetical that the resulting v remains nonnegative would reassure readers scanning the pseudocode.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity; GTH accuracy claims rest on classical subtraction-free analysis, independent perturbation theory, and external numerical benchmarks.

full rationale

The paper's central claims concern Lapack-style implementations (unblocked/recursive/blocked GTH kernels for left/right triplets, plus derived solvers for systems, inverses, Schur complements, SVD, square roots, and MAREs) that deliver componentwise accuracy when a usable triplet is supplied. These rest on the well-known subtraction-free property of the GTH algorithm (cited to Grassmann–Taksar–Heyman and Alfa–Xue–Ye) together with the componentwise perturbation bounds collected in Theorem 3 (proved in the appendix from the all-minors matrix-tree lemma of Chaiken and related literature). No parameters are fitted to data and then re-presented as predictions; the numerical experiments of §5 generate independent test matrices (finite-element, random triplets) and compare against high-precision (vpa) or exact references. Self-citations (e.g., [5] for the square-root iteration, [14] for multilinear PageRank, [16] for fluid queues) supply algorithmic building blocks whose correctness is independently established; they are not load-bearing uniqueness theorems that close a circular loop. The only acknowledged limitation—the need for an a-priori triplet—is stated explicitly (Figure 1, §2.2) and does not render the algorithmic claims circular. The derivation chain is therefore self-contained against external theory and benchmarks.

Assumptions & free parameters 2 free parameters · 3 assumptions · 0 invented entities

The central accuracy claims rest on classical M-matrix theory (nonnegativity of inverses, existence of triplets under irreducibility or nonsingularity) plus the known subtraction-free property of GTH. The only free parameters are algorithmic safety factors and block sizes chosen by the user; no new physical or mathematical entities are postulated.

free parameters (2)
  • safety factor γ (default 1.5) = 1.5
    Used in square-root and related algorithms to keep diagonal subtractions ‘safe’; chosen by hand to bound relative-error amplification by 1/(γ-1).
  • block size n_b
    User-supplied parameter controlling the blocked GTH variants; affects performance but not correctness.
assumptions (3)
  • standard math An M-matrix possesses a right (resp. left) triplet if and only if, in its Frobenius normal form, every singular diagonal block has zero off-block entries in its row (resp. column).
    Theorem 1, taken from the literature and used throughout to guarantee that the input objects are well-defined.
  • standard math Componentwise relative perturbations of a triplet of size ε produce componentwise relative perturbations of the inverse of size O(n)ε (first-order).
    Theorem 3 / Appendix; underpins the claim that GTH is forward-stable in the componentwise sense.
  • domain assumption The GTH update never performs a subtraction of two quantities of the same sign when a triplet is maintained.
    Classical property of the GTH algorithm (Alfa–Xue–Ye et al.); the whole toolbox relies on it.

how reviews work

0 comments
Cite this review

Pith. "Pith review of The mmatrix toolbox: componentwise accurate algorithms for M-matrices with triplet representation." pith.science (2026). https://pith.science/paper/YKJUBLMH

@misc{pith2026260705437,
  author       = {Pith},
  title        = {Pith review of: The mmatrix toolbox: componentwise accurate algorithms for M-matrices with triplet representation},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/YKJUBLMH}},
  note         = {Machine review of arXiv:2607.05437}
}
read the original abstract

We introduce the mmatrix toolbox, a Matlab software package for componentwise accurate computations with M-matrices described through left or right triplet representations. The core of the toolbox is a Fortran implementation, in the Lapack style, of the unblocked, recursive, and blocked versions of the GTH algorithm and its applications for the accurate computation of the solution of linear systems with M-matrix coefficient and nonnegative right-hand side, the LU factorization of an M-matrix and its inverse. These algorithms avoid subtractive cancellation; this property ensures high componentwise accuracy even for ill-conditioned problems. The toolbox contains also accurate algorithms for related problems, such as computing the Schur complement, the singular values, the square root of an M-matrix, and the solution of nonsymmetric algebraic Riccati equations associated with M-matrices. The Matlab interface is based on an object-oriented implementation allowing one to use standard Matlab operations on M-matrices with triplet representation.

Figures

Figures reproduced from arXiv: 2607.05437 by the authors.

Figure 1
Figure 1. Conditioning of transformations between various quantities, for an M-matrix [PITH_FULL_IMAGE:figures/full_fig_p005_1.png] view at source ↗
Figure 2
Figure 2. shows that all GTH-based implementations preserve near machine-precision accuracy in both the normwise and componentwise errors, while the standard LU fac￾torization exhibits large componentwise inaccuracies as the problem becomes more ill￾conditioned [PITH_FULL_IMAGE:figures/full_fig_p018_2.png] view at source ↗
Figure 3
Figure 3. Example 5: Results for matrix inversion with [PITH_FULL_IMAGE:figures/full_fig_p020_3.png] view at source ↗
Figures from the paper (2 more)
Figure 4
Figure 4. Figure 4: Example 6: Singular values for random M-matrices of size [PITH_FULL_IMAGE:figures/full_fig_p021_4.png]
Figure 5
Figure 5. Figure 5: Example 8: Accuracy comparison for sqrtm. The package is easy-to-use due to the Matlab interface, efficient due to the Fortran implementation of the core routines, and flexible enough to be extended with future algorithms. A On Theorem 3 A result that can be used to pr…

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

42 extracted references · 6 canonical work pages

  1. [1]

    , TITLE =

    Berman, Abraham and Plemmons, Robert J. , TITLE =. 1994 , PAGES =

    1994

  2. [2]

    Linear Algebra Appl

    Computing the singular value decomposition with high relative accuracy , author=. Linear Algebra Appl. , volume=. 1999 , DOI =

    1999

  3. [3]

    Demmel, James and Koev, Plamen , fjournal=. Numer. Math. , volume=. 2004 , DOI =

    2004

  4. [4]

    Linear Algebra Appl

    Kuo, I Wen , TITLE =. Linear Algebra Appl. , FJOURNAL =. 1977 , NUMBER =

    1977

  5. [5]

    and Iannazzo, Bruno and Meini, Beatrice and Meng, Jie , year=

    Bini, Dario A. and Iannazzo, Bruno and Meini, Beatrice and Meng, Jie , year=. Component-wise accurate computation of the square root of an. 2605.21679 , archivePrefix=

    arXiv

  6. [6]

    Regenerative analysis and steady state distributions for

    Grassmann, Winfried K and Taksar, Michael I and Heyman, Daniel P , fjournal=. Regenerative analysis and steady state distributions for. Oper. Res. , volume=. 1985 , publisher=

    1985

  7. [7]

    Accurate solutions of

    Xue, Jungong and Xu, Shufang and Li, Ren-Cang , FJournal =. Accurate solutions of. Numer. Math. , volume=. 2012 , publisher=

    2012

  8. [8]

    Royal Society Open Science , volume=

    Comparison of social structures within cities of very different sizes , author=. Royal Society Open Science , volume=. 2016 , publisher=

    2016

Show all 42 references
  1. [9]

    2006 , URL =

    Pajek datasets , author=. 2006 , URL =

    2006

  2. [10]

    Numerical computing with

    Moler, Cleve B , year=. Numerical computing with

  3. [11]

    Multilinear

    Gleich, David F and Lim, Lek-Heng and Yu, Yongyang , FJournal=. Multilinear. SIAM J. Matrix Anal. Appl. , volume=. 2015 , publisher=

    2015

  4. [12]

    Proceedings of the 2015 SIAM International Conference on Data Mining , pages=

    Tensor spectral clustering for partitioning higher-order network structures , author=. Proceedings of the 2015 SIAM International Conference on Data Mining , pages=. 2015 , organization=

    2015

  5. [13]

    and Gleich, David F

    Benson, Austin R. and Gleich, David F. and Lim, Lek-Heng , title =. SIAM Rev. , issn =. 2017 , language =. doi:10.1137/16M1074023 , keywords =

    2017 doi

  6. [14]

    The University of

    Davis, Timothy A and Hu, Yifan , Fjournal=. The University of. ACM Trans. Math. Softw. , volume=. 2011 , publisher=

    2011

  7. [15]

    Statistical Science , volume=

    The mixture transition distribution model for high-order Markov chains and non-Gaussian time series , author=. Statistical Science , volume=. 2002 , publisher=

    2002

  8. [16]

    Estimation and modelling repeated patterns in high order

    Raftery, Adrian and Tavar. Estimation and modelling repeated patterns in high order. Journal of the Royal Statistical Society Series C: Applied Statistics , volume=. 1994 , publisher=

    1994

  9. [17]

    Linear and Multilinear Algebra , volume=

    On the limiting probability distribution of a transition probability tensor , author=. Linear and Multilinear Algebra , volume=. 2014 , publisher=

    2014

  10. [18]

    Iterative

    Saad, Yousef , year=. Iterative

  11. [19]

    Varga, Richard S. , year=. Matrix

  12. [20]

    Introduction to

    Stoer, Josef and Bulirsch, Roland and Bartels, R and Gautschi, Walter and Witzgall, Christoph , volume=. Introduction to. 1980 , publisher=

    1980

  13. [21]

    Gleich, David F and Yu, Yongyang , title=

  14. [22]

    Entrywise perturbation theory for diagonally dominant

    Alfa, Attahiru Sule and Xue, Jungong and Ye, Qiang , FJournal =. Entrywise perturbation theory for diagonally dominant. Numer. Math. , volume=. 2002 , publisher=

    2002

  15. [23]

    Accurate computation of the smallest eigenvalue of a diagonally dominant

    Alfa, Attahiru Sule and Xue, Jungong and Ye, Qiang , fjournal=. Accurate computation of the smallest eigenvalue of a diagonally dominant. Math. Comp. , volume=. doi:10.1090/S0025-5718-01-01325-4 , year=

    doi

  16. [24]

    SIAM Journal on Mathematics of Data Science , volume=

    Ergodicity coefficients for higher-order stochastic processes , author=. SIAM Journal on Mathematics of Data Science , volume=. 2020 , publisher=

    2020

  17. [25]

    Perron-based algorithms for the multilinear

    Meini, Beatrice and Poloni, Federico , FJournal=. Perron-based algorithms for the multilinear. Numer. Linear Algebra Appl. , volume=. 2018 , publisher=

    2018

  18. [26]

    Linear Algebra Appl

    Poloni, Federico , Title =. Linear Algebra Appl. , ISSN =. 2013 , Language =. doi:10.1016/j.laa.2011.05.036 , Keywords =

    2013 doi

  19. [27]

    Xue, Jungong and Xu, Shufang and Li, Ren-Cang , Title =. Numer. Math. , ISSN =. 2012 , DOI =

    2012

  20. [28]

    Extrapolation methods for fixed-point multilinear

    Cipolla, Stefano and Redivo-Zaglia, Michela and Tudisco, Francesco , FJournal=. Extrapolation methods for fixed-point multilinear. Numer. Linear Algebra Appl. , volume=. 2020 , publisher=

    2020

  21. [29]

    Anderson accelerated fixed-point iteration for multilinear

    Lai, Fuqi and Li, Wen and Peng, Xiaofei and Chen, Yannan , FJournal=. Anderson accelerated fixed-point iteration for multilinear. Numer. Linear Algebra Appl. , volume=. 2023 , publisher=

    2023

  22. [30]

    Relative-error bounds for the

    O'Cinneide, Colm Art , FJournal =. Relative-error bounds for the. Numer. Math. , volume=. 1996 , publisher=

    1996

  23. [31]

    Abdesselam, Abdelmalek , FJournal=. The. Adv. Appl. Math. , volume=. 2004 , publisher=

    2004

  24. [32]

    A combinatorial proof of the all minors matrix tree theorem , author=. SIAM J. Alg. Disc. Meth. , volume=. 1982 , publisher=

    1982

  25. [33]

    Newton's method in floating point arithmetic and iterative refinement of generalized eigenvalue problems , fjournal =

    Tisseur, Fran. Newton's method in floating point arithmetic and iterative refinement of generalized eigenvalue problems , fjournal =. SIAM J. Matrix Anal. Appl. , issn =. 2001 , language =. doi:10.1137/S0895479899359837 , keywords =

    2001 doi

  26. [34]

    , title =

    Higham, Nicholas J. , title =. 2002 , publisher =. doi:10.1137/1.9780898718027 , keywords =

    2002 doi

  27. [35]

    and Poloni, Federico , title =

    Nguyen, Giang T. and Poloni, Federico , title =. Numer. Math. , issn =. 2015 , language =. doi:10.1007/s00211-014-0675-4 , keywords =

    2015 doi

  28. [36]

    and Sun, Ji-guang , title =

    Stewart, Gilbert W. and Sun, Ji-guang , title =. 1990 , publisher =

    1990

  29. [37]

    Accuracy and componentwise accuracy in multilinear

    Mehdi Najafi Kalyani and Federico Poloni , year=. Accuracy and componentwise accuracy in multilinear. 2506.18356 , archivePrefix=

    arXiv

  30. [38]

    Linear Algebra Appl

    Guo, Chun-Hua , title =. Linear Algebra Appl. , issn =. 2013 , language =. doi:10.1016/j.laa.2013.08.018 , keywords =

    2013 doi

  31. [39]

    Linear Algebra Appl

    Guo, Chun-Hua and Lu, Di , title =. Linear Algebra Appl. , issn =. 2016 , language =. doi:10.1016/j.laa.2015.11.024 , keywords =

    2016 doi

  32. [40]

    and Iannazzo, Bruno and Meini, Beatrice , TITLE =

    Bini, Dario A. and Iannazzo, Bruno and Meini, Beatrice , TITLE =. 2012 , PAGES =

    2012

  33. [41]

    2018 , publisher =

    Huang, Tsung-Ming and Li, Ren-Cang and Lin, Wen-Wei , title =. 2018 , publisher =. doi:10.1137/1.9781611975369 , keywords =

    2018 doi

  34. [42]

    and Prokop, Harald and Ramachandran, Sridhar , title =

    Frigo, Matteo and Leiserson, Charles E. and Prokop, Harald and Ramachandran, Sridhar , title =. Proceedings of the 40th Annual Symposium on Foundations of Computer Science , series =. 1999 , pages =

    1999

Pith tools

Reviewed July 12, 2026 · model on record in the stance chip above.