An improved model and algorithm were developed for the numerical solution of three-dimensional, laminar natural convection in enclosures. Illustrative calculations were carried out for heating from below in a cubical box, a long channel with a square cross section and in the region between infinite, horizontal plates. The results are in good agreement with prior experimental and theoretical results. For the infinite flat plates the three-dimensional model correctly produces a two-dimensional solution. For the cube the solution produces different stable solutions depending on the initial conditions. A fluid-particle path is shown to be a good method of illustrating the three-dimensional motion. In the cube and in each cubical cell in the channel this streak-line was found to consist of a pair of double helices.