\|_{L^\infty(0,T; H^{-1}(\Omega))} \le C \left (\Delta t + \left (\ln \left (\frac{1}{h}\right )\right )^{\frac{2m+3}{2m+2}} h\right )
\label{rulla-1}$$ by Rulla and Walkington [@Rulla1996], $$\| u_h - u \|_{L^2(0,T; L^2(\Omega))} \le C h^{\frac{m^2+6m+8}{6m^2+14m+8}}
\quad \text{when}\quad \Delta t = \mathcal{O}(h^{\frac{5 m+4}{2m}})
\label{ebmeyer-1}$$ by Ebmeyer [@Ebmeyer-1998], and $$\| u_h - u \|_{L^{m+2}(0,T; L^{m+2}(\Omega))} \le
C \left ( \Delta t^{\frac 1 2} + h + h^{\frac{1}{m+1} \left (\frac{d m}{2m+4}+1\right )} \right )^{\frac{1}{m+2}}
\label{wei-1}$$ by Wei and Lefton [@Wei1999], where $d$ is the space dimension. It is remarked that these estimates are obtained for quasi-uniform meshes. The convergence rate is first order at best and decreases with $m$. Some of these estimates are shown to be optimal in the corresponding norm in lieu of the known regularity of the solution of PME. Moreover, it is worth mentioning that Ebmeyer and Liu [@Ebmeyer-2008] obtain error estimates in quasi-norm and Emmrich and [Š]{}i[š]{}ka [@Emmrich2012] prove that a Galerkin finite element – backward Euler approximation converges to the weak solution of PME using the theory of monotone operators. More recently, Duque et al. [@Duque2013] establish $L^{1+\max (\gamma/2)}$ error bounds for the approximation of a general order continuous Galerkin in space and a general order discontinuous Galerkin in time for 