Richards' equation (RE) is often used to model flow in unsaturated porous media. This model captures physical effects, such as sharp fronts in fluid pressures and saturations, which are present in more complex models of multiphase flow. The numerical solution of RE is difficult not only because of these physical effects but also because of the mathematical problems that arise in dealing with the nonlinearities. The method of lines has been shown to be very effective for solving RE in one space dimension. When solving RE in two space dimensions, the nonlinear equations that must be solved in an implicit time-stepping approach usually become impractical to solve using a direct linear solver (e.g., LU decomposition). In this work, we show how Newton-iterative methods can be applied in this context; these methods use iterative linear solvers, such as the Krylov methods GMRES and Bi-CGSTAB to solve the linear systems associated with Newton's method. We present a theorem on the convergence of inexact Newton methods, and based on this theorem, an adaptive method that selects an appropriate linear tolerance. Numerical results show the method to be effective compared with an existing approach.