In this paper, we propose a weakly compressible model for numerical simulation of the incompressible viscous flows. This model asymptotically approximates the incompressible Navier-Stokes equations when the mach number tends to zero. The main advantage of this model is that its numerical discretization avoids any Poisson solver, thus is very attractive for problems with complicated geometries. This model is discretized by high-order center differences in space and the so-called 'I-stable' method for time. The linear stability region of an I-stable method contains part of the imaginary axis. When solving a system of convection-diffusion equations with a small viscosity, the I-stable method allows a very large cell Reynolds number, thus is particularly suitable for the simulation of fluid flows with large Reynolds numbers. Numerical experiments illustrate the efficiency and robustness of this approach. (C) 2001 Elsevier Science B.V. All rights reserved.