We have studied single-mode mean field equations associated with variable viscosity convection for both steady state and time-dependent situations. Steady state mean field solutions can be obtained by treating the governing mean field equations as two coupled fourth-order differential systems, which are solved as a series of coupled two-point boundary value problems with underrelaxation. Steady solutions can be achieved in which there can exist viscosity contrasts exceeding 108 and interior Rayleigh numbers of 0(109). For comparison with two-dimensional solutions we have employed results derived from two finite element methods. Both temperature-dependent and temperature- and pressure-dependent viscosity with the surface viscosity fixed have been studied. Our results show that the differences between the two approaches grow with both viscosity contrast and convective vigor. The steady state Nusselt numbers and interior temperatures are greater in the case of mean field solutions. The power law index β governing the relationship between the interior Rayleigh number and the Nusselt number is larger for the mean field. Comparison of time-dependent solutions shows that one can monitor rather faithfully the evolution of th averaged interior temperature and surface heat flow over a long time scale with the mean field method. Initially, thermal instabilities originating from the mean field boundary layers are found to be correlated with two-dimensional boundary layer instabilities and are much more violent in character and cause large oscillations of the surface heat flow. The time scales associated with the secular variations of the Urey number, representing the ratio of the heat production and the surface heat loss, agree well between the two approaches. These results suggest that the mean field equations may serve as an efficient vehicle for studying planetary convection with complicated physics.