As a consequence of the Taylor–Proudman balance, a balance between the pressure, Coriolis and buoyancy forces in the radial and latitudinal momentum equations (that is expected to be amply satisfied in the lower solar convection zone), the superadiabatic gradient is determined by the rotation law and by an unspecified function of r, say, S′Ω(r), where r is the radial coordinate. If the rotation law and S′Ω(r) are known, then the solution of the energy equation, performed in this paper in the framework of the MLΩ formalism, leads to a knowledge of the Reynolds stresses, convective fluxes, and meridional motions. The MLΩ-formalism is an extension of the mixing length theory to rotating convection zones, and the calculations also involve the azimuthal momentum equation, from which an expression for the meridional motions in terms of the Reynolds stresses can be derived. The meridional motions are expanded as U r(r,θ)=P 2(cosθ)ψ2(r)/r 2ρ+P 4(cosθ)ψ4(r)/r 2 ρ+..., and a corresponding equation for U θ(r,θ). Here θ is the polar angle, ρ is the density, and P 2(cosθ), P 4(cosθ) are Legendre polynomials. A good approximation to the meridional motion is obtained by setting ψ4(r)=−Hψ2(r) with H≈−1.6, a constant. The value of ψ2(r) is negative, i.e., the P 2 flow rises at the equator and sinks at the poles. For the value of H obtained in the numerical calculations, the meridional motions have a narrow countercell at the poles, and the convective flux has a relative maximum at the poles, a minimum at mid latitudes and a larger maximum at the equator. Both results are in agreement with the observations.