The plan is to divide the talk in three distinct but related parts. First, the question of asymptotic stability for equilibria of delay differential equations is addressed numerically. The proposed method, based on the discretization of the infinitesimal generator of the solution operator semigroup via pseudospectral differentiation, allows to approximate the stability determining eigenvalue with spectral accuracy. Hence it is fast and suitable for robust analysis. Second, the numerical scheme is extended for investigating the stability of steady states of population dynamics, where the study of the associated transcendental characteristic equations is often too difficult to be approached analytically. The fruitful interplay between theoretical and numerical analysis is highlighted through examples taken from age- and physiologically-structured models, as well as delayed epidemics. Third, recent advances in the numerical stability analysis of delay systems are illustrated, showing how equilibria (characteristic roots), periodic orbits (Floquet multipliers) and chaotic motion (Lyapunov exponents) can be faced under the same discretization framework. Examples arising in the populations context are discussed which demand for adapting such treatment.