The supercooled large droplet (SLD) environments are sharply different from those of traditional supercooled droplets with small sizes, lying in the fact that the SLD environment has a wider range of droplet sizes and complex droplet size distribution forms, which makes ice shapes complex and has brought an unprecedented challenge to airline safety. The existing numerical icing simulation methods can only simulate traces of single-size droplets, lacking the ability to simulate traces of the bimodal droplet size distributions and the corresponding ice shapes in the SLD environment precisely. In order to simulate the characteristics of the SLD droplet size distributions precisely and efficiently, a mixed Eulerian-Lagrangian sampling method based on Rosin-Rammler distribution function was proposed. Water droplets were collected via this method, and heat transfer and ice accretion were calculated with the icing model and the water film model, through which a precise and efficient numerical icing simulation of the SLD droplet size distributions was realized. Furthermore, research was done via 2.5D cases on the impact of the SLD droplet size distributions on ice shape characteristics and ice shape spatial randomness. The results show that droplet size distributions do affect ice shapes, and ice shape characteristics and ice shape randomness are closely related to the MVD and the distribution variance.