  We present an algorithm for the solution of a simultaneous space-time
discretization of linear parabolic evolution equations with a symmetric
differential operator in space. Building on earlier work, we recast this
discretization into a Schur-complement equation whose solution is a
quasi-optimal approximation to the weak solution of the equation at hand.
Choosing a tensor-product discretization, we arrive at a remarkably simple
linear system. Using wavelets in time and standard finite elements in space, we
solve the resulting system in linear complexity on a single processor, and in
polylogarithmic complexity when parallelized in both space and time. We
complement these theoretical findings with large-scale parallel computations
showing the effectiveness of the method.
