Graduate Course in Numerical Analysis

SEMINAR SERIES ANNOUNCEMENT

A seminar series will be held at the Katholieke Universiteit Leuven and the Universite Catholique de Louvain on the topic of

NUMERICAL METHODS FOR SYSTEMS AND CONTROL

This seminar series represents the graduate courses "Special Topics in Linear Algebra", P. Van Dooren (UCL) and "Numerical Analysis", R. Cools (KUL).

PROGRAM

The course comprises a total of 15 hours of lectures, and consists of one afternoon of three one-hour lectures, given on December 14, 1998 at UCL and subsequently a sequence of several two-hour seminars given at both universities in the second semester (February and March 1999). A detailed list of subjects and speakers is given below.
The subjects treated include:

• Identification
• Structured matrix problems
• Robustness and distance problems
• DAE, DDE, and ODE techniques

REGISTRATION

Participation is free of charge. For organisational purposes, announcements of any last minute changes to the program, and the distribution of the abstracts we do require registration. Please send your name, address and email to either adhemar.bultheel at cs.kuleuven.be or vandooren at anma.ucl.ac.be before December 10, 1998.

VENUE

The lectures at the U.C.L. will take place at Auditoire Euler, Av. G. Lemaitre, 1348 Louvain-la-Neuve (see instructions " how to get to CESAME") except for the lectures of 14/12/98 that take place in Croix du Sud : Room 11 (see instructions " how to get there).

The lectures at the K.U.Leuven will be held at the at the campus Celestijnenlaan 200 and 300, 3001 Heverlee (see instructions "how to get to the lecture rooms").
C300.02.030 = Celestijnenlaan 300, 3001 Heverlee, building Mechanica, 2nd floor, room 030
C200.C.01.D = Celestijnenlaan 200, 3001 Heverlee, building C, 1st floor, room D

The two-hour lectures in 1999 start at 14:00 and are interrupted by a break around 15:00.

PROGRAM DETAILS

Mo December 14, 1998 (at the UCL, Croix du Sud, Room 11) (3 times 1 hour)

• 14:00 Introduction
• 14:15 Gene Golub (Stanford University): Reconstruction of a Polygon from its Moments
• 15:15 Break
• 15:45 Adhemar Bultheel (K.U.Leuven): Rational approximation in systems and control
• 16:45 Paul Van Dooren (U.C.L.): The Schur decomposition in systems and control

Mo February 8, 1999 (at the K.U.Leuven, C300.02.030)

• Bart de Moor (K.U.Leuven): Subspace based methods for identification

Tu February 16, 1999 (at the U.C.L., Aud. Euler)

• Paul Van Dooren (U.C.L.): Stability radii of dynamical systems

Mo February 22, 1999 (at the K.U.Leuven, C300.02.030)

• Marc Van Barel/Peter Kravanja (K.U.Leuven): Structured matrix problems

Tu March 9, 1999 (at the U.C.L., Aud. Euler)

• Volker Mehrmann (T.U. Chemnitz): DAE and ODE techniques

Mo March 15, 1999 (at the K.U.Leuven, C200.C.01.D)

• Dirk Roose and Tatyana Luzyanina (K.U.Leuven): Numerical methods for delay differential equations : time integration, stability analysis, computation of periodic solutions

Tu March 23, 1999 (at the U.C.L., Aud. Euler)

• Sabine Van Huffel (K.U.Leuven): Total least squares and identification

FURTHER INFORMATION

For any additional information please contact:

ABSTRACTS OF THE TALKS AND COURSE MATERIAL

Mo December 14, 1998 (at the UCL, audit. Euler) (3 times 1 hour)

• Gene Golub (Stanford University): Reconstruction of a Polygon from its Moments

  Computation of certain kinds of numerical quadratures on polygonal regions of the plane and the reconstruction of these regions from their moments can be viewed as dual problems. In fact, this is a consequence of a little-known result of Motzkin and Schoenberg. In this talk, we discuss this result and address the inverse problem of (uniquely) reconstructing a polygonal region in the complex plane from a finite number of its complex moments. Algorithms have been developed for polygon reconstruction from moments and have been applied to tomographic image reconstruction problems. The numerical computations involved in the algorithm can be very ill-conditioned. We have managed to improve the algorithms used, and to recognize when the problem will be ill-conditioned. Some numerical results will be given.
  ( the paper in compressed postscript)

• Adhemar Bultheel (K.U.Leuven): Rational approximation in systems and control

  Rational approximation has always played an important role in linear time-invariant systems theory and signal processing. Some classical examples are linear prediction of time series, system identification, realization theory, H-infinity control etc. The mathematical equivalent corresponds to moment problems, orthogonal polynomials, interpolation problems in the complex plane, continued fractions, fast algorithms for structured problems in linear algebra, etc. We give a brief introduction to some of these topics and sketch recent evolutions.
  (transparencies in compressed postscript)

• Paul Van Dooren (U.C.L.): The Schur decomposition in systems and control

  Canonical forms of state-space models of multivariable linear systems have often been proposed for solving certain analysis and design problems encountered in linear system theory. We show that for many problems one can as well make use of so-called condensed forms, which can be obtained under orthogonal state-space transformations. The Schur form is one of the most useful ones in this context and is recommended for problems of pole placement, optimal control, and for solving Sylvester, Lyapunov and Riccati equations. We also present an extension of this decompositions for a sequence of matrices which we call the periodic Schur decomposition.
  (text in compressed postscript)

• Bart de Moor (K.U.Leuven): Subspace based methods for identification

  Mathematical models of dynamical systems are used for analysis, simulation, prediction, optimization, monitoring, fault detection, signal modelling, filtering and detection, training and control. Engineering solutions in telecommunications, biomedical signal processing, industrial process control etc. are increasingly based on mathematical models of systems. For instance, many problems in biomedical and telecommunications signal processing can be phrased as a model-based filtering and/or detection problem. The design of industrial process control systems, especially for multivariable ones, is critically based upon a reliable model of the process at hand.

  It goes without saying that state-space models are good engineering models in many of these applications.

  In this tutorial, we give a survey of so-called subspace identification techniques for multivariable linear time-invariant models (LTI) of the form

  x(k+1) = A x(k) + B u(k) + w(k) ,
  y(k) = C x(k) + D u(k) + v(k) ,

  where k is the discrete time index, u(k) and y(k) are the (observed, measured) inputs respectively outputs of the system, x(k) is the (unknown) systems state, w(k) and v(k) are unobserved white noise inputs (process noise, resp. measurement noise). The identification problem now consists of finding the real model matrices A, B, C and D and the noise covariance matrices from observations of the inputs u(k) and outputs y(k) only, when a sufficient amount of data is available ( k * infinity).

  Since the beginning of the nineties, there has been an increasing interest (witnessed by numerous papers and books) in so-called subspace identification methods, in which geometrical, algebraic and numerical concepts and algorithms are elegantly intertwined, leading to extremely user-friendly software for obtaining LTI models from input-output data.

  It is the purpose of this tutorial to

  • introduce the audience to the basic geometrical, algebraic and numerical ingredients of subspace methods;
  • organise the conceptual understanding of the building blocks of subspace algorithms, so that the user can "understand" the different variants that have been described in the literature;
  • demonstrate that subspace methods lead to elegant numerically reliable algorithms that have proved useful in many applications and provide him/her with information about software implementations;
  • provide the attending researcher with a state-of-the-art review of the most relevant milestones in the development of subspace methods as well as survey recently obtained results and current research activities in subspace identification;

  Subspace methods derive their user-friendliness and robustness from the fact that

  • There is no need for a user-defined parametrization, the determination of which is a notoriously difficult problem for MIMO systems; In a subspace method, typically, one of the only critical parameters to be determined is the number of states; Said in other words, subspace methods lead to fully (over-)parametrized models (but the overparametrization does not lead to additional difficulties). As an additional feature, there is no need to parametrize the initial state (as is the case with identification methods that are based on input-output models);
  • Subspace methods have led to new conceptual insights in (linear) system identification, such as the insight on how to obtain a (Kalman filter) state estimate directly from inputs and outputs, or on how to obtain reduced models "directly" (i.e. without calculating the "larger" model first);
  • In subspace methods, there are typically three basic conceptual steps:
  • First an estimate of the (Kalman filter) state sequence and/or the column space of the observability matrix is obtained; The number of states derives from the rank of certain oblique and orthogonal projections of subspaces generated by the data;
  • Next the model matrices are obtained by solving a least squares problem;
  • The noise covariance matrices are estimated from the least squares residuals;

  There is a large number of possible variations on and within these basic steps, each of which leads to a variant that has been described in the literature.

  • Each of these steps can be robustly and efficiently implemented, resulting in the fact that subspace algorithms are "nice", consisting of well understood basic building blocks such as the QR and the singular value decomposition, hence avoiding problems of ill convergence as in non-linear iterative optimization algorithms; This 2-hour tutorial is aimed at providing both a conceptual insight in subspace identification methods as well as a survey of more recent developments. The intended audience consists of
  • PhD students and researchers active in system identification;
  • Regular or occasional users of system identification tools (signal processing and/or control engineers), both from academia or industry, who want to get acquainted with
  • The basic ideas of subspace identification techniques;
  • The user-friendliness and robustness of subspace software;

  Basic reference

  Peter Van Overschee, Bart De Moor.
  Subspace Identification for Linear Systems: Theory, Implementation and Applications. Kluwer Academic Publishers, 1996, 254 pp. ISBN 0-7923-9717-7 (http://www.wkap.nl) Matlab floppy disk with subspace identification m.files included.

  Outline of the tutorial

  Part I: Subspace identification for deterministic and stochastic systems (3/4 hour)

  • Basic ingredients of subspace methods
  • Geometry:
  • Row and column spaces of matrices, dimension;
  • Orthogonal and oblique projections;
  • System theory:
  • Block Hankel matrices with inputs and outputs;
  • Realization algorithms;
  • The Kalman filter state directly from inputs and outputs;
  • Numerical issues
  • QR-decomposition;
  • Rank deficiency;
  • Singular Value decomposition;
  • Subspace identification of deterministic systems
  • The state as intersection between past and future;
  • Calculating the model matrices;
  • Subspace identification of stochastic systems
  • Backward and forward innovation models;
  • The Kalman filter state estimate as an orthogonal projection from the future onto the past;
  • Canonical correlations;
  • The issue of positive realness;

  Part II: Subspace identification for combined systems (3/4 hour)

  • The Kalman filter state estimate as an oblique projection;
  • Computing the model matrices and noise covariances;
  • Weighted oblique projections and special cases (MOESP, N4SID, CVA, others described in literature,...)
  • Direct model reduction; Enns's conjecture;
  • User choices;

  Part III: Implementations, applications, extensions and open problems (1/2 hour)

  • Implementations:
    • How to numerically implement robust reliable software for subspace identificaton; Numerical issues; large problems (large number of inputs, outputs, states); Updating issues; Graphical User Interfaces;
    • Survey of subspace identification codes, both the commercially available ones (Matlab, Xmath, Adaptics,...) as the ones available in public domain (SLICOT/NICONET,...);
  • Applications
    • We will illustrate properties and results with applications from process industry, mechatronics and signal processing from datasets publically available at http://www.esat.kuleuven.be/sista/daisy.
  • Extensions Review of extensions to periodic systems, bilinear system identification, system identification for nonlinear system (e.g. Wiener/Hammerstein), frequency domain versions, closed-loop subspace identification, including a priori information (e.g. enforcing stability,...), blind identification, stochastic identification with non-stationary inputs;
  • Recent results and open problems Here we give a short (literature) survey with respect to the current research on subspace identification (including some recently obtained statistical, system theoretic (algebraic and positive degree), genericity (input conditions) and numerical results); We also formulate the most challenging open problems;

  (the transparancies in compressed postscript)

• Paul Van Dooren (U. C.L.): Stability radii of dynamical systems

  The basic stability radius problem can be defined as follows. Let A be a stable matrix with all its eigenvalues in the stability region S (typically the open left half plane or the open unit disc). The stability radius of A is the norm of smallest (complex) perturbation D that causes at least one eigenvalue of A+D to become unstable. This problem has extensively been studied and there exists an analytic formulation for r_C as well as a quadratically convergent algorithm to compute it. We discuss several variants and extensions of this standard problem to the generalized eigenvalue problem, polynomial matrices, periodic eigenvalue problems, and structured matrices. We discuss the use of different norms as well as structured perturbations and link this problem to the concept of pseudospectra. We will also talk about computational aspects.

  (two papers: paper 1, paper 2 in compressed postscript)

• Marc Van Barel/Peter Kravanja (K.U.Leuven): Structured matrix problems

  Part I: Discrete linearized least squares rational approximation on the unit circle

  Given the abscissae z_i, i = 1,2,...,m on the unit circle in the complex plane and the corresponding function values f_i. We look for a rational function n(z)/d(z) such that n(z_i)/d(z_i) » f_i, i=1,2,...,m. More precisely, we look for the minimum of å_i=1^m w_i | r_i|² with r_i:= f_i d(z_i) - n(z_i) and w_i some positive weight. In this talk, we show how this problem can be solved in an efficient and accurate way.

  (transparancies part 1 in compressed postscript)
  

  References

  1. M. Van Barel and A. Bultheel Discrete linearized least squares approximation on the unit circle . J. Comput. Appl. Math., 50:545-563, 1994.
  2. M. Van Barel and A. Bultheel Orthogonal polynomial vectors and least squares approximation for a discrete inner product . Electron. Trans. Numer. Anal., 3:1-23, 1995.
  3. A. Bultheel and M. Van Barel, Vector orthogonal polynomials and least squares approximation. SIAM J. Matrix Anal. Appl., 16(3):863-885, 1995.

  Part II: A stabilized superfast solver for indefinite Toeplitz systems

  We present a stabilized superfast solver for indefinite Toeplitz systems Tx=b. An explicit formula for T^-1 is expressed in such a way that the matrix-vector product T^-1b can be calculated via FFTs and Hadamard products. This inversion formula involves certain polynomials that can be computed by solving two linearized rational interpolation problems on the unit circle. The heart of our Toeplitz solver is a superfast algorithm to solve these interpolation problems. This algorithm is stabilized via pivoting, iterative improvement, downdating, and by giving ``difficult'' interpolation points an adequate treatment. A Fortran 90 implementation of our algorithm is available. Large scale numerical examples will illustrate the effectiveness of our approach.
  ( text part 2 in compressed postscript and transparancies part 2 in compressed postscript )

• Volker Mehrmann (T.U. Chemnitz): DAE and ODE techniques

  We discuss differential-algebraic systems and control problems in descriptor form. Based on a new canonical form we present the analysis as well as the numerical techniques for the numerical integration of DAEs as well as the analysis and numerical treatment of control problems for DAEs. The general approach allows to study over- and underdetermined systems and thus we can use this approach also to study control problems in a behaviour setting.

  We also discuss an the corresponding software packages for the solution of these problems as well as their extensions to desriptor control problems.

  References

  1. P. Kunkel und V. Mehrmann, A New Class of Discretization Methods for the Solution of Linear Differential Algebraic Equations with Variable Coefficients, SIAM Journal Numerical Analysis, Vol. 33, No.5, October 1996, pp. 1941-1961.
  2. P. Kunkel, V. Mehrmann, W. Rath und J. Weickert, GELDA: A Software Package for the Solution of General Linear Differential Algebraic equations, SIAM Journal Scientific Computing, Vol. 18, 1997, pp. 115 - 138.
  3. P. Kunkel und V. Mehrmann, Regular Solutions of Nonlinear Differential-Algebraic Equations and their Numerical Determination, Numerische Mathematik, Vol. 79, 1998, pp. 581--600.
  4. P. Kunkel, V. Mehrmann und W. Rath, Analysis and Numerical Solution of Control Problems in Descriptor Form Preprint SFB393/99-01, Sonderforschungsbereich 393, `Numerische Simulation auf massiv parallelen Rechnern', Fakultät für Mathematik, TU Chemnitz.

• Dirk Roose and Tatyana Luzyanina (K.U.Leuven): Numerical methods for delay differential equations : time integration, stability analysis, computation of periodic solutions

  In a delay differential equation, not only the present state of a physical system is taken into account, but also its past history. Different ways to incorporate the history can be considered: at one or more time-points in the past (discrete delays) or as an integral over the model's past state (distributed delays). Also other types of delay are possible, e.g. a state-dependent delay. If the history is present not only through the state variable but also appears through the time derivative, the model is called a neutral differential equation.

  In this lecture, we will give an overview of numerical methods for differential equations with delay. We will first discuss the stability analysis of steady st ate solutions, which requires the solution of a nonlinear eigenvalue problem. W e then describe time integration methods, which are adaptations of time integrat ors for ordinary differential equations, and we give information about available software packages. Finally, we discuss numerical methods for the the computation and the stability analysis of parameter dependent nonlinear delay and neutral differential equations, based on a shooting approach or on collocation.

  We will illustrate the numerical methods with practical examples, incl. results from the analysis of a closed loop system with a nonlinear control law, in which a delay in the feedback loop is taken into account.

  References

  1. C.T.H. Baker and C.A.H. Paul and D.R. Willé Issues in the numerical solution of evolutionary delay differential equations . Advances in Computational Mathematics, vol 3 (3), 171-196, 1995. Also Numerical Analysis Report, 248, Univ. of Manchester, Dept. mathematics, April 1994.
  2. K. Engelborghs, D. Roose Numerical computation of stability and detection of Hopf bifurcations of steady state solutions of delay differential equations . Report TW274 Dept. Computer Sci. K.U.Leuven, 1998. (To appear in Adv. Comput. Math.)
  3. T. Luzyanina, K. Engelborghs, K. Lust, D. Roose Computation, continuation and bifurcation analysis of periodic solutions of delay differential equations . International Journal of Bifurcation and Chaos vol. 7, number 11, p. 2547-2560, 1997. Also Report TW252, Dept. Computer Sci. K.U.Leuven, 1997.

• Sabine Van Huffel (K.U.Leuven): Total least squares and identification

  The Total Least Squares (TLS) approach is a popular method in data fitting problems where we have to solve an overdetermined system of linear equations A x » b. When the noise on the entries of [A b] is Gaussian i.i.d. white noise of equal variance, TLS yields a Maximum Likelihood (ML) estimate for x, DA and Db such that (A + DA)x = b+Db. After a description of its basic principle, the algebraic and statistical properties of the TLS technique are treated, as well as its computational aspects and merits in practical applications.

  Furthermore, generalizations of the basic TLS principle are presented. In particular,the above-mentioned statistical condition is not satisfied in many system identification and signal processing applications where we have to deal with structured (Hankel, Toeplitz, sparse matrices,...) data matrices [A b] with i.i.d. Gausssian white noise of equal variance on the different entries. A necessary condition for ML estimates in that case is that the structure of [DA Db] equals that of [A b]. Therefore several formulations and corresponding solution methods have been developed, that allow to fix the structure of [DA Db], while minimizing ||[DA Db]||_F.
  We present these so-called Structured TLS approaches, and discuss the properties (convergence rate and accuracy) of the algorithms. Finally, the application of the Structured TLS method to various identification is studied in which the computed correction matrix applied to A or [A b] keeps the same Toeplitz structure as the data matrix A or [A b], respectively. In particular, the Structured TLS method is compared with the LS and TLS methods in deconvolution, transfer function modeling and linear prediction problems, and shown to improve the accuracy of the parameter estimates by a factor of 2 to 40 at any signal-to-noise ratio.
  (the paper in compressed postscript)