On the accuracy of classic numerical schemes for modeling flow in saturated heterogeneous formations