-
Algorithms for constrained optimal transport
Authors:
Martin Corless,
Anthony Quinn,
Sarah Boufelja,
Robert Shorten
Abstract:
We derive iterative scaling algorithms of the Sinkhorn-Knopp (SK) type for constrained optimal transport. The constraints are in the form of prior-imposed zeroes in the transport plan. Based on classical Bregman arguments, we prove asymptotic convergence of our algorithms to a unique optimal solution. New insights obtained from the convergence proof are highlighted. An example from electrical vehi…
▽ More
We derive iterative scaling algorithms of the Sinkhorn-Knopp (SK) type for constrained optimal transport. The constraints are in the form of prior-imposed zeroes in the transport plan. Based on classical Bregman arguments, we prove asymptotic convergence of our algorithms to a unique optimal solution. New insights obtained from the convergence proof are highlighted. An example from electrical vehicle charging in a smart city context is outlined, in which the prior zero-constraints prevent energy from being transported from some providers to some vehicles.
△ Less
Submitted 16 February, 2024;
originally announced April 2024.
-
Teaching Linear Algebra in a Mechanized Mathematical Environment
Authors:
Robert M. Corless,
David J. Jeffrey,
Azar Shakoori
Abstract:
This paper outlines our ideas on how to teach linear algebra in a mechanized mathematical environment, and discusses some of our reasons for thinking that this is a better way to teach linear algebra than the ``old fashioned way''. We discuss some technological tools such as Maple, Matlab, Python, and Jupyter Notebooks, and some choices of topics that are especially suited to teaching with these t…
▽ More
This paper outlines our ideas on how to teach linear algebra in a mechanized mathematical environment, and discusses some of our reasons for thinking that this is a better way to teach linear algebra than the ``old fashioned way''. We discuss some technological tools such as Maple, Matlab, Python, and Jupyter Notebooks, and some choices of topics that are especially suited to teaching with these tools. The discussion is informed by our experience over the past thirty or more years teaching at various levels, especially at the University of Western Ontario.
△ Less
Submitted 31 May, 2023;
originally announced June 2023.
-
Blendstrings: an environment for computing with smooth functions
Authors:
Robert M. Corless
Abstract:
A "blendstring" is a piecewise polynomial interpolant with high-degree two-point Hermite interpolational polynomials on each piece, analogous to a cubic spline. Blendstrings are smoother and can be more accurate than cubic splines, and can be used to represent smooth functions on a line segment or polygonal path in the complex plane. I sketch some properties of blendstrings, including efficient me…
▽ More
A "blendstring" is a piecewise polynomial interpolant with high-degree two-point Hermite interpolational polynomials on each piece, analogous to a cubic spline. Blendstrings are smoother and can be more accurate than cubic splines, and can be used to represent smooth functions on a line segment or polygonal path in the complex plane. I sketch some properties of blendstrings, including efficient methods for evaluation, differentiation, and integration, as well as a prototype Maple implementation. Blendstrings can be differentiated and integrated exactly and can be combined algebraically. I also show applications of blendstrings to solving differential equations and computing Mathieu functions and generalized Mathieu eigenfunctions.
△ Less
Submitted 18 May, 2023;
originally announced May 2023.
-
Elliptic cross sections in blood flow regulation
Authors:
Chris Brimacombe,
Robert M. Corless,
Mair Zamir
Abstract:
Arterial deformations arise in blood flow when surrounding tissue invades the space available for a blood vessel to maintain its circular cross section, the most immediate effects being a reduction in blood flow and redistribution of shear stress. Here we consider deformations from circular to elliptic cross sections. Solution of this problem in steady flow is fairly straightforward. The focus in…
▽ More
Arterial deformations arise in blood flow when surrounding tissue invades the space available for a blood vessel to maintain its circular cross section, the most immediate effects being a reduction in blood flow and redistribution of shear stress. Here we consider deformations from circular to elliptic cross sections. Solution of this problem in steady flow is fairly straightforward. The focus in the present paper is on pulsatile flow where the change from circular to elliptic cross sections is associated with a transition in the character of the equations governing the flow from Bessel to Mathieu equations. The study of this problem has been hampered in the past because of difficulties involved in the solution of the governing equations. In the present study we describe methods we have used to overcome some of these difficulties and present a comprehensive set of results based on these methods. In particular, vessel deformation is examined under two different conditions relevant to blood flow regulation: (i) keeping cross sectional area constant and (ii) keeping cross sectional circumference constant. The results provide an important context for the mechanism of neurovascular control of blood flow under the pathological conditions of vessel deformation.
△ Less
Submitted 19 January, 2023;
originally announced April 2023.
-
Fully Probabilistic Design for Optimal Transport
Authors:
Sarah Boufelja Y.,
Anthony Quinn,
Martin Corless,
Robert Shorten
Abstract:
The goal of this paper is to introduce a new theoretical framework for Optimal Transport (OT), using the terminology and techniques of Fully Probabilistic Design (FPD). Optimal Transport is the canonical method for comparing probability measures and has been successfully applied in a wide range of areas (computer vision Rubner et al. [2004], computer graphics Solomon et al. [2015], natural languag…
▽ More
The goal of this paper is to introduce a new theoretical framework for Optimal Transport (OT), using the terminology and techniques of Fully Probabilistic Design (FPD). Optimal Transport is the canonical method for comparing probability measures and has been successfully applied in a wide range of areas (computer vision Rubner et al. [2004], computer graphics Solomon et al. [2015], natural language processing Kusner et al. [2015], etc.). However, we argue that the current OT framework suffers from two shortcomings: first, it is hard to induce generic constraints and probabilistic knowledge in the OT problem; second, the current formalism does not address the question of uncertainty in the marginals, lacking therefore the mechanisms to design robust solutions. By viewing the OT problem as the optimal design of a probability density function with marginal constraints, we prove that OT is an instance of the more generic FPD framework. In this new setting, we can furnish the OT framework with the necessary mechanisms for processing probabilistic constraints and deriving uncertainty quantifiers, hence establishing a new extended framework, called FPD-OT. Our main contribution in this paper is to establish the connection between OT and FPD, providing new theoretical insights for both. This will lay the foundations for the application of FPD-OT in a subsequent work, notably in processing more sophisticated knowledge constraints, as well as in designing robust solutions in the case of uncertain marginals.
△ Less
Submitted 19 December, 2022;
originally announced December 2022.
-
Perturbation Theory
Authors:
Nicolas Fillion,
Robert M. Corless
Abstract:
This article aims to explain essential elements of perturbation theory and their conceptual underpinnings. It is not meant as a summary of popular perturbation methods, though some illustrative examples are given to underline the main methodological insights and concerns. We also give brief explications of the mathematical notions of limit, continuity, differentiability, convergence, and divergenc…
▽ More
This article aims to explain essential elements of perturbation theory and their conceptual underpinnings. It is not meant as a summary of popular perturbation methods, though some illustrative examples are given to underline the main methodological insights and concerns. We also give brief explications of the mathematical notions of limit, continuity, differentiability, convergence, and divergence, which provide the necessary foundation.
△ Less
Submitted 14 December, 2022;
originally announced December 2022.
-
Bohemian Matrix Geometry
Authors:
Robert M. Corless,
George Labahn,
Dan Piponi,
Leili Rafiee Sevyeri
Abstract:
A Bohemian matrix family is a set of matrices all of whose entries are drawn from a fixed, usually discrete and hence bounded, subset of a field of characteristic zero. Originally these were integers -- hence the name, from the acronym BOunded HEight Matrix of Integers (BOHEMI) -- but other kinds of entries are also interesting. Some kinds of questions about Bohemian matrices can be answered by nu…
▽ More
A Bohemian matrix family is a set of matrices all of whose entries are drawn from a fixed, usually discrete and hence bounded, subset of a field of characteristic zero. Originally these were integers -- hence the name, from the acronym BOunded HEight Matrix of Integers (BOHEMI) -- but other kinds of entries are also interesting. Some kinds of questions about Bohemian matrices can be answered by numerical computation, but sometimes exact computation is better. In this paper we explore some Bohemian families (symmetric, upper Hessenberg, or Toeplitz) computationally, and answer some open questions posed about the distributions of eigenvalue densities.
△ Less
Submitted 26 April, 2022; v1 submitted 15 February, 2022;
originally announced February 2022.
-
Computational Discovery with Newton Fractals, Bohemian Matrices, & Mandelbrot Polynomials
Authors:
Neil J. Calkin,
Eunice Y. S. Chan,
Robert M. Corless
Abstract:
The authors have been using a largely algebraic form of ``computational discovery'' in various undergraduate classes at their respective institutions for some decades now to teach pure mathematics, applied mathematics, and computational mathematics. This paper describes what we mean by ``computational discovery,'' what good it does for the students, and some specific techniques that we used.
The authors have been using a largely algebraic form of ``computational discovery'' in various undergraduate classes at their respective institutions for some decades now to teach pure mathematics, applied mathematics, and computational mathematics. This paper describes what we mean by ``computational discovery,'' what good it does for the students, and some specific techniques that we used.
△ Less
Submitted 8 September, 2021;
originally announced September 2021.
-
A Fractal Eigenvector
Authors:
Neil J. Calkin,
Eunice Y. S. Chan,
Robert M. Corless,
David J. Jeffrey,
Piers W. Lawrence
Abstract:
The recursively-constructed family of Mandelbrot matrices $M_n$ for $n=1$, $2$, $\ldots$ have nonnegative entries (indeed just $0$ and $1$, so each $M_n$ can be called a binary matrix) and have eigenvalues whose negatives $-λ= c$ give periodic orbits under the Mandelbrot iteration, namely $z_k = z_{k-1}^2+c$ with $z_0=0$, and are thus contained in the Mandelbrot set. By the Perron--Frobenius theor…
▽ More
The recursively-constructed family of Mandelbrot matrices $M_n$ for $n=1$, $2$, $\ldots$ have nonnegative entries (indeed just $0$ and $1$, so each $M_n$ can be called a binary matrix) and have eigenvalues whose negatives $-λ= c$ give periodic orbits under the Mandelbrot iteration, namely $z_k = z_{k-1}^2+c$ with $z_0=0$, and are thus contained in the Mandelbrot set. By the Perron--Frobenius theorem, the matrices $M_n$ have a dominant real positive eigenvalue, which we call $ρ_n$. This article examines the eigenvector belonging to that dominant eigenvalue and its fractal-like structure, and similarly examines (with less success) the dominant singular vectors of $M_n$ from the singular value decomposition.
△ Less
Submitted 29 March, 2021;
originally announced April 2021.
-
Approximate GCD in Lagrange bases
Authors:
Leili Rafiee Sevyeri,
Robert M. Corless
Abstract:
For a pair of polynomials with real or complex coefficients, given in any particular basis, the problem of finding their GCD is known to be ill-posed. An answer is still desired for many applications, however. Hence, looking for a GCD of so-called approximate polynomials where this term explicitly denotes small uncertainties in the coefficients has received significant attention in the field of hy…
▽ More
For a pair of polynomials with real or complex coefficients, given in any particular basis, the problem of finding their GCD is known to be ill-posed. An answer is still desired for many applications, however. Hence, looking for a GCD of so-called approximate polynomials where this term explicitly denotes small uncertainties in the coefficients has received significant attention in the field of hybrid symbolic-numeric computation. In this paper we give an algorithm, based on one of Victor Ya. Pan, to find an approximate GCD for a pair of approximate polynomials given in a Lagrange basis. More precisely, we suppose that these polynomials are given by their approximate values at distinct known points. We first find each of their roots by using a Lagrange basis companion matrix for each polynomial, cluster the roots of each polynomial to identify multiple roots, and then "marry" the two polynomials to find their GCD. At no point do we change to the monomial basis, thus preserving the good conditioning properties of the original Lagrange basis. We discuss advantages and drawbacks of this method. The computational cost is dominated by the rootfinding step; unless special-purpose eigenvalue algorithms are used, the cost is cubic in the degrees of the polynomials. In principle, this cost could be reduced but we do not do so here.
△ Less
Submitted 25 March, 2021;
originally announced March 2021.
-
Memoir on the vibratory movement of an elliptical membrane
Authors:
Robert H. C. Moir,
Robert M. Corless
Abstract:
Translation of Émile Mathieu's "Mémoire sur le mouvement vibratoire d'une membrane de forme elliptique", Journal de mathématiques pures et appliquées 2e série, tome 13 (1868), p. 137-203.
Translation of Émile Mathieu's "Mémoire sur le mouvement vibratoire d'une membrane de forme elliptique", Journal de mathématiques pures et appliquées 2e série, tome 13 (1868), p. 137-203.
△ Less
Submitted 3 March, 2021;
originally announced March 2021.
-
Equivalences for Linearizations of Matrix Polynomials
Authors:
Robert M. Corless,
Leili Rafiee Sevyeri,
B. David Saunders
Abstract:
One useful standard method to compute eigenvalues of matrix polynomials ${\bf P}(z) \in \mathbb{C}^{n\times n}[z]$ of degree at most $\ell$ in $z$ (denoted of grade $\ell$, for short) is to first transform ${\bf P}(z)$ to an equivalent linear matrix polynomial ${\bf L}(z)=z{\bf B}-{\bf A}$, called a companion pencil, where ${\bf A}$ and ${\bf B}$ are usually of larger dimension than ${\bf P}(z)$ b…
▽ More
One useful standard method to compute eigenvalues of matrix polynomials ${\bf P}(z) \in \mathbb{C}^{n\times n}[z]$ of degree at most $\ell$ in $z$ (denoted of grade $\ell$, for short) is to first transform ${\bf P}(z)$ to an equivalent linear matrix polynomial ${\bf L}(z)=z{\bf B}-{\bf A}$, called a companion pencil, where ${\bf A}$ and ${\bf B}$ are usually of larger dimension than ${\bf P}(z)$ but ${\bf L}(z)$ is now only of grade $1$ in $z$. The eigenvalues and eigenvectors of ${\bf L}(z)$ can be computed numerically by, for instance, the QZ algorithm. The eigenvectors of ${\bf P}(z)$, including those for infinite eigenvalues, can also be recovered from eigenvectors of ${\bf L}(z)$ if ${\bf L}(z)$ is what is called a "strong linearization" of ${\bf P}(z)$.
In this paper we show how to use algorithms for computing the Hermite Normal Form of a companion matrix for a scalar polynomial to direct the discovery of unimodular matrix polynomial cofactors ${\bf E}(z)$ and ${\bf F}(z)$ which, via the equation ${\bf E}(z){\bf L}(z){\bf F}(z) = \mathrm{diag}( {\bf P}(z), {\bf I}_n, \ldots, {\bf I}_n)$, explicitly show the equivalence of ${\bf P}(z)$ and ${\bf L}(z)$. By this method we give new explicit constructions for several linearizations using different polynomial bases. We contrast these new unimodular pairs with those constructed by strict equivalence, some of which are also new to this paper. We discuss the limitations of this experimental, computational discovery method of finding unimodular cofactors.
△ Less
Submitted 18 February, 2021;
originally announced February 2021.
-
A Note on Order and Index Reduction for Descriptor Systems
Authors:
Martin Corless,
Robert Shorten
Abstract:
We present order reduction results for linear time invariant descriptor systems. Results are given for both forced and unforced systems as well methods for constructing the reduced order systems. Our results establish a precise connection between classical and new results on this topic, and lead to an elementary construction of quasi-Weierstrass forms for a descriptor system. Examples are given to…
▽ More
We present order reduction results for linear time invariant descriptor systems. Results are given for both forced and unforced systems as well methods for constructing the reduced order systems. Our results establish a precise connection between classical and new results on this topic, and lead to an elementary construction of quasi-Weierstrass forms for a descriptor system. Examples are given to illustrate the usefulness of our results.
△ Less
Submitted 18 January, 2021;
originally announced January 2021.
-
Chaos Game Representation
Authors:
Eunice Y. S. Chan,
Robert M. Corless
Abstract:
The chaos game representation (CGR) is an interesting method to visualize one-dimensional sequences. In this paper, we show how to construct a chaos game representation. The applications mentioned here are biological, in which CGR was able to uncover patterns in DNA or proteins that were previously unknown. We also show how CGR might be introduced in the classroom, either in a modelling course or…
▽ More
The chaos game representation (CGR) is an interesting method to visualize one-dimensional sequences. In this paper, we show how to construct a chaos game representation. The applications mentioned here are biological, in which CGR was able to uncover patterns in DNA or proteins that were previously unknown. We also show how CGR might be introduced in the classroom, either in a modelling course or in a dynamical systems course. Some sequences that are tested are taken from the Online Encyclopedia of Integer Sequences, and others are taken from sequences that arose mainly from a course in experimental mathematics.
△ Less
Submitted 15 December, 2020;
originally announced December 2020.
-
Teaching Programming for Mathematical Scientists
Authors:
Jack Betteridge,
Eunice Y. S. Chan,
Robert M. Corless,
James H. Davenport,
James Grant
Abstract:
Over the past thirty years or so the authors have been teaching various programming for mathematics courses at our respective Universities, as well as incorporating computer algebra and numerical computation into traditional mathematics courses. These activities are, in some important ways, natural precursors to the use of Artificial Intelligence in Mathematics Education. This paper reflects on so…
▽ More
Over the past thirty years or so the authors have been teaching various programming for mathematics courses at our respective Universities, as well as incorporating computer algebra and numerical computation into traditional mathematics courses. These activities are, in some important ways, natural precursors to the use of Artificial Intelligence in Mathematics Education. This paper reflects on some of our course designs and experiences and is therefore a mix of theory and practice. Underlying both is a clear recognition of the value of computer programming for mathematics education. We use this theory and practice to suggest good techniques for and to raise questions about the use of AI in Mathematics Education.
△ Less
Submitted 30 October, 2020;
originally announced October 2020.
-
Integrals of functions containing parameters
Authors:
Robert M. Corless,
David J. Jeffrey,
David R. Stoutemyer
Abstract:
This paper offers what seems at first to be a minor technical correction to the current practice of computing indefinite integrals, and introduces the idea of a "Kahanian constant of integration". However, the total impact of this minor correction is potentially large because the current practice is taught early at the university level and to very many students---most of whom do not go on to becom…
▽ More
This paper offers what seems at first to be a minor technical correction to the current practice of computing indefinite integrals, and introduces the idea of a "Kahanian constant of integration". However, the total impact of this minor correction is potentially large because the current practice is taught early at the university level and to very many students---most of whom do not go on to become mathematics majors. Moreover, computer algebra systems have become widespread, including good free ones, some of which are available for smartphones. Most current computer algebra systems apply current textbook rules and amplify the effects of fundamental "minor" errors such as the error in continuity that we address in this article. So in practice, the correction we present is important.
△ Less
Submitted 17 September, 2020;
originally announced September 2020.
-
Pure tone modes for a 5:3 elliptic drum
Authors:
Robert M. Corless
Abstract:
The paper exhibits several standing modes of a 5:3 elliptic drum computed using Mathieu functions. To match the boundary conditions, I used Newton's method on the appropriate modified Mathieu equation using the Squire-Trapp formula for computing derivatives. I tabulate the requisite values of the parameter $q$ for these low-frequency modes.
The paper exhibits several standing modes of a 5:3 elliptic drum computed using Mathieu functions. To match the boundary conditions, I used Newton's method on the appropriate modified Mathieu equation using the Squire-Trapp formula for computing derivatives. I tabulate the requisite values of the parameter $q$ for these low-frequency modes.
△ Less
Submitted 16 August, 2020;
originally announced August 2020.
-
Computation and applications of Mathieu functions: A historical perspective
Authors:
Chris Brimacombe,
Robert M. Corless,
Mair Zamir
Abstract:
Mathieu functions of period $π$ or $2π$, also called elliptic cylinder functions, were introduced in 1868 by Émile Mathieu together with so-called modified Mathieu functions, in order to help understand the vibrations of an elastic membrane set in a fixed elliptical hoop. These functions still occur frequently in applications today: our interest, for instance, was stimulated by a problem of pulsat…
▽ More
Mathieu functions of period $π$ or $2π$, also called elliptic cylinder functions, were introduced in 1868 by Émile Mathieu together with so-called modified Mathieu functions, in order to help understand the vibrations of an elastic membrane set in a fixed elliptical hoop. These functions still occur frequently in applications today: our interest, for instance, was stimulated by a problem of pulsatile blood flow in a blood vessel compressed into an elliptical cross-section. This paper surveys and recapitulates the historical development of the theory and methods of computation for Mathieu functions and modified Mathieu functions and identifies some gaps in current software capability, particularly to do with double eigenvalues of the Mathieu equation. We demonstrate how to compute Puiseux expansions of the Mathieu eigenvalues about such double eigenvalues, and give methods to compute the generalized eigenfunctions that arise there. In examining Mathieu's original contribution, we bring out that his use of anti-secularity predates that of Lindstedt. For interest, we also provide short biographies of some of the major mathematical researchers involved in the history of the Mathieu functions: Émile Mathieu, Sir Edmund Whittaker, Edward Ince, and Gertrude Blanch.
△ Less
Submitted 30 June, 2021; v1 submitted 4 August, 2020;
originally announced August 2020.
-
Inverse Cubic Iteration
Authors:
Robert M. Corless
Abstract:
There are thousands of papers on rootfinding for nonlinear scalar equations. Here is one more, to
talk about an apparently new method, which I call ``Inverse Cubic Iteration'' (ICI) in analogy to the Inverse Quadratic Iteration in Richard Brent's zeroin method. The possibly new method is based on a cubic blend of tangent-line approximations for the inverse function. We rewrite this iteration for…
▽ More
There are thousands of papers on rootfinding for nonlinear scalar equations. Here is one more, to
talk about an apparently new method, which I call ``Inverse Cubic Iteration'' (ICI) in analogy to the Inverse Quadratic Iteration in Richard Brent's zeroin method. The possibly new method is based on a cubic blend of tangent-line approximations for the inverse function. We rewrite this iteration for numerical stability as an average of two Newton steps and a secant step: only one new function evaluation and derivative evaluation is needed for each step. The total cost of the method is therefore only trivially more than Newton's method, and we will see that it has order $1+\sqrt{3} = 2.732...$, thus ensuring that to achieve a given accuracy it usually takes fewer steps than Newton's method while using essentially the same effort per step.
△ Less
Submitted 13 July, 2020;
originally announced July 2020.
-
Blends in Maple
Authors:
Robert M. Corless,
Erik Postma
Abstract:
A blend of two Taylor series for the same smooth real- or complex-valued function of a single variable can be useful for approximation. We use an explicit formula for a two-point Hermite interpolational polynomial to construct such blends. We show a robust Maple implementation that can stably and efficiently evaluate blends using linear-cost Horner form, evaluate their derivatives to arbitrary ord…
▽ More
A blend of two Taylor series for the same smooth real- or complex-valued function of a single variable can be useful for approximation. We use an explicit formula for a two-point Hermite interpolational polynomial to construct such blends. We show a robust Maple implementation that can stably and efficiently evaluate blends using linear-cost Horner form, evaluate their derivatives to arbitrary order at the same time, or integrate a blend exactly. The implementation is suited for use with evalhf. We provide a top-level user interface and efficient module exports for programmatic use. This work was presented at the Maple Conference 2020. See www.maplesoft.com/mapleconference
△ Less
Submitted 27 November, 2020; v1 submitted 9 July, 2020;
originally announced July 2020.
-
Approximate GCD in a Bernstein basis
Authors:
Robert M. Corless,
Leili Rafiee Sevyeri
Abstract:
We adapt Victor Y. Pan's root-based algorithm for finding approximate GCD to the case where the polynomials are expressed in Bernstein bases. We use the numerically stable companion pencil of Gudbjörn F. Jónsson to compute the roots, and the Hopcroft-Karp bipartite matching method to find the degree of the approximate GCD. We offer some refinements to improve the process.
We adapt Victor Y. Pan's root-based algorithm for finding approximate GCD to the case where the polynomials are expressed in Bernstein bases. We use the numerically stable companion pencil of Gudbjörn F. Jónsson to compute the roots, and the Hopcroft-Karp bipartite matching method to find the degree of the approximate GCD. We offer some refinements to improve the process.
△ Less
Submitted 4 October, 2019;
originally announced October 2019.
-
Upper Hessenberg and Toeplitz Bohemians
Authors:
Eunice Y. S. Chan,
Robert M. Corless,
Laureano Gonzalez-Vega,
J. Rafael Sendra,
Juana Sendra
Abstract:
We look at Bohemians, specifically those with population $\{-1, 0, {+1}\}$ and sometimes $\{0,1,i,-1,-i\}$. More, we specialize the matrices to be upper Hessenberg Bohemian. From there, focusing on only those matrices whose characteristic polynomials have maximal height allows us to explicitly identify these polynomials and give useful bounds on their height, and conjecture an accurate asymptotic…
▽ More
We look at Bohemians, specifically those with population $\{-1, 0, {+1}\}$ and sometimes $\{0,1,i,-1,-i\}$. More, we specialize the matrices to be upper Hessenberg Bohemian. From there, focusing on only those matrices whose characteristic polynomials have maximal height allows us to explicitly identify these polynomials and give useful bounds on their height, and conjecture an accurate asymptotic formula. The lower bound for the maximal characteristic height is exponential in the order of the matrix; in contrast, the height of the matrices remains constant. We give theorems about the numbers of normal matrices and the numbers of stable matrices in these families.
△ Less
Submitted 23 July, 2019;
originally announced July 2019.
-
L2 Observers for a Class of Nonlinear Systems with Unknown Inputs
Authors:
Martin Corless,
Ankush Chakrabarty
Abstract:
We consider the problem of estimating the state and unknown input for a large class of nonlinear systems subject to unknown exogenous inputs. The exogenous inputs themselves are modeled as being generated by a nonlinear system subject to unknown inputs. The nonlinearities considered in this work are characterized by multiplier matrices that include many commonly encountered nonlinearities. We obta…
▽ More
We consider the problem of estimating the state and unknown input for a large class of nonlinear systems subject to unknown exogenous inputs. The exogenous inputs themselves are modeled as being generated by a nonlinear system subject to unknown inputs. The nonlinearities considered in this work are characterized by multiplier matrices that include many commonly encountered nonlinearities. We obtain a linear matrix inequality (LMI), that, if feasible, provides the gains for an observer which results in certified L2 performance of the error dynamics associated with the observer. We also present conditions which guarantee that the L2 norm of the error can be made arbitrarily small and investigate conditions for feasibility of the proposed LMIs.
△ Less
Submitted 21 February, 2019;
originally announced February 2019.
-
Bohemian Upper Hessenberg Toeplitz Matrices
Authors:
Eunice Y. S. Chan,
Robert M. Corless,
Laureano Gonzalez-Vega,
J. Rafael Sendra,
Juana Sendra,
Steven E. Thornton
Abstract:
We look at Bohemian matrices, specifically those with entries from $\{-1, 0, {+1}\}$. More, we specialize the matrices to be upper Hessenberg, with subdiagonal entries $1$. Even more, we consider Toeplitz matrices of this kind. Many properties remain after these specializations, some of which surprised us. Focusing on only those matrices whose characteristic polynomials have maximal height allows…
▽ More
We look at Bohemian matrices, specifically those with entries from $\{-1, 0, {+1}\}$. More, we specialize the matrices to be upper Hessenberg, with subdiagonal entries $1$. Even more, we consider Toeplitz matrices of this kind. Many properties remain after these specializations, some of which surprised us. Focusing on only those matrices whose characteristic polynomials have maximal height allows us to explicitly identify these polynomials and give a lower bound on their height. This bound is exponential in the order of the matrix.
△ Less
Submitted 27 September, 2018;
originally announced September 2018.
-
Bohemian Upper Hessenberg Matrices
Authors:
Eunice Y. S. Chan,
Robert M. Corless,
Laureano Gonzalez-Vega,
J. Rafael Sendra,
Juana Sendra,
Steven E. Thornton
Abstract:
We look at Bohemian matrices, specifically those with entries from $\{-1, 0, {+1}\}$. More, we specialize the matrices to be upper Hessenberg, with subdiagonal entries $\pm1$. Many properties remain after these specializations, some of which surprised us. We find two recursive formulae for the characteristic polynomials of upper Hessenberg matrices. Focusing on only those matrices whose characteri…
▽ More
We look at Bohemian matrices, specifically those with entries from $\{-1, 0, {+1}\}$. More, we specialize the matrices to be upper Hessenberg, with subdiagonal entries $\pm1$. Many properties remain after these specializations, some of which surprised us. We find two recursive formulae for the characteristic polynomials of upper Hessenberg matrices. Focusing on only those matrices whose characteristic polynomials have maximal height allows us to explicitly identify these polynomials and give a lower bound on their height. This bound is exponential in the order of the matrix. We count stable matrices, normal matrices, and neutral matrices, and tabulate the results of our experiments. We prove a theorem about the only possible kinds of normal matrices amongst a specific family of Bohemian upper Hessenberg matrices.
△ Less
Submitted 27 September, 2018;
originally announced September 2018.
-
Differentiation Matrices for Univariate Polynomials
Authors:
Amirhossein Amiraslani,
Robert M. Corless,
Madhusoodan Gunasingham
Abstract:
We collect here elementary properties of differentiation matrices for univariate polynomials expressed in various bases, including orthogonal polynomial bases and non-degree-graded bases such as Bernstein bases and Lagrange \& Hermite interpolational bases.
We collect here elementary properties of differentiation matrices for univariate polynomials expressed in various bases, including orthogonal polynomial bases and non-degree-graded bases such as Bernstein bases and Lagrange \& Hermite interpolational bases.
△ Less
Submitted 15 September, 2018;
originally announced September 2018.
-
Revisiting Gilbert Strang's "A Chaotic Search for $i$"
Authors:
Ao Li,
Robert M. Corless
Abstract:
In the paper "A Chaotic Search for $i$"~(\cite{strang1991chaotic}), Strang completely explained the behaviour of Newton's method when using real initial guesses on $f(x) = x^{2}+1$, which has only a pair of complex roots $\pm i$. He explored an exact symbolic formula for the iteration, namely $x_{n}=\cot{ \left( 2^{n} θ_{0} \right) }$, which is valid in exact arithmetic. In this paper, we extend t…
▽ More
In the paper "A Chaotic Search for $i$"~(\cite{strang1991chaotic}), Strang completely explained the behaviour of Newton's method when using real initial guesses on $f(x) = x^{2}+1$, which has only a pair of complex roots $\pm i$. He explored an exact symbolic formula for the iteration, namely $x_{n}=\cot{ \left( 2^{n} θ_{0} \right) }$, which is valid in exact arithmetic. In this paper, we extend this to to $k^{th}$ order Householder methods, which include Halley's method, and to the secant method. Two formulae, $x_{n}=\cot{ \left( θ_{n-1}+θ_{n-2} \right) }$ with $θ_{n-1}=\mathrm{arccot}{\left(x_{n-1}\right)}$ and $θ_{n-2}=\mathrm{arccot}{\left(x_{n-2}\right)}$, and $x_{n}=\cot{ \left( (k+1)^{n} θ_{0} \right) }$ with $θ_{0} = \mathrm{arccot}(x_{0})$, are provided. The asymptotic behaviour and periodic character are illustrated by experimental computation. We show that other methods (Schröder iterations of the first kind) are generally not so simple. We also explain an old method that can be used to allow Maple's \textsl{Fractals[Newton]} package to visualize general one-step iterations by disguising them as Newton iterations.
△ Less
Submitted 8 August, 2018;
originally announced August 2018.
-
Algebraic Linearizations of Matrix Polynomials
Authors:
Eunice Y. S. Chan,
Robert M. Corless,
Laureano Gonzalez-Vega,
J. Rafael Sendra,
Juana Sendra
Abstract:
We show how to construct linearizations of matrix polynomials $z\mathbf{a}(z)\mathbf{d}_0 + \mathbf{c}_0$, $\mathbf{a}(z)\mathbf{b}(z)$, $\mathbf{a}(z) + \mathbf{b}(z)$ (when $\mathrm{deg}\left(\mathbf{b}(z)\right) < \mathrm{deg}\left(\mathbf{a}(z)\right)$), and $z\mathbf{a}(z)\mathbf{d}_0\mathbf{b}(z) + \mathbf{c_0}$ from linearizations of the component parts, $\mathbf{a}(z)$ and $\mathbf{b}(z)$.…
▽ More
We show how to construct linearizations of matrix polynomials $z\mathbf{a}(z)\mathbf{d}_0 + \mathbf{c}_0$, $\mathbf{a}(z)\mathbf{b}(z)$, $\mathbf{a}(z) + \mathbf{b}(z)$ (when $\mathrm{deg}\left(\mathbf{b}(z)\right) < \mathrm{deg}\left(\mathbf{a}(z)\right)$), and $z\mathbf{a}(z)\mathbf{d}_0\mathbf{b}(z) + \mathbf{c_0}$ from linearizations of the component parts, $\mathbf{a}(z)$ and $\mathbf{b}(z)$. This allows the extension to matrix polynomials of a new companion matrix construction.
△ Less
Submitted 29 May, 2018;
originally announced May 2018.
-
Optimal Residuals and the Dahlquist Test Problem
Authors:
Robert M. Corless,
C. Yalcin Kaya,
Robert H. C. Moir
Abstract:
We show how to compute the \emph{optimal relative backward error} for the numerical solution of the Dahlquist test problem by one-step methods. This is an example of a general approach that uses results from optimal control theory to compute optimal residuals, but elementary methods can also be used here because the problem is so simple. This analysis produces new insight into the numerical soluti…
▽ More
We show how to compute the \emph{optimal relative backward error} for the numerical solution of the Dahlquist test problem by one-step methods. This is an example of a general approach that uses results from optimal control theory to compute optimal residuals, but elementary methods can also be used here because the problem is so simple. This analysis produces new insight into the numerical solution of stiff problems.
△ Less
Submitted 28 May, 2018;
originally announced May 2018.
-
Compact Finite Differences and Cubic Splines
Authors:
Robert M. Corless
Abstract:
In this paper I uncover and explain---using contour integrals and residues---a connection between cubic splines and a popular compact finite difference formula. The connection is that on a uniform mesh the simplest Padé scheme for generating fourth-order accurate compact finite differences gives \textsl{exactly} the derivatives at the interior nodes needed to guarantee twice-continuous differentia…
▽ More
In this paper I uncover and explain---using contour integrals and residues---a connection between cubic splines and a popular compact finite difference formula. The connection is that on a uniform mesh the simplest Padé scheme for generating fourth-order accurate compact finite differences gives \textsl{exactly} the derivatives at the interior nodes needed to guarantee twice-continuous differentiability for cubic splines. %I found this connection surprising, because the two problems being solved are different. I also introduce an apparently new spline-like interpolant that I call a compact cubic interpolant; this is similar to one introduced in 1972 by Swartz and Varga, but has higher order accuracy at the edges. I argue that for mildly nonuniform meshes the compact cubic approach offers some potential advantages, and even for uniform meshes offers a simple way to treat the edge conditions, relieving the user of the burden of deciding to use one of the three standard options: free (natural), complete (clamped), or "not-a-knot" conditions. Finally, I establish that the matrices defining the compact cubic splines (equivalently, the fourth-order compact finite difference formulæ) are positive definite, and in fact totally nonnegative, if all mesh widths are the same sign.
△ Less
Submitted 22 November, 2019; v1 submitted 19 May, 2018;
originally announced May 2018.
-
Generalized Standard Triples for Algebraic Linearizations of Matrix Polynomials
Authors:
Eunice Y. S. Chan,
Robert M. Corless,
Leili Rafiee Sevyeri
Abstract:
We define \emph{generalized standard triples} $\mathbf{X}$, $\mathbf{Y}$, and $L(z) = z\mathbf{C}_{1} - \mathbf{C}_{0}$, where $L(z)$ is a linearization of a regular matrix polynomial $\mathbf{P}(z) \in \mathbb{C}^{n \times n}[z]$, in order to use the representation $\mathbf{X}(z \mathbf{C}_{1}~-~\mathbf{C}_{0})^{-1}\mathbf{Y}~=~\mathbf{P}^{-1}(z)$ which holds except when $z$ is an eigenvalue of…
▽ More
We define \emph{generalized standard triples} $\mathbf{X}$, $\mathbf{Y}$, and $L(z) = z\mathbf{C}_{1} - \mathbf{C}_{0}$, where $L(z)$ is a linearization of a regular matrix polynomial $\mathbf{P}(z) \in \mathbb{C}^{n \times n}[z]$, in order to use the representation $\mathbf{X}(z \mathbf{C}_{1}~-~\mathbf{C}_{0})^{-1}\mathbf{Y}~=~\mathbf{P}^{-1}(z)$ which holds except when $z$ is an eigenvalue of $\mathbf{P}$. This representation can be used in constructing so-called \emph{algebraic linearizations} for matrix polynomials of the form $\mathbf{H}(z) = z \mathbf{A}(z)\mathbf{B}(z) + \mathbf{C} \in \mathbb{C}^{n \times n}[z]$ from generalized standard triples of $\mathbf{A}(z)$ and $\mathbf{B}(z)$. This can be done even if $\mathbf{A}(z)$ and $\mathbf{B}(z)$ are expressed in differing polynomial bases. Our main theorem is that $\mathbf{X}$ can be expressed using the coefficients of the expression $1 = \sum_{k=0}^\ell e_k φ_k(z)$ in terms of the relevant polynomial basis. For convenience, we tabulate generalized standard triples for orthogonal polynomial bases, the monomial basis, and Newton interpolational bases; for the Bernstein basis; for Lagrange interpolational bases; and for Hermite interpolational bases. We account for the possibility of common similarity transformations.
△ Less
Submitted 25 March, 2021; v1 submitted 11 May, 2018;
originally announced May 2018.
-
Optimal Solution of Linear Ordinary Differential Equations by Conjugate Gradient Method
Authors:
Wenqiang Yang,
Wenyuan Wu,
Robert M. Corless
Abstract:
Solving initial value problems and boundary value problems of Linear Ordinary Differential Equations (ODEs) plays an important role in many applications. There are various numerical methods and solvers to obtain approximate solutions represented by points. However, few work about optimal solution to minimize the residual can be found in the literatures. In this paper, we first use Hermit cubic spl…
▽ More
Solving initial value problems and boundary value problems of Linear Ordinary Differential Equations (ODEs) plays an important role in many applications. There are various numerical methods and solvers to obtain approximate solutions represented by points. However, few work about optimal solution to minimize the residual can be found in the literatures. In this paper, we first use Hermit cubic spline interpolation at mesh points to represent the solution, then we define the residual error as the square of the L2 norm of the residual obtained by substituting the interpolation solution back to ODEs. Thus, solving ODEs is reduced to an optimization problem in curtain solution space which can be solved by conjugate gradient method with taking advantages of sparsity of the corresponding matrix. The examples of IVP and BVP in the paper show that this method can find a solution with smaller global error without additional mesh points.
△ Less
Submitted 20 May, 2018; v1 submitted 2 May, 2018;
originally announced May 2018.
-
The Runge Example for Interpolation and Wilkinson's Examples for Rootfinding
Authors:
Robert M. Corless,
Leili Rafiee Sevyeri
Abstract:
We look at two classical examples in the theory of numerical analysis, namely the Runge example for interpolation and Wilkinson's example (actually two examples) for rootfinding. We use the modern theory of backward error analysis and conditioning, as instigated and popularized by Wilkinson, but refined by Farouki and Rajan. By this means, we arrive at a satisfactory explanation of the puzzling ph…
▽ More
We look at two classical examples in the theory of numerical analysis, namely the Runge example for interpolation and Wilkinson's example (actually two examples) for rootfinding. We use the modern theory of backward error analysis and conditioning, as instigated and popularized by Wilkinson, but refined by Farouki and Rajan. By this means, we arrive at a satisfactory explanation of the puzzling phenomena encountered by students when they try to fit polynomials to numerical data, or when they try to use numerical rootfinding to find polynomial zeros. Computer algebra, with its controlled, arbitrary precision, plays an important didactic role.
△ Less
Submitted 23 April, 2018;
originally announced April 2018.
-
Stirling's Original Asymptotic Series from a Formula like one of Binet's and its Evaluation by Sequence Acceleration
Authors:
Robert M. Corless,
Leili Rafiee Sevyeri
Abstract:
We give an apparently new proof of Stirling's original asymptotic formula for the behavior of $\ln z!$ for large $z$. Stirling's original formula is not the formula widely known as "Stirling's formula", which was actually due to De Moivre. We also show by experiment that this old formula is quite effective for numerical evaluation of $\ln z!$ over $\mathbb{C}$, when coupled with the sequence accel…
▽ More
We give an apparently new proof of Stirling's original asymptotic formula for the behavior of $\ln z!$ for large $z$. Stirling's original formula is not the formula widely known as "Stirling's formula", which was actually due to De Moivre. We also show by experiment that this old formula is quite effective for numerical evaluation of $\ln z!$ over $\mathbb{C}$, when coupled with the sequence acceleration method known as Levin's $u$-transform. As an homage to Stirling, who apparently used inverse symbolic computation to identify the constant term in his formula, we do the same in our proof.
△ Less
Submitted 6 May, 2019; v1 submitted 14 April, 2018;
originally announced April 2018.
-
Observer-Based Controllers for Incrementally Quadratic Nonlinear Systems with Disturbances
Authors:
Xiangru Xu,
Behcet Acikmese,
Martin J. Corless
Abstract:
Robust global stabilization of nonlinear systems by observer-based feedback controllers is a challenging task. This article investigates the problem of designing observer-based stabilizing controllers for incrementally quadratic nonlinear systems with external disturbances. The nonlinearities considered in the system model satisfy the incremental quadratic constraints, which are characterized by i…
▽ More
Robust global stabilization of nonlinear systems by observer-based feedback controllers is a challenging task. This article investigates the problem of designing observer-based stabilizing controllers for incrementally quadratic nonlinear systems with external disturbances. The nonlinearities considered in the system model satisfy the incremental quadratic constraints, which are characterized by incremental multiplier matrices and encompass many common nonlinearities. The simultaneous search for the observer and the controller gain matrices is formulated as a feasibility problem of linear matrix inequalities, for two parameterizations (i.e., the block diagonal parameterization and the block anti-triangular parameterization) of the incremental multiplier matrices, respectively. The closed-loop system implementing the observer-based feedback controller is proven to be input-to-state stable with respect to external disturbances. Using the proposed continuous-time observer-based controllers, event-triggered controllers with time regularization are constructed for globally Lipschitz systems, such that the closed-loop system is Zeno-free and input-to-state practically stable.
△ Less
Submitted 3 June, 2020; v1 submitted 24 February, 2018;
originally announced February 2018.
-
A random walk through experimental mathematics
Authors:
Eunice Y. S. Chan,
Robert M. Corless
Abstract:
We describe our adventures in creating a new first-year course in Experimental Mathematics that uses active learning. We used a state-of-the-art facility, called The Western Active Learning Space, and got the students to "drive the spaceship" (at least a little bit). This paper describes some of our techniques for pedagogy, some of the vignettes of experimental mathematics that we used, and some o…
▽ More
We describe our adventures in creating a new first-year course in Experimental Mathematics that uses active learning. We used a state-of-the-art facility, called The Western Active Learning Space, and got the students to "drive the spaceship" (at least a little bit). This paper describes some of our techniques for pedagogy, some of the vignettes of experimental mathematics that we used, and some of the outcomes. EYSC was a student in the simultaneously-taught senior sister course "Open Problems in Experimental Mathematics" the first time it was taught and an unofficial co-instructor the second time. Jon Borwein attended the Project Presentation Day (the second time) and gave thoughtful feedback to each student. This paper is dedicated to his memory.
△ Less
Submitted 7 February, 2018; v1 submitted 16 January, 2018;
originally announced January 2018.
-
Minimal height companion matrices for Euclid polynomials
Authors:
Eunice Y. S. Chan,
Robert M. Corless
Abstract:
We define Euclid polynomials $E_{k+1}(λ) = E_{k}(λ)\left(E_{k}(λ) - 1\right) + 1$ and $E_{1}(λ) = λ+ 1$ in analogy to Euclid numbers $e_k = E_{k}(1)$. We show how to construct companion matrices $\mathbb{E}_k$, so $E_k(λ) = \operatorname{det}\left(λ\mathbf{I} - \mathbb{E}_{k}\right)$, of height 1 (and thus of minimal height over all integer companion matrices for $E_{k}(λ)$). We prove various prop…
▽ More
We define Euclid polynomials $E_{k+1}(λ) = E_{k}(λ)\left(E_{k}(λ) - 1\right) + 1$ and $E_{1}(λ) = λ+ 1$ in analogy to Euclid numbers $e_k = E_{k}(1)$. We show how to construct companion matrices $\mathbb{E}_k$, so $E_k(λ) = \operatorname{det}\left(λ\mathbf{I} - \mathbb{E}_{k}\right)$, of height 1 (and thus of minimal height over all integer companion matrices for $E_{k}(λ)$). We prove various properties of these objects, and give experimental confirmation of some unproved properties.
△ Less
Submitted 12 December, 2017;
originally announced December 2017.
-
Symbolic-Numeric Integration of Rational Functions
Authors:
Robert M. Corless,
Robert H. C. Moir,
Marc Moreno Maza,
Ning Xie
Abstract:
We consider the problem of symbolic-numeric integration of symbolic functions, focusing on rational functions. Using a hybrid method allows the stable yet efficient computation of symbolic antiderivatives while avoiding issues of ill-conditioning to which numerical methods are susceptible. We propose two alternative methods for exact input that compute the rational part of the integral using Hermi…
▽ More
We consider the problem of symbolic-numeric integration of symbolic functions, focusing on rational functions. Using a hybrid method allows the stable yet efficient computation of symbolic antiderivatives while avoiding issues of ill-conditioning to which numerical methods are susceptible. We propose two alternative methods for exact input that compute the rational part of the integral using Hermite reduction and then compute the transcendental part two different ways using a combination of exact integration and efficient numerical computation of roots. The symbolic computation is done within BPAS, or Basic Polynomial Algebra Subprograms, which is a highly optimized environment for polynomial computation on parallel architectures, while the numerical computation is done using the highly optimized multiprecision rootfinding package MPSolve. We show that both methods are forward and backward stable in a structured sense and away from singularities tolerance proportionality is achieved by adjusting the precision of the rootfinding tasks.
△ Less
Submitted 25 October, 2018; v1 submitted 5 December, 2017;
originally announced December 2017.
-
A Distributed Algorithm for Least Square Solutions of Linear Equations
Authors:
Xuan Wang,
Jingqiu Zhou,
Shaoshuai Mou,
Martin J. Corless
Abstract:
A distributed discrete-time algorithm is proposed for multi-agent networks to achieve a common least squares solution of a group of linear equations, in which each agent only knows some of the equations and is only able to receive information from its nearby neighbors. For fixed, connected, and undirected networks, the proposed discrete-time algorithm results in each agents solution estimate to co…
▽ More
A distributed discrete-time algorithm is proposed for multi-agent networks to achieve a common least squares solution of a group of linear equations, in which each agent only knows some of the equations and is only able to receive information from its nearby neighbors. For fixed, connected, and undirected networks, the proposed discrete-time algorithm results in each agents solution estimate to converging exponentially fast to the same least squares solution. Moreover, the convergence does not require careful choices of time-varying small step sizes.
△ Less
Submitted 28 September, 2017;
originally announced September 2017.
-
On the SPRification of linear descriptor systems via output feedback
Authors:
Martin Corless,
Ezra Zeheb,
Robert Shorten
Abstract:
We consider input-output systems in descriptor form and ask when such systems can be rendered SPR (strictly positive real) via output feedback. Time and frequency domain conditions are given to determine when and how this is possible. In addition, a synthesis procedure for controller design is also derived. Together, the results provide a complete answer to when a linear descriptor system can be m…
▽ More
We consider input-output systems in descriptor form and ask when such systems can be rendered SPR (strictly positive real) via output feedback. Time and frequency domain conditions are given to determine when and how this is possible. In addition, a synthesis procedure for controller design is also derived. Together, the results provide a complete answer to when a linear descriptor system can be made SPR via output feedback, and give a recipe for design of the feedback controller when it exists. Simple examples are given to illustrate our results and to demonstrate their efficacy.
△ Less
Submitted 6 April, 2017;
originally announced April 2017.
-
Stability tests for a class of switched descriptor systems with non-homogenous indices
Authors:
Shravan Sajja,
Martin Corless,
Ezra Zeheb,
Robert Shorten
Abstract:
In this paper we derive stability conditions for a switched system where switching occurs between linear descriptor systems of different indices. In particular, our results can be used to analyse the stability of the important case when switching between a standard system and an index one descriptor system, and systems where switching occurs between an index one and and an index two descriptor sys…
▽ More
In this paper we derive stability conditions for a switched system where switching occurs between linear descriptor systems of different indices. In particular, our results can be used to analyse the stability of the important case when switching between a standard system and an index one descriptor system, and systems where switching occurs between an index one and and an index two descriptor system. Examples are given to illustrate the use of our results.
△ Less
Submitted 18 May, 2017; v1 submitted 4 April, 2017;
originally announced April 2017.
-
Gamma and Factorial in the Monthly
Authors:
Jonathan M. Borwein,
Robert M. Corless
Abstract:
The Monthly has published roughly fifty papers on the $Γ$ function or Stirling's formula. We survey those papers (discussing only our favourites in any detail) and place them in the context of the larger mathematical literature on $Γ$.
The Monthly has published roughly fifty papers on the $Γ$ function or Stirling's formula. We survey those papers (discussing only our favourites in any detail) and place them in the context of the larger mathematical literature on $Γ$.
△ Less
Submitted 15 March, 2017;
originally announced March 2017.
-
Backward Error Analysis for Perturbation Methods
Authors:
Robert M. Corless,
Nicolas Fillion
Abstract:
We demonstrate via several examples how the backward error viewpoint can be used in the analysis of solutions obtained by perturbation methods. We show that this viewpoint is quite general and offers several important advantages. Perhaps the most important is that backward error analysis can be used to demonstrate the validity of the solution, however obtained and by whichever method. This include…
▽ More
We demonstrate via several examples how the backward error viewpoint can be used in the analysis of solutions obtained by perturbation methods. We show that this viewpoint is quite general and offers several important advantages. Perhaps the most important is that backward error analysis can be used to demonstrate the validity of the solution, however obtained and by whichever method. This includes a nontrivial safeguard against slips, blunders, or bugs in the original computation. We also demonstrate its utility in deciding when to truncate an asymptotic series, improving on the well-known rule of thumb indicating truncation just prior to the smallest term. We also give an example of elimination of spurious secular terms even when genuine secularity is present in the equation. We give short expositions of several well-known perturbation methods together with computer implementations (as scripts that can be modified). We also give a generic backward error based method that is equivalent to iteration (but we believe useful as an organizational viewpoint) for regular perturbation.
△ Less
Submitted 5 September, 2016;
originally announced September 2016.
-
Narayana, Mandelbrot, and A New Kind of Companion Matrix
Authors:
Eunice Y. S. Chan,
Robert M. Corless
Abstract:
We demonstrate a new kind of companion matrix, for polynomials of the form $c(λ) = λa(λ)b(λ) + c_0$ where upper Hessenberg companions are known for the polynomials $a(λ)$ and $b(λ)$. This construction can generate companion matrices with smaller entries than the Fiedler or Frobenius forms. This generalizes Piers Lawrence's Mandelbrot companion matrix. We motivate the construction by use of Narayan…
▽ More
We demonstrate a new kind of companion matrix, for polynomials of the form $c(λ) = λa(λ)b(λ) + c_0$ where upper Hessenberg companions are known for the polynomials $a(λ)$ and $b(λ)$. This construction can generate companion matrices with smaller entries than the Fiedler or Frobenius forms. This generalizes Piers Lawrence's Mandelbrot companion matrix. We motivate the construction by use of Narayana-Mandelbrot polynomials, which are also new to this paper.
△ Less
Submitted 20 June, 2016;
originally announced June 2016.
-
On the stability and convergence of a class of consensus systems with a nonlinear input
Authors:
Mingming Liu,
Fabian Wirth,
Martin Corless,
Robert Shorten
Abstract:
We consider a class of consensus systems driven by a nonlinear input. Such systems arise in a class of IoT applications. Our objective in this paper is to determine conditions under which a certain partially distributed system converges to a Lur'e-like scalar system, and to provide a rigorous proof of its stability. Conditions are derived for the non-uniform convergence and stability of such a sys…
▽ More
We consider a class of consensus systems driven by a nonlinear input. Such systems arise in a class of IoT applications. Our objective in this paper is to determine conditions under which a certain partially distributed system converges to a Lur'e-like scalar system, and to provide a rigorous proof of its stability. Conditions are derived for the non-uniform convergence and stability of such a system and an example is given of a speed advisory system where such a system arises in real engineering practice.
△ Less
Submitted 11 July, 2017; v1 submitted 3 June, 2015;
originally announced June 2015.
-
Nonhomogeneous Place-Dependent Markov Chains, Unsynchronised AIMD, and Network Utility Maximization
Authors:
Fabian Wirth,
Sonja Stuedli,
Jia Yuan Yu,
Martin Corless,
Robert Shorten
Abstract:
We present a solution of a class of network utility maximization (NUM) problems using minimal communication. The constraints of the problem are inspired less by TCP-like congestion control but by problems in the area of internet of things and related areas in which the need arises to bring the behavior of a large group of agents to a social optimum. The approach uses only intermittent feedback, no…
▽ More
We present a solution of a class of network utility maximization (NUM) problems using minimal communication. The constraints of the problem are inspired less by TCP-like congestion control but by problems in the area of internet of things and related areas in which the need arises to bring the behavior of a large group of agents to a social optimum. The approach uses only intermittent feedback, no inter-agent communication, and no common clock.
The proposed algorithm is a combination of the classical AIMD algorithm in conjunction with a simple probabilistic rule for the agents to respond to a capacity signal. This leads to a nonhomogeneous Markov chain and we show almost sure convergence of this chain to the social optimum.
△ Less
Submitted 19 March, 2016; v1 submitted 20 April, 2014;
originally announced April 2014.
-
Stieltjes, Poisson and other integral representations for functions of Lambert $W$
Authors:
German A. Kalugin,
David J. Jeffrey,
Robert M. Corless
Abstract:
We show that many functions containing $W$ are Stieltjes functions. Explicit Stieltjes integrals are given for functions $1/W(z)$, $W(z)/z$, and others. We also prove a generalization of a conjecture of Jackson, Procacci & Sokal. Integral representations of $W$ and related functions are also given which are associated with the properties of their being Pick or Bernstein functions. Representations…
▽ More
We show that many functions containing $W$ are Stieltjes functions. Explicit Stieltjes integrals are given for functions $1/W(z)$, $W(z)/z$, and others. We also prove a generalization of a conjecture of Jackson, Procacci & Sokal. Integral representations of $W$ and related functions are also given which are associated with the properties of their being Pick or Bernstein functions. Representations based on Poisson and Burniston--Siewert integrals are given as well.
△ Less
Submitted 27 March, 2011;
originally announced March 2011.
-
On the inversion of $y^αe^y$ in terms of associated Stirling numbers
Authors:
David J. Jeffrey,
Robert M. Corless,
David E. G. Hare,
Donald E. Knuth
Abstract:
The function $y=Φ_α(x)$, the solution of $y^αe^y=x$ for $x$ and $y$ large enough, has a series expansion in terms of $\ln x$ and $\ln\ln x$, with coefficients given in terms of Stirling cycle numbers. It is shown that this expansion converges for $x>(αe)^α$ for $α\ge 1$. It is also shown that new expansions can be obtained for $Φ_α$ in terms of associated Stirling numbers. The new expansions con…
▽ More
The function $y=Φ_α(x)$, the solution of $y^αe^y=x$ for $x$ and $y$ large enough, has a series expansion in terms of $\ln x$ and $\ln\ln x$, with coefficients given in terms of Stirling cycle numbers. It is shown that this expansion converges for $x>(αe)^α$ for $α\ge 1$. It is also shown that new expansions can be obtained for $Φ_α$ in terms of associated Stirling numbers. The new expansions converge more rapidly and on a larger domain.
△ Less
Submitted 30 November, 1995;
originally announced December 1995.