A new, computationally efficient approach for modelling finite amplitude ultrasound propagation is described. The model is able to simulate nonlinear distortion of CW and pulsed excitations from non-axisymmetric sources in tissue. We have developed a second order operator splitting approach, enabling the effects of diffraction, nonlinearity, and absorption to be calculated separately over relatively large incremental distances using a fractional step marching scheme. A computationally efficient angular spectrum algorithm has also been developed to calculate the diffractive propagation from non-axisymmetric, non-separable sources. Results of our model have shown close agreement with published data. Moreover, our approach may offer computational savings compared with existing models. Indeed, with our algorithm it should be possible to simulate the nonlinear propagation of sound beams from realistic medical ultrasound scanners, and perhaps to investigate ways to improve the design of tissue harmonic imaging systems.