In high precision laser ranging, simulation of intersatellite ultra-long distance propagation plays a crucial role. The lasers’ propagation distance between satellites takes hundreds to millions of kilometers; due to wavefront aberration and the diffraction effect in the detection approach, the lasers’ propagation cannot be directly simulated using the classical Gaussian beam tracing approach. A mode expansion method is used in this study to simulate the ultra-long distance intersatellite diffraction propagation. The simulation accuracy is confirmed by comparing with the analytical approaches and Fraunhofer diffraction integral. The mode order influence and ratio of the initial Gaussian beam waist radius to the aperture radius on the simulation accuracy are examined. Thus, this numerical simulation approach with both computation accuracy and efficiency is offered for the high precision simulation laser ranging in space.