In this paper, we develop a high-order time-stepping scheme for solving the nonlinear diffusion-wave equation with the Caputo fractional time derivative of order \(\alpha ~(1<\alpha <2)\) . Based on smooth initial data and assumptions regarding the nonlinear source term, we discuss the regularity estimates of the solution. Linear interpolation and quadratic interpolation are employed to approximate the integral in the representation of the solution. Then a multistep Mittag–Leffler integrator is proposed for solving the nonlinear diffusion-wave equation. Additionally, two methods are proposed for computing the Mittag–Leffler function in time discretization. The error estimates of the numerical scheme are rigorously proven. The convergence order is \(O(\tau ^{\alpha +1})\) , where \(\tau \) denotes the temporal grid size. Numerical examples are presented to demonstrate the theoretical results of the proposed scheme.
