Submerged impingement jets are widely used in erosion/corrosion experiments as it is easy to control jet standoff distance, jet angle, and flow velocities in experiments. In addition to experiments, typically computational fluid dynamics (CFD) technique has been used to simulate slurry flow in this geometry to investigate erosion process and develop erosion models or equations. The traditional CFD simulations of erosion in this geometry use the Reynolds-averaged Navier–Stokes (RANS) equations with turbulence models. By using this technique, time-averaged fluid flow is revealed, and thus, time-averaged erosion rate can be obtained by tracking particles in the fluid flow field. However, this seemingly simple flow displays unsteady flow structures in the stagnation zone of the flow field and its effects on the erosion process were previously unclear. In this study, large eddy simulation (LES) is used to simulate unsteady fluid flow in different impingement jets in an Eulerian scheme. Then, transient particle tracking is performed in a Lagrangian scheme. Particles are injected randomly at the inlet plane and tracked to simulate unsteady erosion that occurs on the target surface. Finally, an erosion equation is used to calculate solid particle erosion rates. The LES Eulerian–Lagrangian erosion modeling is further validated by available experimental data for fluid velocities and an erosion profile. The results show that the accuracy of erosion prediction of small particles is improved significantly by using the LES method. In addition, the unsteady particle motion and erosion process can be revealed by using this method.