The role that the classical special functions, such as the Airy, Bessel, Hermite, Legendre and hypergeometric functions, started to play in the 19th century, has now been greatly expanded by the Painlevé functions. Increasingly, as nonlinear science develops, people are finding that the solutions to an extraordinarily broad array of scientific problems, from neutron scattering theory, to partial differential equations, to transportation problems, to combinatorics, etc., can be expressed in terms of Painlevé transcendents. Much can be, and has been, proved regarding the algebraic and asymptotic properties of Painlevé transcendents. Here the role of integral representations and the classical steepest descent method in deriving precise asymptotics and connection formulae for the classical special functions is played, and expanded, by a Riemann-Hilbert representation of the Painlevé equations. The Riemann-Hilbert method is based on the observation that the Painlevé equations describe the isomonodromy deformations of certain systems of linear differential equations with rational coefficients, so solving a Painlevé equation is equivalent to solving an inverse monodromy problem. However on the other hand, very little is known, beyond some ad hoc calculations, about the numerical solution of the Painlevé equations. Writing useful software for nonlinear equations such as the Painlevé equations presents many challenges, conceptual, philosophical and technical. Without the help of linearity, it is not at all clear how to select a broad enough class of “representative problems”.