Summary: Blood flow in arterial systems is described by the three-dimensional Navier-Stokes equations within a time-dependent spatial domain that accounts for the viscoelasticity of the arterial walls. These equations are simplified by assuming cylindrical geometry, axially symmetric flow, and hydrostatic equilibrium in the radial direction. An efficient semi-implicit method is formulated in such a fashion that numerical stability is obtained at a minimal computational cost. The resulting computer model is relatively simple, robust, accurate, and extremely efficient. These features are illustrated on nontrivial test cases where the exact analytical solution is known and by an example of a realistic flow through a complex arterial system.