Numerical models are frequently used in order to quantify the effect of the spatial heterogeneity in the formation hydraulic properties on field-scale solute transport. Because of practical constraints it is common to discretize the flow domain into relatively large numerical cells, with characteristic length scales, ℓ1, ℓ2, ℓ3, comparable with those of the formation heterogeneity, I1, I2, I3, leading to a loss in the velocity variability, and, concurrently, in solute spreading. We analyzed block-effective dispersion coefficients, D'ij (i, j = 1, 2, 3), required to compensate for this loss in three-dimensional, variably saturated heterogeneous formations under conditions of steady state, gravity-dominated, unsaturated flow. The results that the principal components of D'ij are controlled by the ratio ℓ' = ℓ2/I2 = ℓ3/I3, that they may reach their ergodic limits if ℓ' is sufficiently large, and that their tendency to their asymptotic limits with increasing scaled travel time, t', is inversely related to ℓ', are in agreement with previous results for saturated flow. New findings suggest that under unsaturated flow conditions, both the time required to reach the asymptotic limits of D'ij and the length scale needed to reach the ergodic limits of D'ij depend on soil stratification, soil texture, and mean soil water saturation and differ for the longitudinal and the transverse components of D'ij.