A Cartesian grid method has been developed for simulating two-dimensional unsteady, viscous, incompressible flows with complex immersed boundaries. A finitevolume method based on a second-order accurate central-difference scheme is used in conjunction with a two-step fractional-step procedure. The key aspects that need to be considered in developing such a solver are imposition of boundary conditions on the immersed boundaries and accurate discretization of the governing equation in cells that are cut by these boundaries. A new interpolation procedure is presented which allows systematic development of a spatial discretization scheme that preserves the second-order spatial accuracy of the underlying solver. The presence of immersed boundaries alters the conditioning of the linear operators and this can slow down the iterative solution of these equations. The convergence is accelerated by using a preconditioned conjugate gradient method where the preconditioner takes advantage of the structured nature of the underlying mesh. The accuracy and fidelity of the solver is validated by simulating a number of canonical flows and the ability of the solver to simulate flows with very complicated immersed boundaries is demonstrated.