For simulations involving highly complex geometries, as they occur in many fields of science and engineering, the process of generating a high-quality grid is extremely time-consuming. Specifically, for flows containing moving boundaries, the considerable advantage of the immersed interface/boundary method becomes evident. While good results, both qualitative and quantitative, have been obtained, most current schemes rely on low-order corrections in the vicinity of boundary/interface. The objective of the paper is to present a newly developed immersed interface method for solving the incompressible Navier-Stokes equations with moving boundaries in primitive variable formulation on non-staggered grids. This time-explicit immersed interface method is based on a local Taylor series expansion at irregular grid points whereby numerical stability is enforced through a local stability condition. A matrix stability analysis of the full discretization matrix was used to rigorously analyze the impact of different boundary discretizations on the overall stability of the numerical scheme. To validate the immersed incompressible Navier-Stokes solver, two different physically relevant flows are simulated. The two-dimensional flow around a circular cylinder placed in a uniform free-stream is simulated in the range ReD = 20 to ReD = 200 covering the range of steady and unsteady flow regimes. Computations of Tollmien-Schlichting waves in a two-dimensional channel flow are presented for a sub-critical and super-critical case, challenging the near wall accuracy of current schemes. Finally, an analytical test function approach demonstrates the full capability of the present immersed Navier-Stokes solver with moving boundaries.