Author Archives: Nick Higham

The Contribution of Dr. J. H. Wilkinson to Numerical Analysis

Posted on by Leave a comment
Wilkinson_033.jpg

President HRH Duke of Edinburgh presenting honorary fellowship of The Institute of Mathematics and its Applications to James Wilkinson in 1977 (© The IMA).

The title of this post is the same as that of a symposium organized by Michael J. D. Powell and the Institute of Mathematics and its Applications (IMA) at the Royal Society in London on July 6th, 1977. The meeting commemorated the election of James Hardy Wilkinson to an Honorary Fellowship of the IMA.

The proceedings of the meeting were published by the IMA in a 91-page A5 booklet. As far as I am aware, few copies of the booklet survive and its contents have not previously been made available online. I am grateful to David Youdan, Executive Director of the IMA, for giving me permission to provide here a scan of the booklet. It is timely to do so, because this year marks the 100th anniversary of the birth of Wilkinson.

Here are the individual chapters, with comments from Mike Powell’s preface in quotes.

  • About Jim Wilkinson, with a Commemorative Snippet on Backward Error Analysis, L. Fox (Oxford University Computing Laboratory). “Leslie Fox describes many of Jim Wilkinson’s achievements that have not been published before and he exposes the accuracy of some ill-conditioned least squares calculations.”
  • Inverse Iteration, Newton’s Method, and Non-linear Eigenvalue Problems, M. R. Osborne (Australian National University). “Mike Osborne unifies the convergence properties of a main class of iterative methods for calculating eigenvalues.”
  • A New Look at Error Analysis, C. W. Clenshaw (University of Lancaster). “Charles Clenshaw develops an idea, due to Frank Olver, for treating the accumulation of errors in floating point arithmetic.”
  • A Problem in Numerical Linear Algebra, J. H. Wilkinson (National Physical Laboratory). “Jim Wilkinson shows the relevance in practice of the equivalence of repeated matrix eigenvalues, the ill-conditioning of the matrix eigenvector calculation, and the orthogonality of left and right hand eigenvectors that have a common eigenvalue.”

For more information about Wilkinson, see this web page that Sven Hammarling and I have created.

ima78_cover.jpg

Wilkinson Quotes

Posted on by Leave a comment
Wilkinson_021.jpg

by Sven Hammarling and Nick Higham

We collect here some quotes from the work of Jim Wilkinson. These reflect his unique perspective as a mathematician who was involved in designing and building one of the first digital computers and who subsequently developed and analyzed a variety of numerical algorithms

We have arranged the quotes under the following headings:

Program libraries | Floating-point arithmetic | Rounding error analysis | Conditioning | Backward error analysis | Polynomials | Interaction on Pilot ACE | Communication avoidance | Linear algebra on Pilot ACE

Program libraries

Since the programming is likely to be the main bottleneck in the use of an electronic computer we have given a good deal of thought to the preparation of standard routines of considerable generality for the more important processes involved in computation. By this means we hope to reduce the time taken to code up large-scale computing problems, by building them up, as it were, from prefabricated units. [W48, p. 286]

In spite of the self-contained nature of the linear algebra field, experience has shown that even here the preparation of a fully tested set of algorithms is a far greater task than had been anticipated. [W71a, p. v]

Floating-point arithmetic

At a time when the arithmetic provided on modern computers is often so disappointing, it is salutary to recall that the subroutines included provision for accumulating inner products in double-precision floating-point arithmetic and all rounding was immaculate! [W80, p. 105]

Rounding error analysis

The two main classes of rounding error analysis are not, as my audience might imagine, `backwards’ and `forwards’, but rather `one’s own’ and `other people’s’. One’s own is, of course, a model of lucidity; that of others serves only to obscure the essential simplicity of the matter in hand. [W85, p. 5]

In general, the statistical distribution of the rounding errors will reduce considerably the function of n occurring in the relative errors. We might expect in each case that this function should be replaced by something which is no bigger than its square root and is usually appreciably smaller. [W61, p. 38]

