The present study models the multi-material topology optimization problems as the multi-valued integer programming (MVIP) or named as combinatorial optimization. By extending classical convex analysis and convex programming to discrete point-set functions, the discrete convex analysis and discrete steepest descent (DSD) algorithm are introduced. To overcome combinatorial complexity of the DSD algorithm, we employ the sequential approximate integer programming (SAIP) to explicitly and linearly approximate the implicit objective and constraint functions. Considering the multiple potential changed directions for multi-valued design variables, the random discrete steepest descent (RDSD) algorithm is proposed, where a random strategy is implemented to select a definitive direction of change. To analytically calculate multi-material discrete variable sensitivities, topological derivatives with material contrast is applied. In all, the MVIP is finally transferred as the linear 0–1 programming that can be efficiently solved by the canonical relaxation algorithm (CRA). Explicit nonlinear examples demonstrate that the RDSD algorithm owns nearly three orders of magnitude improvement compared with the commercial software (GUROBI). The proposed approach, without using any continuous variable relaxation and interpolation penalization schemes, successfully solves the minimum compliance problem, strength-related problem, and frequency-related optimization problems. Given the algorithm efficiency, mathematical generality and merits over other algorithms, the proposed RDSD algorithm is meaningful for other structural and topology optimization problems involving multi-valued discrete design variables.