Most of existing methods for solving the discrete-velocity Bhatnagar-Gross-Krook (BGK) model of the Boltzmann equation are explicit or semi-implicit, both have stability constraints on the time step size that leads to a large number of time step for steady state problems. In this talk, a fully implicit scheme will be discussed, together with a nonlinearly preconditioned algorithm for the solution of the large sparse nonlinear system of equations arising at each time step. We show by numerical experiments that the new method is robust for high Reynolds number flows, has no CFL constraints, and is scalable on a supercomputer with thousands of processors. This is a joint work with J. Huang and C. Yang.