A new coupling method, immersed moving boundary–discrete element method–lattice Boltzmann method (IMB-DEM-LBM), is proposed to simulate particle flow for application in soil mechanics or coastal engineering. In this study, LBM fluid is simulated on the regular Eulerian grid and Lagrangian particle motion is governed by DEM while IMB couples the two algorithms. The new method is promising and robust as it resolves numerical instability near the particle boundary caused by mesh distortion in the conventional grid method. In IMB, the interface lattice solid fraction determines the distribution function ratio of non-equilibrium bounce back and Bhatnagar-Gross-Krook (BGK) collision. The non-equilibrium bounce back at moving boundary results in the fluid momentum change and contributes to the hydrodynamic force on particle. For numerical stability, this paper introduces the hydrodynamic force calculation concept from IB (immersed boundary method) to IMB, and at the same time, proposes a new solid fraction calculation method for sphere that divides the intersection into simple sector and triangle, as well as calculates the intersection area by vector. With this method, approximate inaccuracy is overcome while complicated integration is avoided.