For me, then, the primary purpose of the rounding error analysis was insight. [W86, p. 197]

Conditioning

The system was mildly ill-conditioned, though we were not so free with such terms of abuse in those days, [W71, p. 144]

Backward error analysis

“You have been solving these damn problems better than I can pose them.” Sir Edward Bullard, Director NPL, in a remark to Wilkinson (mid 1950s) [W85, p. 11]

I first used backward error analysis in connection with simple programs for computing zeros of polynomials soon after the PILOT ACE came into use. [W85, p. 8]

There does seem to be some misunderstanding about the purpose of an a priori backward error analysis. All too often, too much attention is paid to the precise error bound that has been established. The main purpose of such an analysis is either to establish the essential numerical stability of an algorithm or to show why it is unstable and in doing so to expose what sort of change is necessary to make it stable. The precise error bound is not of great importance. [W74, p. 356]

The great stability of unitary transformations in numerical analysis springs from the fact that both the \ell_2-norm and the Frobenius norm are unitarily invariant. This means in practice that even when rounding errors are made, no substantial growth takes place in the norms of the successive transformed matrices. [W65, p. 77]

Although backward analysis is a perfectly straightforward concept there is strong evidence that a training in classical mathematics leaves one unprepared to adopt it. … I have even detected a note of moral disapproval in the attitude of many to its use and there is a tendency to seek a forward error analysis even when a backward error analysis has been spectacularly successful. [W85, p. 5]

Polynomials

The Fundamental Theorem of Algebra asserts that every polynomial equation over the complex field has a root. It is almost beneath the dignity of such a majestic theorem to mention that in fact it has precisely n roots. [W84, p. 21]

The cosy relationship that mathematicians enjoyed with polynomials suffered a severe setback in the early fifties when electronic computers came into general use. Speaking for myself I regard it as the most traumatic experience in my career as a numerical analyst. [W84, p. 3]

Interaction on Pilot ACE

Since the use of the punched-card equipment required the use of an operator, it encouraged user participation generally, and this was a distinctive feature of Pilot ACE operation. For example, various methods of accelerating the convergence of matrix iterative processes were left under the control of operators, and the skill with which these stratagems were used by young women with no more than high school mathematics qualifications was most impressive. Speaking for myself I gained a great deal of experience from user participation, and it was this that led to my own conversion to backward error analysis. [W80, p. 112]

Communication avoidance

Since all machines have stores of finite size often divided up into high speed and auxiliary sections, storage considerations often have a vitally important part to play. [W55, p. 188]

Linear algebra on Pilot ACE

An interesting feature of these codes is that they make a very intensive use of subroutines; the addition of two vectors, multiplication of a vector by a scalar, inner products, etc, are all coded in this way. [W80, p. 105]

From 1946–1948 a great deal of quite detailed coding was done.… The subroutines for floating-point arithmetic were … produced by Alway and myself in 1947 … They were almost certainly the earliest floating-point subroutines. [W80, pp. 104–105]

References

[W48]   J. H.Wilkinson, The Automatic Computing Engine at the National Physical Laboratory, Proc. Roy. Soc. London Ser. A 195, 285-286, 1948.

[W55]   J. H. Wilkinson, The use of iterative methods for finding the latent roots and vectors of matrices, Mathematical Tables and Other Aids to Computation 9, 184-191, 1955.

[W65]   J. H. Wilkinson, Error Analysis of Transformations Based on the Use of Matrices of the Form I-2ww^H, pages 77-101, in Louis Rall, ed., Error in Digital Computation, vol. 2, Wiley, 1965.

[W61]   J. H. Wilkinson. Error analysis of direct methods of matrix inversion. J. ACM, 8:281-330, 1961.

[W71]   J. H. Wilkinson. Some comments from a numerical analyst. J. ACM, 18:137–147, 1971. (The 1970 A. M. Turing lecture).

[W71a]   J. H. Wilkinson and C. Reinsch, eds, Linear Algebra, II, Springer, 1971.

[W74]   J. H. Wilkinson, Numerical linear algebra on digital computers, IMA Bull. 10, 354-356, 1974.

