  In this paper, we propose a parallel-in-time algorithm for approximately
solving parabolic equations. In particular, we apply the $k$-step backward
differentiation formula, and then develop an iterative solver by using the
waveform relaxation technique. Each resulting iteration represents a
periodic-like system, which could be further solved in parallel by using the
diagonalization technique. The convergence of the waveform relaxation iteration
is theoretically examined by using the generating function method. The approach
we established in this paper extends the existing argument of single-step
methods in Gander and Wu [Numer. Math., 143 (2019), pp. 489--527] to general
BDF methods up to order six. The argument could be further applied to the
time-fractional subdiffusion equation, whose discretization shares common
properties of the standard BDF methods, because of the nonlocality of the
fractional differential operator. Illustrative numerical results are presented
to complement the theoretical analysis.
