Abstract The non-Newtonian fluid flow with a free surface occurring during a plane channel filling in the gravity field has been simulated numerically. The mathematical statement of the problem is formulated on the basis of the motion equations, continuity equation, and natural boundary conditions on the free surface with an application of the Herschel-Bulkley rheological model. A traditional mathematical model singularity on the three-phase contact line is eliminated using a slip condition. A numerical algorithm based on the finite-difference method is developed for solving the problem. Regularization of the rheological equation has been carried out using a shock-capturing method for the flow with unyielded regions. A parametric investigation of the flow kinematics and free surface behavior in terms of the governing parameters has been implemented. The flow structures distinguished by the presence of unyielded regions have been demonstrated depending on relation of the viscous, gravity, and plastic forces in the flow.