Ever increasing computational resources make high-resolution mesoscale simulations (~1 km or less grid-cell size) possible for a wide range of applications including: transport and dispersion, renewable energy forecasting, wildland fire prediction, etc. However, reduced grid-cell sizes and consequent violation of the the assumption of horizontal homogeneity require extending turbulence parameterizations to three spatial dimensions. We have therefore developed a three-dimensional planetary boundary layer (3D PBL) parameterization based on the work of Mellor and Yamada (1982). To determine the scales at which the three-dimensional parameterization would be needed we carried out a nested, high-resolution, turbulence resolving, large-eddy simulation (LES) of flow over complex terrain based on the Wind Forecast Improvement 2 (WFIP2) field study. The WFIP2 field study took place in the Columbia River Gorge area during 2016 and early 2017. We focus on selected cases when physical phenomena of significance for forecasting in complex terrain such as mountain waves, topographic wakes and gap flows were observed. LES results are first validated using observations from the field study and then used to compute turbulent stresses and fluxes as well as their divergence at different spatial scales. Analysis of turbulence quantities computed at different scales provides an estimate of a scale below which 3D PBL parameterization should be used. Furthermore, together with filed study observations, validated LES results are used to assess the new 3D PBL parameterization.