[W80]   J. H. Wilkinson, Turing’s work at the National Physical Laboratory and the construction of Pilot ACE, DEUCE, and ACE, pages 101-114, in N. Metropolis, J. Howlett and G.-C. Rota, eds, A History of Computing in the Twentieth Century: A Collection of Essays, Academic Press, 1980.

[W84]   James Wilkinson, The Perfidious Polynomial, in G. H. Golub, ed., Studies in Numerical Analysis, Mathematical Association of America, Washington, D.C., 24, pp 1-28, 1984.

[W85]   J. H. Wilkinson. The state of the art in error analysis. NAG Newsletter, 2/85:5–28, 1985. (Invited lecture for the NAG 1984 Annual General Meeting).

[W86]   J. H. Wilkinson, Error Analysis Revisited, IMA Bull. 22, 192-200, 1986

Wilkinson and Backward Error Analysis

Posted on by Leave a comment

by Sven Hammarling and Nick Higham

It is often thought that Jim Wilkinson developed backward error analysis because of his early involvement in solving systems of linear equations. In his 1970 Turing lecture [5] he described an experience, during world war II at the Armament Research Department, of solving a system of twelve linear equations on a desk computer, using Gaussian elimination. (He doesn’t say how long it took, but it must surely have been several days.) The coefficients were of order unity and, using ten decimal digit computation, he found that the coefficients of the reduced equation determining x_{12} had four leading zeros, so he felt that the solutions could surely have no more than six correct figures. As a check on his calculations, he then computed the residuals and to his surprise the left hand sides agreed with the right hand sides to the full ten figures.

Wilkinson_slide2_cropped.jpg

A slide of Wilkinson’s describing the solution of the system of 18 equations.

After the war Wilkinson joined the Mathematics Division at the National Physical Laboratory. Soon after his arrival, a system of eighteen equations were given to the Mathematics Division. This required a joint effort, which was manned by Leslie Fox, Eric Goodwin, Alan Turing and Wilkinson. Again, the solution was somewhat ill conditioned, as revealed by the final reduced equation, but again in computing the residuals the right and left hand sides agreed to full accuracy. Incidentally, Wilkinson and his colleagues used iterative refinement, which convinced them that the first solution had been accurate to six figures.

These experiences did not straightforwardly lead Wilkinson to develop backward error analysis. In [7] he says that he first used backward error analysis in connection with simple programs for computing zeros of polynomials soon after PILOT ACE came into use, specifically a program for evaluating a polynomial by nested multiplication and a program for carrying out polynomial deflation. But he nevertheless did not recognise backward error analysis as a general tool. Wilkinson explains

It is natural to ask why I did not immediately set about using this type of error analysis as a general purpose tool. In retrospect it seems amazing that I did not try it on Gaussian elimination and on various eigenvalue algorithms in which I was keenly interested at the time. The truth is that it did not occur to me for one moment to do so.

His explicit recognition of a tool that he decided to call backward error analysis soon came through his experience of solving eigenvalue problems on PILOT ACE. He states further in [7]:

Because of the small storage capacity of the PILOT ACE virtually the only algorithm that could be used for dealing with large unsymmetric eigenvalue problems was the power method supplemented by various techniques for accelerating convergence. After each eigenvalue/eigenvector was determined this pair was removed by deflation. Now at that time deflation was generally held to be extremely unstable and accordingly I used it at first with great trepidation. However, it soon became evident that it was being remarkably effective.

As with the linear equation problem, Wilkinson computed residuals, r = A\hat{x} - \hat{\lambda} \hat{x}, where \hat{x} and \hat{\lambda} are the computed values, with \hat{x} normalised so that \hat{x}^T \hat{x} = 1, and, even after many deflations, he found that the residuals were remarkably small. He then realised that

(A - r \hat{x}^T) \hat{x} = \hat{\lambda} \hat{x} \;\;\; \mbox{exactly}.

