We introduce a novel variational framework for the regularized reconstruction of time-resolved volumetric flow fields. Our objective functional takes the physical characteristics of the underlying flow into account in both the spatial and the temporal domains. For an efficient minimization of the objective functional, we apply a proximal-splitting algorithm and perform parallel computations. To demonstrate the utility of our variational method, we first denoise a simulated flow-field in the human aorta and show that our method outperforms spatial-only regularization in terms of signal-to-noise ratio (SNR). We then apply the scheme to a real 3D+time phase-contrast MRI dataset and obtain high-quality visualizations.