This paper derives a new operational matrix of the variable-order (VO) time fractional partial derivative involved in anomalous diffusion for shifted Chebyshev polynomials. We then develop an accurate numerical algorithm to solve the 1 + 1 and 2 + 1 VO and constant-order fractional diffusion equation with Dirichlet conditions. The contraction of the present method is based on shifted Chebyshev collocation procedure in combination with the derived shifted Chebyshev operational matrix. The main advantage of the proposed method is to investigate a global approximation for spatial and temporal discretizations, and it reduces such problems to those of solving a system of algebraic equations, which greatly simplifies the solution process. In addition, we analyze the convergence of the present method graphically. Finally, comparisons between the algorithm derived in this paper and the existing algorithms are given, which show that our numerical schemes exhibit better performances than the existing ones.