We revisit an approximate stochastic dynamic programming method that we proposed earlier for the optimization of multireservoir problems. The method exploits the convexity properties of the value function to sample the reservoir level space based on the local curvature of the value function, which is estimated by the difference between a lower and an upper bounds (error bound). Unlike the previous approach where the state space was exhaustively partitioned into full dimensional simplices whose vertices formed a discrete grid over which the value function was approximated, here we propose instead a new randomized approach for selecting the grid points from a small number of randomly sampled simplices from which an error bound is estimated. Results of numerical experiments on three literature test problems and simulated midterm reservoir optimization problems illustrate the advantages of the randomized approach which can solve models of higher dimensions than with the exhausitive approach.