pith. sign in

arxiv: 2604.23306 · v1 · submitted 2026-04-25 · 🧮 math.NA · cs.NA

Summation-by-parts operators for general function spaces: optimal nodes

Pith reviewed 2026-05-08 07:27 UTC · model grok-4.3

classification 🧮 math.NA cs.NA
keywords summation-by-partsGauss-Lobatto quadratureoptimal nodesgeneral function spacesnumerical methodsboundary conditionsquadrature rulesinitial boundary value problems
0
0 comments X

The pith

Generalized Gauss-Lobatto quadrature gives the smallest summation-by-parts operators for any chosen function space.

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

The paper shows that the classical optimality of Gauss-Lobatto nodes for polynomial-based summation-by-parts operators carries over to arbitrary function spaces. For any prescribed space the associated generalized Gauss-Lobatto quadrature supplies the unique set of nodes and weights that produces an SBP operator of minimal dimension while retaining the boundary closure property. This matters because smaller operators reduce the cost of high-order discretizations of initial-boundary-value problems without sacrificing stability. The authors supply a stable algorithm to compute the rules and demonstrate their accuracy on both polynomial and non-polynomial spaces.

Core claim

For any prescribed function space an associated generalised Gauss-Lobatto quadrature provides the optimal nodes and weights for the SBP formulation, yielding an operator whose dimension is minimal among all possible choices.

What carries the argument

The generalised Gauss-Lobatto quadrature rule, defined so that it integrates exactly the products of functions from the chosen space while forcing nodes at the endpoints.

If this is right

  • The minimal-dimension SBP operator is obtained directly from the quadrature nodes and weights without further search or optimization.
  • The same construction works for non-polynomial spaces such as trigonometric or exponential bases.
  • The resulting operators can be inserted into existing SBP-based codes for initial-boundary-value problems with no change in the stability proof.
  • Computational cost scales with the dimension of the space rather than with an arbitrary choice of grid points.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

  • Users could design custom function spaces adapted to known solution features and obtain the corresponding optimal nodes automatically.
  • The approach may reduce the number of degrees of freedom needed in high-order simulations of problems with localized features.
  • Extension to time-dependent or nonlinear problems would require only that the spatial operator retain the SBP property at each step.

Load-bearing premise

A generalized Gauss-Lobatto quadrature exists, can be computed stably, and produces an SBP operator of strictly minimal dimension for the given function space.

What would settle it

A concrete counter-example: a function space together with a set of nodes and weights that produces a valid closed SBP operator whose size is smaller than the size obtained from the corresponding generalized Gauss-Lobatto rule.

read the original abstract

Gauss-Lobatto quadrature nodes and weights are optimal for closed summation-by-parts (SBP) formulations based on polynomial approximation spaces in the sense that for a prescribed function space they yield an SBP operator of minimal dimension. We show that the same principle extends to general (possibly non-polynomial) function spaces: an associated generalised Gauss-Lobatto quadrature provides the optimal nodes and weights for the SBP formulation. We present an algorithm for computing these quadrature rules, demonstrate their accuracy and efficiency across a range of function spaces, and illustrate their use in solving initial boundary value problems.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

3 major / 2 minor

Summary. The manuscript claims that Gauss-Lobatto quadrature nodes and weights are optimal for closed summation-by-parts (SBP) operators on polynomial spaces in the sense of yielding minimal-dimension operators for a prescribed space, and that this optimality principle extends to general (possibly non-polynomial) function spaces via an associated generalized Gauss-Lobatto quadrature. It presents an algorithm for computing the nodes and weights, demonstrates accuracy and efficiency on selected spaces, and illustrates application to initial-boundary-value problems.

Significance. If the result holds, the work would extend a standard optimality result for SBP discretizations beyond polynomials, enabling minimal-dimension operators for specialized approximation spaces in numerical PDEs. The algorithmic contribution and numerical demonstrations on non-polynomial examples are concrete strengths that support practical utility.

