The weak maximum principle of finite element methods for parabolic equations is proved for both semi-discretization in space and fully discrete methods with k-step backward differentiation formulae for \(k = 1,\ldots ,6\) , on a two-dimensional general polygonal domain or a three-dimensional convex polyhedral domain. The semi-discrete result is established via a dyadic decomposition argument and local energy estimates in which the nonsmoothness of the domain can be handled. The fully discrete result for multistep backward differentiation formulae is proved by utilizing the solution representation via the discrete Laplace transform and the resolvent estimates, which are inspired by the analysis of convolutional quadrature for parabolic and fractional-order partial differential equations.
