The injection molding process is used in the fabrication of a large variety of plastic articles. A numerical simulation of the filling stage along the thickness direction is proposed by combining the free surface boundary condition with the relevant governing equations. The mathematical model is based on the equations of continuity, momentum and energy along with inelastic power-law model and relevant boundary conditions. Due to the significant implications for microstructure development in the pro duct, the fountain effect at the advancing free surface is explicitly taken into consideration in the simulation. The model yields data on free surface shape as well as velocity, pressure, temperature and shear stress distributions within the mold cavity. The rearrangement of the velocity and temperature profiles in the vicinity of the melt front is considered in detail.