and this led him directly to the backward error analysis since, if we put E = -r \hat{x}^T, then \hat{\lambda} and \hat{x} are an exact eigenvalue and eigenvector of the matrix A + E. He now recognised that the process could be widely used and this, of course, led to his 1963 book Rounding Errors in Algebraic Processes [3] and, soon after, to The Algebraic Eigenvalue Problem [4].

It should be noted that Wilkinson did not claim to be the first to perform a backward error analysis. He attributes the first analysis to von Neumann and Goldstine in their 1947 paper [2], which, as Wilkinson said in his von Neumann prize paper “is not exactly bedtime reading” [6]. Wilkinson also gives great credit to Givens for his backward error analysis of orthogonal tridiagonalisation in his, sadly, unpublished technical report [1].

References

[1]   W. Givens. Numerical computation of the characteristic values of a real symmetric matrix. Technical Report ORNL-1574, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37831, USA, 1954.

[2]   J. von Neumann and H. H. Goldstine. Numerical inverting of matrices of high order. Bull. Amer. Math. Soc., 53:1021–1099, 1947.

[3]   J. H. Wilkinson. Rounding Errors in Algebraic Processes. Notes on Applied Science, No.32. HMSO, London, UK, 1963. (Also published by Prentice-Hall, Englewood Cliffs, NJ, USA, 1964. Reprinted by Dover Publications, New York, 1994).

[4]   J. H.Wilkinson. The Algebraic Eigenvalue Problem. Oxford University Press, Oxford, UK, 1965.

[5]   J. H. Wilkinson. Some comments from a numerical analyst. J. ACM, 18:137–147, 1971. (The 1970 A. M. Turing lecture).

[6]   J. H. Wilkinson. Modern error analysis. SIAM Review, 13:548–568, 1971. (The 1970 von Neumann lecture).

[7]   J. H. Wilkinson. The state of the art in error analysis. NAG Newsletter, 2/85:5–28, 1985. (Invited lecture for the NAG 1984 Annual General Meeting).

Computing the Wave-Kernel Matrix Functions

Posted on by Leave a comment

The wave-kernel functions \cosh{\sqrt{A}} and \mathrm{sinhc}{\sqrt{A}} arise in the solution of second order differential equations such as u''(t) - Au(t) = b(t) with initial conditions at t=0. Here, A is an arbitrary square matrix and \mathrm{sinhc}{z} = \sinh(z)/z. The square root in these formulas is illusory, as both functions can be expressed as power series in A, so there are no questions about existence of the functions.

How can these functions be computed efficiently? In Computing the Wave-Kernel Matrix Functions (SIAM J. Sci. Comput., 2018) Prashanth Nadukandi and I develop an algorithm based on Padé approximation and the use of double angle formulas. The amount of scaling and the degree of the Padé approximant are chosen to minimize the computational cost subject to achieving backward stability for \cosh{\sqrt{A}} in exact arithmetic.

In the derivation we show that the backward error of any approximation to \cosh{\sqrt{A}} can be explicitly expressed in terms of a hypergeometric function. To bound the backward error we derive and exploit a new bound for \|A^k\|^{1/k} in terms of the norms of lower powers of A; this bound is sharper than one previously obtained by Al-Mohy and Higham.

Numerical experiments show that the algorithm behaves in a forward stable manner in floating-point arithmetic and is superior in this respect to the general purpose Schur–Parlett algorithm applied to these functions.

principal_domain_cosh_sqrt.jpg

The fundamental regions of the

function cosh(sqrt(z)), needed for the backward error analysis

underlying the algorithm.

Celebrating the Centenary of James H. Wilkinson’s Birth

Posted on by Leave a comment

by Sven Hammarling and Nick Higham

Wilkinson_035.jpg

September 27, 2019 is the 100th anniversary of the birth of James Hardy Wilkinson—the renowned numerical analyst who died in 1986. We are marking this special anniversary year in several ways:

The tag wilkinson lists all the posts in this series.

Numerical Algorithms for High-Performance Computational Science, Royal Society, London, April 8-9, 2019

Posted on by Leave a comment

This 2-day scientific discussion meeting is being held at the Royal Society, London, UK, and is organized by Jack Dongarra, Laura Grigori and Nick Higham.