major comments (3)
  1. [Abstract and §2] Abstract and opening of §2: the claim that a generalized Gauss-Lobatto quadrature 'provides the optimal nodes and weights for the SBP formulation' for arbitrary prescribed function spaces is not supported by a general existence result; the analysis establishes conditional optimality (when the quadrature exists) but reduces the construction to a nonlinear system whose solvability is space-dependent and is not proven in general.
  2. [§3] §3 (algorithm description): the procedure for solving the nonlinear system for nodes/weights does not include analysis of existence, uniqueness, or stability conditions (e.g., positive weights and forced boundary nodes); this is load-bearing because the central optimality claim for general spaces rests on the algorithm always succeeding.
  3. [Numerical results] Numerical results section: demonstrations are confined to selected examples where the quadrature exists and is stable; no counter-examples, failure cases, or bounds on the class of admissible spaces are provided, leaving the generality of the optimality claim unverified.
minor comments (2)
  1. [§2] Notation for the generalized quadrature weights and the SBP operator matrix could be made more explicit when first introduced to avoid confusion with the classical polynomial case.
  2. [Numerical results] Figure captions in the examples section would benefit from listing the precise function spaces and polynomial degrees (or equivalent) used in each test.

Simulated Author's Rebuttal

3 responses · 0 unresolved

We thank the referee for the thorough and constructive report. The comments correctly identify areas where the manuscript's claims require clarification regarding scope and supporting analysis. We address each major comment point by point below, with revisions planned to improve precision without altering the core contributions.

read point-by-point responses
  1. Referee: [Abstract and §2] Abstract and opening of §2: the claim that a generalized Gauss-Lobatto quadrature 'provides the optimal nodes and weights for the SBP formulation' for arbitrary prescribed function spaces is not supported by a general existence result; the analysis establishes conditional optimality (when the quadrature exists) but reduces the construction to a nonlinear system whose solvability is space-dependent and is not proven in general.

    Authors: We agree that the manuscript does not establish a general existence result for the generalized Gauss-Lobatto quadrature across all possible function spaces. The optimality result is indeed conditional: when such a quadrature exists and satisfies the required moment conditions, it yields the minimal-dimension closed SBP operator for the prescribed space. The construction is reduced to a nonlinear algebraic system whose solvability is space-dependent. In the revised version we will modify the abstract and the opening paragraphs of §2 to state this conditional character explicitly, replacing the current phrasing with language that emphasizes existence as a prerequisite. This revision will align the claims more precisely with the analysis provided. revision: yes

  2. Referee: [§3] §3 (algorithm description): the procedure for solving the nonlinear system for nodes/weights does not include analysis of existence, uniqueness, or stability conditions (e.g., positive weights and forced boundary nodes); this is load-bearing because the central optimality claim for general spaces rests on the algorithm always succeeding.

    Authors: The referee is correct that §3 presents a practical numerical procedure without a supporting theoretical analysis of existence, uniqueness, or stability of the nonlinear solver. We will add a dedicated subsection to §3 that discusses the numerical behavior we have observed, including the enforcement of boundary nodes, the requirement of positive weights, and practical strategies (such as initial guesses and regularization) that promote reliable convergence in the examples considered. While a complete theoretical characterization of solvability for arbitrary spaces lies beyond the present scope, the added discussion will make the load-bearing assumptions transparent and will note that the algorithm is intended for spaces where the underlying quadrature problem is well-posed. revision: partial

  3. Referee: [Numerical results] Numerical results section: demonstrations are confined to selected examples where the quadrature exists and is stable; no counter-examples, failure cases, or bounds on the class of admissible spaces are provided, leaving the generality of the optimality claim unverified.

    Authors: We acknowledge that the numerical section focuses on successful cases. In the revision we will enlarge the set of examples to include additional function spaces that test the limits of the approach (for instance, spaces with near-singular behavior or higher effective degrees). We will also insert a short discussion of observed limitations, including situations in which the nonlinear system may fail to produce a real solution with positive weights, together with practical diagnostics that indicate when the method is likely to succeed. These additions will provide readers with clearer guidance on the range of admissible spaces while preserving the concrete demonstrations already present. revision: yes

Circularity Check

0 steps flagged

No significant circularity detected in derivation chain

full rationale

