SUMMARY The ‘pQononic lattice solid’ (PLS) approach based on simulating the movement, interaction and scattering of quasi-particles on a discrete lattice has been developed recently to model wave propagation through heterogeneous media. In the initial version, the quasi-particles carried pressure and the differential equation describing the transport of quasi-particle number densities on the discrete lattice was solved using a finite-difference scheme. The significant problems and the limitations of this method are that only compressional waves are considered, the convergence is slow in homogeneous regions, the results are inaccurate when impedance contrasts are large and small lattice spacings are required to specify sharp interfaces (cf the classical finite-difference solution to the wave equation). To address the last three points, we develop an improved PLS approach by interpolation (PLSI) to directly model the behaviour of the quasi-particle number densities on the discrete lattice rather than solving the corresponding macroscopic transport equation using finite differences. This involves simulating the three microscopic processes describing the behaviour of quasi-particles: movement along the links between lattice nodes, scattering by medium heterogeneities and interaction between quasi-particles arriving at lattice nodes. In the movement step, quasi-particle number densities are moved from the nodes along the links by an amount cAr where c is the quasi-particle speed and Ar is the time step as if there were no heterogeneity in the link. Number densities are then interpolated to the locations of interfaces between parts of links with different properties and scattering is taken into account using the known 1-D reflection and transmission coefficients. Theoretical analysis demonstrates that, in the macroscopic limit, the PLSI models the acoustic wave equation for N-D heterogeneous media (N = 1,2, 3). A 2-D numerical example illustrates the PLSI approach. The PLSI is comparable with the lattice Boltzmann lattice gas appproach in that no finite-difference errors are present but models wave phenomena in heterogeneous media rather than fluid flow and acoustic waves in a constant bulk modulus gas. Because the PLSI can handle sharp interfaces at any location (cf classical finite difference methods have difficulty handling sharp interfaces), it may enable numerical experiments to be conducted of wave propagation through complex fractured and porous rocks to study the causes of anisotropy and attenuation including the effect of non-linear solid-fluid interactions. This would require the extension of the approach to model shear waves in addition to compressional waves.