The programme consists of invited talks and contributed posters. The invited speakers are

  1. Guillaume Aupy, INRIA Bordeaux
  2. Erin Carson, Charles University, Prague
  3. George Constantinides, Imperial College
  4. Steve Furber, FRS, University of Manchester
  5. Mike Heroux, Sandia National Laboratories
  6. Tony Hey, Science and Technology Facilities Council
  7. David Keyes, King Abdullah University of Science and Technology, Saudi Arabia
  8. Doug Kothe, Director of the DOE Exascale Computing Project (ECP)
  9. Satoshi Matsuoka, RIKEN Center for Computational Science, Tokyo
  10. Tim Palmer, FRS, Oxford University
  11. Jack Poulson, Hodge Star Scientific Computing, Toronto
  12. Anna Scaife, University of Manchester
  13. John Shalf, Lawrence Berkeley National Laboratory
  14. Rick Stevens, Argonne National Laboratory
  15. Michela Taufer, University of Delaware
  16. Kathy Yelick, Lawrence Berkeley National Laboratory,

The deadline for submission of poster abstracts is Friday March 1, 2019 (extended from Monday 4 February 2019).

Attendance is free, but places are limited and advance registration is essential.

For more information, including the programme and registration, see https://royalsociety.org/science-events-and-lectures/2019/04/high-performance-computing/

Numerical Linear Algebra Group Activities 2018

Posted on by Leave a comment

The Numerical Linear Algebra Group, many of whom are in the October 2018 photo below, had a busy year. This post summarizes what we got up to. Publications are not included here, but many of them can be found on MIMS EPrints under the category Numerical Analysis.

As well as the new group members mentioned below, Stephanie Lai joined the group in September as Executive Assistant to Nick Higham. Stephanie, who holds an MMath degree from the University of Manchester, helps to maintain this website, which launched this autumn. Martin Lotz left the group in August 2018 and is now an Associate Professor at The University of Warwick.

DSC00120.ARW

Software

Together with Jack Dongarra’s team at the University of Tennessee the group is one of the two partners involved in the development of PLASMA: Parallel Linear Algebra Software for Multicore Architectures.

PhD students Weijian Zhang, Steven Elsworth and Jonathan Deakin continued to develop Etymo—a technology for developing knowledge graphs, including Etymo Scholar, a discovery engine for AI research.

We make our research codes available as open source, principally on GitHub; see the repositories of Connolly, Fasi, Higham, Liu, Pranesh, Relton, Sego, Tisseur, Zhang, Zounon.

We also put MATLAB software on MATLAB Central File Exchange and on our own web sites, e.g., the Rational Krylov Toolbox (RKToolbox).

PhD Students

We welcomed new PhD students Michael Connolly, Xiaobo Liu, and Conor Rogers.

Weijian Zhang successfully defended his PhD on Evolving Graphs and similarity-based Graphs with Applications in October 2018.

Postdoctoral Research Associates (PDRAs)

Theo Mary joined us in January 2018 to work jointly on the ICONIC project.

Mawussi Zounon, who was working on the Parallel Numerical Linear Algebra for Extreme Scale Systems (NLAFET) project, moved over in November 2018 to the position of KTP Associate in a three-year Knowledge Transfer Partnership with NAG.

Jakub Sistek returned to the Czech Academy of Sciences after the completion of the project Programming Model INTERoperability ToWards Exascale (INTERTWinE).

Prashanth Nadukandi completed his Marie Skłodowska-Curie Individual Fellowship in September 2018 and is now a Technical Advisor in Advanced Mathematics at the Repsol Technology Lab, Móstoles, Spain.

Maksims Abalenkovs completed his work on the EPSRC project SERT: Scale-free, Energy-aware, Resilient and Transparent Adaptation of CSE Applications to Mega-core Systems.

Presentations at Conferences and Workshops

UK and Republic of Ireland Section of SIAM Annual Meeting, National Oceanography Centre, Southampton, January 11, 2018: Blanchard, Nadukandi.

