  For linear parabolic initial-boundary value problems with self-adjoint,
time-homogeneous elliptic spatial operator in divergence form with
Lipschitz-continuous coefficients, and for incompatible, time-analytic forcing
term in polygonal/polyhedral domains $D$, we prove time-analyticity of
solutions. Temporal analyticity is quantified in terms of weighted, analytic
function classes, for data with finite, low spatial regularity and without
boundary compatibility. Leveraging this result, we prove exponential
convergence of a conforming, semi-discrete $hp$-time-stepping approach. We
combine this semi-discretization in time with first-order, so-called
"$h$-version" Lagrangian Finite Elements with corner-refinements in space into
a tensor-product, conforming discretization of a space-time formulation. We
prove that, under appropriate corner- and corner-edge mesh-refinement of $D$,
error vs. number of degrees of freedom in space-time behaves essentially (up to
logarithmic terms), to what standard FEM provide for one elliptic boundary
value problem solve in $D$. We focus on two-dimensional spatial domains and
comment on the one- and the three-dimensional case.
