Based on the principle that the volume of solid particles carried out of or into the stratum by muddy water equals the change in pore volume, a formula for the variation of porosity (void ratio) with time during the seepage process of muddy water is derived. Based on Darcy’s permeability coefficient formula, the relationship between the permeability coefficient and the initial porosity (initial void ratio) as well as the initial permeability coefficient is established, and a formula for the temporal variation of permeability coefficient during the muddy water seepage process is developed. According to the critical initiation condition for piping, where the seepage drag force acting on free soil particles is balanced by their effective weight, formulas for calculating the critical seepage velocity and critical hydraulic gradient for piping initiation are derived, and a numerical method based on finite element software is proposed to iteratively solve the initiation and evolution process of unsteady seepage piping. The proposed method is validated using a piping failure case in the filter layer of a medium-sized reservoir dam in Shandong Province. The results show that the method can accurately reproduce the evolution of piping failure in the filter layer. The permeability coefficient of the filter layer increases sharply within a short period after particle loss is completed, leading to a rapid increase in hydraulic gradient and triggering a chain failure, eventually resulting in the breakthrough of the clay inclined wall and the formation of a concentrated seepage channel. The calculated piping failure mode agrees well with the measured locations and development patterns of collapse pits in the field.