Waseda Research Institute for Science and Engineering Institute for Mathematical Science Kickoff Meeting, Waseda University, Tokyo, March 6, 2018: Higham.

SIAM Conference on Parallel Processing for Scientific Computing, Waseda University, Tokyo, Japan, March 7-10, 2018: Hammarling, Higham, Mary, Zounon.

Conference on Scientific Computing and Approximation, Purdue University, March 30-31, 2018: Higham.

SIAM Conference on Applied Linear Algebra, May 4-8, 2018, Hong Kong Baptist University: Elsworth, Güttel, Zemaityte.

Structured Matrix Days 2018, May 14-15, 2018, ENS Lyon, France: Mary, Tisseur.

Mathematical Modelling, Numerical Analysis and Scientific Computing, Kácov, Czech Republic, May 27-June 1, 2018: Higham, Tisseur.

International Conference on Approximation and Matrix Functions, May 31 and June 1, 2018, Université de Lille, France: Güttel.

Bohemian Matrices and Applications, University of Manchester, June 20-22, 2018: Higham.

International Supercomputing conference (ISC18), Frankfurt, June 24-28, 2018: Hammarling, Zounon.

10th International Workshop on Parallel Matrix Algorithms and Applications (PMAA18), June 27-29, 2018, ETH Zurich: Zounon.

6th IMA Conference on Numerical Linear Algebra and Optimization, June 27-29 2018, University of Birmingham: Blanchard, Mary, Tisseur.

SIAM Annual Meeting, Portland, USA, July 9-13, 2018: Blanchard, Higham, Tisseur.

JuliaCon 2018, London, August 7-11, 2018: Higham.

Sparse Day Meeting 2018, September 26-27 2018, Toulouse, France: Mary.

Conference and Workshop Organization

Visitors

Recognition and Service

    • Nick Higham gave the Samuel D. Conte Distinguished Lecture at the Department of Computer Science, Purdue University, USA, in March 2018.
  • Massimilano Fasi was presented with the SIAM Student Chapter Certificate of Recognition 2018.

Other Notable Tweets

A New Approach to Probabilistic Rounding Error Analysis

Posted on by Leave a comment

by Nick Higham and Theo Mary

hima18_fig.jpg

Backward error for random linear system of dimension n.

James Wilkinson developed a systematic way of analyzing rounding errors in numerical algorithms in the 1960s—a time when computers typically used double precision arithmetic and a matrix of dimension n = 100 was considered large. As a result, an error bound with a dominant term p(n)u, for p a low degree polynomial and u the machine precision, was acceptably small.

Today, the supercomputers heading the TOP500 list solve linear systems of equations of dimension 10^8 and half precision arithmetic (accurate to about 4 significant decimal digits) is increasingly available in hardware, notably on graphical processing units (GPUs) from AMD and NVIDIA. Traditional rounding error bounds cannot guarantee accurate results when the dimension is very large or the precision is low. Yet computations are often successful in such cases—for example in machine learning algorithms making use of half, or even lower, precision.

This discrepancy between theory and practice stems from the fact that traditional rounding error bounds are worst-case bounds and so tend to be pessimistic. Indeed, while a single floating-point operation incurs a rounding error bounded in modulus by u, the composition of n operations leads to a worst-case error bound with dominant term proportional to nu. But this worst-case error bound is attained only when each rounding error is of maximal magnitude and identical sign, which is very unlikely. Since the beginning of the digital computer era many researchers have modelled rounding errors as random variables in an attempt to obtain better estimates of how the error behaves on average. This line of thinking has led to the well-known rule of thumb, based on informal arguments and assumptions, that constants in rounding error bounds can be replaced by their square roots.

In our EPrint A New Approach to Probabilistic Rounding Error Analysis we make minimal probabilistic assumptions on the rounding errors and make use of a tool from probability theory called a concentration inequality. We show that in several core linear algebra algorithms, including matrix-vector products and LU factorization, the backward error can be bounded with high probability by a relaxed constant proportional to \sqrt{n\log n}u instead of nu. Our analysis provides the first rigorous justification of the rule of thumb.

