The Integral Transform Method is employed in the hybrid numerical‐analytical solution of laminar natural convection inside rectangular enclosures. The hybrid nature of this approach allows for the establishment of benchmark results in the solution of non‐linear partial differential systems, including the coupled set of heat and fluid flow equations that govern the steady natural convection problem under consideration. Therefore, numerical results under user prescribed accuracy are obtained for different values of the governing parameters, and the convergence behaviour of the proposed eigenfunction expansions is illustrated. Critical comparisons against previously reported benchmark solutions are then performed. Copyright © 1999 John Wiley & Sons, Ltd.