Turbulent flow through porous media is of considerable theoretical and practical interest in various engineering disciplines such as chemical, mechanical, nuclear, geological, environmental, petroleum, etc. In the present work, a numerical study of turbulent flow through porous media was carried out. The porous media was considered to be made of a spatially periodic array of infinite square cylinders. The low Re–k–ε–Lam–Bremhorst (LB) turbulent model was employed for the description of the turbulent flow in the porous media. The transport equations were solved only through the flow domains having different porosities. The porosity of the medium (ϕ) and the Reynolds number (Re D), were taken in the range of 0.3⩽ ϕ⩽ 0.84, and 100⩽ Re D⩽ 40,000, respectively. The influence of ϕ and Re D on normalized turbulent kinetic energy,< k> f and dissipation rate,< ε> f was investigated and the simulated data were compared with two different turbulent flow models, namely, large eddy simulation (LES), and v2f model. The predictions of turbulent parameters using low Re k–ε–LB turbulent model were found to be in good agreement with the LES than that with v2f model. The macroscopic pressure variation across the flow domain was also investigated and the predictions of the pressure gradient were found to be in good agreement with the Forchheimer-extended Darcy law. At low Re D (Re D= 100) the low Re–k–ε–LB predictions were deviating from the DNS, LES and v2f predictions. This means that the proposed low Re–k–ε–LB model will hold valid only for flows having Re D> 300. A correlation for macroscopic pressure gradient as a function of porosity is also proposed for high Re D.