- Porous media are the materials containing void space or pores, where fluids and gases can pass through. Unsteady flows in porous media are encountered in many engineering applications and natural problems, such as CO₂ sequestration, high temperature nuclear reactor cooling, high efficiency combustion, chemical reactors, noise reduction on airplane trailing edges, dispersion of pollutants in the hyporheic zones, groundwater remediation, etc. In this work, unsteady inertial and turbulent flows are studied numerically in a triply-periodic unit cell of a face-centered cubic (FCC) lattice. Pore-scale direct numerical simulations (DNS) are performed using a fictitious domain method (FDM). Flows at three different pore Reynolds numbers (Rep = 300, 500, and 1000) are simulated, ranging from unsteady inertial to turbulent flow regimes. First, the flow physics are studied mostly in the Eulerian frame. The characteristics of the flow topology, turbulent kinetic energy (TKE) transport, and anisotropy distribution at different Reynolds numbers are investigated. The mean flow through the unit cell is driven by an axial pressure gradient. Negative TKE production is observed in the flow domain near bead surfaces, attributed to the tortuous flow paths. The turbulence length scales are quantified by the two-point auto correlations to be less than 10% of the bead diameter. Next, the flow field is analyzed using angular Lagrangian multiscale statistics of fluid particles to study their directional change at different time scales. Two power laws are observed for the evolution of the mean absolute angle as a function of time lag demarking the early time and an intermediate, inertial range. The effect of the geometric confinement on the asymptotic behavior of the angular statistics is examined in detail. An asymptotic limit different than that for three-dimensional isotropic turbulence is obtained. Furthermore, a Monte Carlo-based stochastic model is developed to predict the asymptotic curvature angle. Lastly, the method of volume averaging is used to obtain the upscaled continuum balance equations for turbulent flows through porous media. A standard closure model based on a modified Ergun equation is applied and evaluated for turbulent flows using the DNS data. However, the empirically obtained coefficients of the modified Ergun correlation tend to overpredict the pressure drop in the fully turbulent flow regime, suggesting the need for a more accurate model. To overcome such drawback, we propose an algebraic model by upscaling the DNS data. Finally, DNS at pore Reynolds number of 300 is performed at a 3x3x3 matrix of the FCC unit cell to verify the algebraic model.