The shifted Legendre algorithm is used to solve the governing differential equation of variable order cantilever beam in time domain. According to the variable fractional model, the governing equations of viscoelastic polymer cantilever beam are built. The integer order and variable fractional order differential operator matrices are obtained based on the properties of shifted Legendre polynomials. The variable fractional differential equations are converted into algebraic equations and solved in time domain by operator matrix. The convergence analysis confirmed the high efficiency and accuracy of the proposed algorithm for solving differential equations. The deflection and stress of the polymer cantilever beams (PET and HDPE) under various external load (uniformly distributed and harmonic) are studied.