This paper presents a new algorithm for performing model order reduction for Volterra series channel models of high-density digital magnetic recording channels. We employ a set-membership approach to the problem in which a set of consistent modeling solutions bounded by an optimal ellipsoid is first developed for the channel. We then present a new algorithm for finding the minimum number of coefficients in the Volterra series expansion which preserves the accuracy, in a least-squares sense, of the reduced order model in comparison with the original model.