The paper defines optimality of nodes/weights for SBP operators via the independent criterion of minimal operator dimension for a prescribed function space, then shows that an associated generalized Gauss-Lobatto quadrature satisfies this criterion when it exists. This is an extension of a known polynomial-space result rather than a self-referential definition or fitted prediction. No load-bearing self-citation chains, uniqueness theorems imported from the authors' prior work, or ansatzes smuggled via citation are required for the central claim; the work supplies an algorithm and example demonstrations instead. Concerns about general existence for arbitrary spaces pertain to completeness rather than circularity in the presented derivation.

Axiom & Free-Parameter Ledger

0 free parameters · 1 axioms · 0 invented entities

The central claim rests on the existence and computability of a generalized Gauss-Lobatto quadrature for arbitrary spaces that satisfies the SBP summation property.

axioms (1)
  • domain assumption A generalized Gauss-Lobatto quadrature rule exists for the chosen function space and satisfies the discrete summation-by-parts identity.
    Invoked to guarantee that the constructed nodes and weights produce a valid closed SBP operator of minimal dimension.

pith-pipeline@v0.9.0 · 5393 in / 1132 out tokens · 37298 ms · 2026-05-08T07:27:03.748811+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

62 extracted references · 62 canonical work pages · 1 internal anchor

  1. [1]

    Introduction.Summation-by-parts (SBP) formulations [12, 53] combined with weak imposition of boundary conditions ensure semi-discrete stability for well- posed linear initial boundary value problems (IBVPs) [22, 34, 35]. The SBP framework emulates the principle of continuous integration-by-parts (IBP) at the discrete level, which is the first crucial step...

  2. [2]

    Preliminaries. 2.1. FSBP operators.Following Glaubitzet al.[18], we define anF-based summation-by-parts (FSBP) operator as follows. Definition2.1 (FSBP operators).LetF ⊂C 1([a, b])be a finite dimensional function space andx= [x 1, . . . , xn]⊤ ⊂[a, b]be a set ofndistinct nodes. 1 An operator D=P −1Qis anF-based summation-by-parts (FSBP) operatorif i.Df(x)...

  3. [3]

    nearly Gaussian quadratures

    Generalised Gauss quadrature.In the polynomial case, Gauss–Legendre and Gauss–Legendre–Lobatto quadrature rules underpin the construction of SBP op- erators. For more general function spaces, a natural analogue arises: generalised Gauss quadrature (GGQ) rules. Since ann-point quadrature has 2nfree parameters (nodes and weights), it is reasonable to expect...

  4. [4]

    Applications to IBVPs.To illustrate the efficiency of our development, we will consider two examples of IBVPs, both in a multi-domain setting. 4.1. Advection equation.Consider the constant coefficient problem ut +au x = 0,0≤x≤1, t >0, u=g L(t), x= 0, t≥0,(4.1) u=f(x),0≤x≤1, t= 0, wherea >0. The boundary datag L and initial datafare compatible such that a ...

  5. [5]

    Lastly, the temporal integration leads to an energy estimate and hence stability

    We chooseσ L = 0 in equation (4.5) which leads to an energy rate similar to (4.2) but with damping terms. Lastly, the temporal integration leads to an energy estimate and hence stability. Here,||u|| 2 P =u T Puis the discreteL 2-equivalent norm. We will compute and compare approximate solutions to equation (4.1) for the following constructions of the deri...

  6. [6]

    Therefore, (4.12) with (4.13) leads to an energy rate which is similar to (4.9) but with additional damping terms. Next, we compare approximate solutions to equation (4.7) using (4.10)–(4.11) with derivative operators based on exponential basis functions with equispaced and optimised nodes, respectively. We also implement polynomial bases with equispaced ...

  7. [7]

    Conclusion.We have studied the construction of diagonal-norm FSBP oper- ators for general (non-polynomial) approximation spaces by exploiting the equivalence between such operators and positive quadrature rules exact forG= (F F) ′. WhenG forms a Tchebyshev system, the theory of generalised Gauss quadrature implies ex- istence, uniqueness, and positivity o...

  8. [8]

    AI was not used in the generation of any code, numerical results, or figures

    Acknowledgements.Generative AI (ChatGPT) was used to assist in editing the authors’ text for spelling, grammar, and structure in this work. AI was not used in the generation of any code, numerical results, or figures. References. N. HALE, C. HARLEY, P. NCHUPANG, AND J. NORDSTR ¨OM

  9. [9]

    Abgrall, J

    R. Abgrall, J. Nordstr¨ om, P.¨Offner, and S. Tokareva. Analysis of the SBP-SAT stabilization for finite element methods Part I: Linear problems.J. Sci. Comput., 85(2):Paper No. 43, 29, 2020

  10. [10]

    Abgrall, J

    R. Abgrall, J. Nordstr¨ om, P.¨Offner, and S. Tokareva. Analysis of the SBP-SAT stabilization for finite element methods part II: Entropy stability.Commun. Appl. Math. Comput., 5(2):573–595, 2023

  11. [11]

    Berstein

    S. Berstein. Sur les formules de quadrature de Cotes et Tchebycheff. InComptes Rendus (Doklady) de l’Academie des Sciences de l’URSS, volume 14, pages 323– 326, 1937

  12. [12]

    Bremer, Z

    J. Bremer, Z. Gimbutas, and V. Rokhlin. A nonlinear optimization procedure for generalized Gaussian quadratures.SIAM J. Sci. Comput., 32(4):1761–1788, 2010

  13. [13]

    M. H. Carpenter, T. C. Fisher, E. J. Nielsen, and S. H. Frankel. Entropy sta- ble spectral collocation schemes for the Navier-Stokes equations: discontinuous interfaces.SIAM J. Sci. Comput., 36(5):B835–B867, 2014

  14. [14]

    M. H. Carpenter and D. Gottlieb. Spectral methods on arbitrary grids.J. Comput. Phys., 129(1):74–86, 1996

  15. [15]

    M. H. Carpenter, D. Gottlieb, and S. Abarbanel. Time-stable boundary condi- tions for finite-difference schemes solving hyperbolic systems: methodology and application to high-order compact schemes.J. Comput. Phys., 111(2):220–236, 1994

  16. [16]

    M. H. Carpenter, J. Nordstr¨ om, and D. Gottlieb. A stable and conservative interface treatment of arbitrary spatial accuracy.J. Comput. Phys., 148(2):341– 365, 1999

  17. [17]

    Castonguay, D

    P. Castonguay, D. M. Williams, P. E. Vincent, and A. Jameson. Energy stable flux reconstruction schemes for advection-diffusion problems.Comput. Methods Appl. Mech. Engrg., 267:400–417, 2013

  18. [18]

    Cheng, V

    H. Cheng, V. Rokhlin, and N. Yarvin. Nonlinear optimization, quadrature, and interpolation.SIAM J. Optim., 9(4):901–923, 1999. Dedicated to John E. Dennis, Jr., on his 60th birthday

  19. [19]

    D. C. Del Rey Fern´ andez, P. D. Boom, and D. W. Zingg. A generalized frame- work for nodal first derivative summation-by-parts operators.J. Comput. Phys., 266:214–239, 2014

  20. [20]

    D. C. Del Rey Fern´ andez, J. E. Hicken, and D. W. Zingg. Review of summation- by-parts operators with simultaneous approximation terms for the numerical so- lution of partial differential equations.Comput. & Fluids, 95:171–196, 2014

  21. [21]

    T. A. Driscoll, N. Hale, and L. N. Trefethen. Chebfun guide, 2014

  22. [22]

    G. J. Gassner. A skew-symmetric discontinuous Galerkin spectral element dis- cretization and its relation to SBP-SAT finite difference methods.SIAM J. Sci. Comput., 35(3):A1233–A1253, 2013

  23. [23]

    Gautschi

    W. Gautschi. A survey of Gauss-Christoffel quadrature formulae. In P. L. Butzer and F. Feh´ er, editors,E. B. Christoffel: The Influence of His Work on Mathe- matics and the Physical Sciences, pages 72–147. Birkh¨ auser Basel, Basel, 1981

  24. [24]

    Glaubitz, S.-C

    J. Glaubitz, S.-C. Klein, J. Nordstr¨ om, and P. ¨Offner. Multi-dimensional summation-by-parts operators for general function spaces: theory and construc- tion.J. Comput. Phys., 491:Paper No. 112370, 23, 2023

  25. [25]

    Glaubitz, S.-C

    J. Glaubitz, S.-C. Klein, J. Nordstr¨ om, and P. ¨Offner. Summation-by-parts op- erators for general function spaces: the second derivative.J. Comput. Phys., 504:Paper No. 112889, 20, 2024

  26. [26]

    Glaubitz, J

    J. Glaubitz, J. Nordstr¨ om, and P. ¨Offner. Summation-by-parts operators for SBP OPERATORS FOR GENERAL FUNCTION SPACES: OPTIMAL NODES general function spaces.SIAM J. Numer. Anal., 61(2):733–754, 2023

  27. [27]

    Glaubitz, J

    J. Glaubitz, J. Nordstr¨ om, and P.¨Offner. Energy-stable global radial basis func- tion methods on summation-by-parts form.J. Sci. Comput., 98(1):Paper No. 30, 28, 2024

  28. [28]

    Glaubitz, J

    J. Glaubitz, J. Nordstr¨ om, and P. ¨Offner. An optimization-based construction procedure for function space based summation-by-parts operators on arbitrary grids.J. Sci. Comput., 105(83), 2025

  29. [29]

    Gong and J

    J. Gong and J. Nordstr¨ om. Interface procedures for finite difference approxima- tions of the advection-diffusion equation.J. Comput. Appl. Math., 236(5):602– 620, 2011

  30. [30]

    Gustafsson, H.-O

    B. Gustafsson, H.-O. Kreiss, and J. Oliger.Time dependent problems and differ- ence methods. Pure and Applied Mathematics (New York). John Wiley & Sons, Inc., New York, 1995. A Wiley-Interscience Publication

  31. [31]

    N. Hale. Companion code for this paper. https://github.com/nickhale/fsbpopt

  32. [32]

    J. S. Hesthaven and D. Gottlieb. A stable penalty method for the compressible Navier-Stokes equations. I. Open boundary conditions.SIAM J. Sci. Comput., 17(3):579–612, 1996

  33. [33]

    Huybrechs

    D. Huybrechs. Generalizedgauss.jl. https://github.com/daanhb/ GeneralizedGauss.jl

  34. [34]

    Huybrechs

    D. Huybrechs. On the computation of Gaussian quadrature rules for Chebyshev sets of linearly independent functions.SIAM J. Numer. Anal., 60(3):1168–1192, 2022

  35. [35]

    Huybrechs and R

    D. Huybrechs and R. Cools. On generalized Gaussian quadrature rules for singu- lar and nearly singular integrals.SIAM J. Numer. Anal., 47(1):719–739, 2008/09

  36. [36]

    H. T. Huynh. A flux reconstruction approach to high-order schemes including discontinuous Galerkin methods. In18th AIAA Computational Fluid Dynamics Conference, 2007

  37. [37]

    Karlin and W

    S. Karlin and W. J. Studden.Tchebycheff systems: With applications in analysis and statistics, volume Vol. XV ofPure and Applied Mathematics. Interscience Publishers John Wiley & Sons, New York-London-Sydney, 1966

  38. [38]

    D. A. Kopriva, G. J. Gassner, and J. Nordstr¨ om. Stability of discontinuous Galerkin spectral element schemes for wave propagation when the coefficient matrices have jumps.J. Sci. Comput., 88(1):Paper No. 3, 23, 2021

  39. [39]

    D. A. Kopriva, J. Nordstr¨ om, and G. J. Gassner. Error boundedness of discontin- uous Galerkin spectral element approximations of hyperbolic problems.Journal of Scientific Computing, 72(1):314 – 330, 2017

  40. [40]

    B. G. Korenev.Bessel functions and their applications, volume 8 ofAnalytical Methods and Special Functions. Taylor & Francis Group, London, 2002. Trans- lated from the Russian by E. V. Pankratiev

  41. [41]

    M. G. Krein. The ideas of P. L. ˇCebyˇ sev and A. A. Markov in the theory of limiting values of integrals and their further development.Amer. Math. Soc. Transl. (2), 12:1–121, 1959

  42. [42]

    H.-O. Kreiss. Initial boundary value problems for hyperbolic systems.Comm. Pure Appl. Math., 23:277–298, 1970

  43. [43]

    Kreiss and J

    H.-O. Kreiss and J. Lorenz.Initial-boundary value problems and the Navier- Stokes equations, volume 47 ofClassics in Applied Mathematics. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, 2004. Reprint of the 1989 edition

  44. [44]

    Kreiss and G

    H.-O. Kreiss and G. Scherer. Finite element and finite difference methods for hyperbolic partial differential equations. InAspects of Finite Elements in Partial N. HALE, C. HARLEY, P. NCHUPANG, AND J. NORDSTR ¨OM Differential Equations, Academic Press, Inc., New York (NY), 1974

  45. [45]

    Lundquist and J

    T. Lundquist and J. Nordstr¨ om. The SBP-SAT technique for initial value prob- lems.J. Comput. Phys., 270:86–104, 2014

  46. [46]

    Lundquist and J

    T. Lundquist and J. Nordstr¨ om. Efficient fully discrete summation-by-parts schemes for unsteady flow problems.BIT, 56(3):951–966, 2016

  47. [47]

    Lundquist, A

    T. Lundquist, A. R. Winters, and J. Nordstr¨ om. Encapsulated generalized summation-by-parts formulations for curvilinear and non-conforming meshes.J. Comput. Phys., 498:Paper No. 112699, 16, 2024

  48. [48]

    J. Ma, V. Rokhlin, and S. Wandzura. Generalized Gaussian quadrature rules for systems of arbitrary functions.SIAM J. Numer. Anal., 33(3):971–996, 1996

  49. [49]

    A. Markov. On the limiting values of integrals in connection with interpolation. Zap. Imp. Akad. Nauk. Fiz.-Mat. Otd, 6(8):146–230, 1898

  50. [50]

    Nordstr¨ om

    J. Nordstr¨ om. A roadmap to well posed and stable problems in computational physics.J. Sci. Comput., 71(1):365–385, 2017

  51. [51]

    Nordstr¨ om, S

    J. Nordstr¨ om, S. Eriksson, and P. Eliasson. Weak and strong wall boundary procedures and convergence to steady-state of the Navier-Stokes equations.J. Comput. Phys., 231(14):4867–4884, 2012

  52. [52]

    Nordstr¨ om, K

    J. Nordstr¨ om, K. Forsberg, C. Adamsson, and P. Eliasson. Finite volume meth- ods, unstructured meshes and strict stability for hyperbolic problems.Appl. Numer. Math., 45(4):453–473, 2003

  53. [53]

    Nordstr¨ om and C

    J. Nordstr¨ om and C. La Cognata. Energy stable boundary conditions for the nonlinear incompressible Navier-Stokes equations.Math. Comp., 88(316):665– 690, 2019

  54. [54]

    Nordstr¨ om and T

    J. Nordstr¨ om and T. Lundquist. Summation-by-parts in time.J. Comput. Phys., 251:487–499, 2013

  55. [55]

    Nordstr¨ om and T

    J. Nordstr¨ om and T. Lundquist. Summation-by-parts in time: the second deriv- ative.SIAM J. Sci. Comput., 38(3):A1561–A1586, 2016

  56. [56]

    Nordstr¨ om and A

    J. Nordstr¨ om and A. A. Ruggiu. On conservation and stability properties for summation-by-parts schemes.J. Comput. Phys., 344:451–464, 2017

  57. [57]

    P. Olsson. Summation by parts, projections, and stability. I.Math. Comp., 64(211):1035–1065, S23–S26, 1995

  58. [58]

    P. Olsson. Summation by parts, projections, and stability. II.Math. Comp., 64(212):1473–1493, 1995

  59. [59]

    P. J. Roache. Code verification by the method of manufactured solutions.J. Fluids Eng., 124(1):4–10, 2002

  60. [60]

    B. Strand. Summation by parts for finite difference approximations ford/dx.J. Comput. Phys., 110(1):47–67, 1994

  61. [61]

    Sv¨ ard and J

    M. Sv¨ ard and J. Nordstr¨ om. Review of summation-by-parts schemes for initial- boundary-value problems.J. Comput. Phys., 268:17–38, 2014

  62. [62]

    Yarvin and V

    N. Yarvin and V. Rokhlin. Generalized Gaussian quadratures and singular value decompositions of integral operators.SIAM J. Sci. Comput., 20(2):699–718, 1998