Computing the solution of linear systems of equations is invariably the most time consuming task in the numerical solutions of PDEs in many fields of computational science. In this study, we focus on the numerical simulation of cardiovascular hemodynamics with rigid and deformable walls, discretized in space and time through the variational multiscale finite element method. We focus on three approaches: the problem agnostic generalized minimum residual (GMRES) and stabilized bi-conjugate gradient (BICGS) methods, and a recently proposed, problem specific, bi-partitioned (BIPN) method. We also perform a comparative analysis of several preconditioners, including diagonal, block-diagonal, incomplete factorization, multigrid, and resistance based methods. Solver performance and matrix characteristics (diagonal dominance, symmetry, sparsity, bandwidth and spectral properties) are first examined for an idealized cylindrical geometry with physiologic boundary conditions and then successively tested on several patient-specific anatomies representative of realistic cardiovascular simulation problems. Incomplete factorization preconditioners provide the best performance and results in terms of both strong and weak scalability. The BIPN method was found to outperform other methods in patient-specific models with rigid walls. In models with deformable walls, BIPN was outperformed by BICG with diagonal and Incomplete LU preconditioners.