The influence of two-phase flows containing suspension particles, which are common in nature, on internal erosion with coupling effect of clogging remains unclear. This paper presents a three-dimensional coupled discrete element method and computational fluid dynamics (CFD-DEM) analysis of internal erosion considering different concentrations of suspension C (i.e., mass of the suspended particles in unit volume of fluid) in gap-graded granular soils with different fine fraction Fc (i.e., the percentage by mass of the fine particles in the gap-graded sample). The influences of C and Fc on the erosion and clogging behavior of soils are investigated from both the macroscopic and microscopic perspectives. It is found that for gap-graded samples being under-filled with Fc=15%, the suspension flow (i.e., influent fluid with suspending particles) decreases the cumulative eroded fine particle loss and the increasing rate of soil hydraulic conductivity due to clogging at the top of the sample. The degree of clogging is found to jointly be determined by both constriction size distribution and the suspension concentration. Clogging in a local area usually occurs with the formation of the clusters which has a high resistance to the drag force applied by the fluid flow.