This new bound is illustrated in the figure above, where we consider the solution of a linear system Ax = b by LU factorization. The matrix A and vector x have entries from the random uniform [0,1] distribution and b is formed as Ax. We compare the backward error with its traditional worst-case bound and our relaxed probabilistic bound. The figure shows that the probabilistic bound is in very good agreement with the actual backward error and is much smaller than the traditional bound. Moreover, it successfully captures the asymptotic behavior of the error growth, which follows \sqrt{n} rather than n.

The assumptions underlying our analysis—that the rounding errors are independent random variables of mean zero—do not always hold, as we illustrate with examples in the paper. Nevertheless, our experiments show that the bounds do correctly predict the backward error for a selection of real-life matrices from the SuiteSparse collection.

Fast Solution of Linear Systems via GPU Tensor Cores’ FP16 Arithmetic and Iterative Refinement

Posted on by Leave a comment
181112_sc18.jpg

NVIDIA Founder & CEO Jensen Huang talking about the work reported here in his special address at Supercomputing 2018 (8:30 onwards).

Over the last 30 years, hierarchical computer memories, multicore processors and graphical processing units (GPUs) have all necessitated the redesign of numerical linear algebra algorithms, and in doing so have led to algorithmic innovations. Mixed precision arithmetic—a concept going back to the earliest computers, which had the ability to accumulate inner products in extra precision—attracted renewed interest in the late 1990s once Intel chips were able to execute single precision at twice the rate of double precision. Now the increasing availability of low precision arithmetic is offering new opportunities.

In the paper Harnessing GPU Tensor Cores for Fast FP16 Arithmetic to Speed up Mixed-Precision Iterative Refinement Solvers presented at SC18 (the leading supercomputing conference), Azzam Haidar, Stanimire Tomov, Jack Dongarra and Nick Higham show how to exploit the half precision (fp16) arithmetic that is now available in hardware. Whereas fp16 arithmetic can be expected to run at twice the rate of fp32 (single precision) arithmetic, the NVIDIA V100 GPU has tensor cores that can execute half precision at up to eight times the speed of single precision and can deliver the results to single precision accuracy. Developing algorithms that can exploit half precision arithmetic is important both for a workstation connected to a single V100 GPU and the world’s fastest computer (as of November 2018): Summit at Oak Ridge National Laboratory, which contains 27,648 V100 GPUs.

The paper shows that a dense n-by-n double precision linear system Ax = b can be solved using mixed precision iterative refinement at a rate up to four times faster than a highly optimized double precision solver and with a reduction in energy consumption by a factor five.

The key idea is to LU factorize the matrix A in a mix of half precision and single precision then apply iterative refinement. The update equations in the refinement process are solved by an inner GMRES iteration that uses the LU factors as preconditioners. This GMRES-IR algorithm was proposed by Carson and Higham in two (open access) papers in SIAM J. Sci. Comput. (2017 and 2018). In the form used here, the algorithm converges for matrices with condition numbers up to about 10^8. It provides a backward stable, double precision solution while carrying out almost all the flops at lower precision.

Codes implementing this work will be released through the open-source MAGMA library.

htdh18_fig7a.jpg

Lecturer, Senior Lecturer or Reader in Applied Mathematics

Posted on by Leave a comment

The School of Mathematics is seeking mathematical scientists of outstanding ability or potential for appointments at Lecturer, Senior Lecturer or Reader level.

Applicants working in any area of applied mathematics are welcome, particularly those working in areas that complement and enhance the applied mathematics research of the School, which spans numerical linear algebra, uncertainty quantification, dynamical systems, data science, mathematics in the life sciences, industrial mathematics, inverse problems and continuum mechanics.

Applicants with research experience in numerical linear algebra, or at the interfaces with pure mathematics or probability and statistics, are strongly encouraged.

These positions provide an excellent opportunity to join the Numerical Linear Algebra Group.

The closing date is January 11, 2019. For the advert and more details see here.

« Older Entries Recent Entries »