We can calculate the accelaration of a convective fluid element \[\ddot{r} = - g - \frac{1}{\rho} \frac{\partial {\mathrm{P}}}{\partial r},\] which expanded at first order around the equilibrium state leads to an estimate of the total net force (buoyant minus gravitation) \(f=-g \Delta \rho\). Along the distance \(\Delta r\) a convective fluid element experience the force \(f(\Delta r) = - g \Delta \rho (\Delta r)\). Here the variation of \(g\) with height is considered negligible. The work done per unit volume over the distance \(\Delta r\) is therefore \[W (\Delta r) = \int_{0}^{\Delta r} f(x) {{\mathrm d}}x = - g \int_{0}^{\Delta r} \Delta \rho (x) {{\mathrm d}}x = - \frac{1}{2} g \Delta \rho (\Delta r) \Delta r.\] An average of \(W(\Delta r)\) over all possible values of \(\Delta r\) and a choice of a numerical factor (1/4) in exchanging \(\Delta r\) with \(\lambda\) leads to the following expression (\(W\) depends quadratically on \(\Delta r\)): \[\overline{W} (\lambda) = \frac{1}{4} W(\lambda) = - \frac{1}{8} g \Delta \rho (\lambda) \lambda.\] Under the assumption of no dissipative forces, no transfer of kinetic energy to surrounding matter and no heat loss, then one would conclude \(\bar{W}(\lambda) = \bar{1/2\rho v^2}\). However in the MLT one assumes that only half of this work is transformed into kinetic energy. This leads to \[\overline{\frac{1}{2}\rho v^2} \simeq \frac{1}{2}\rho \bar{v}^2 = \frac{1}{2}\overline{W} (\lambda) = - \frac{1}{16} g \Delta \rho (\lambda) \lambda,\] so that the velocity of convective elements at the level where the mixing length is \(\lambda\) can be expressed as \[\label{eq:velocity2} \bar{v}^2 = - \frac{1}{8} g (\Delta \rho (\lambda)/\rho) \lambda.\] Relating \(\Delta \rho\) and \(\Delta T\) requires the equation of state for the gas \(\rho = \rho (\mu, T, {\mathrm{P}})\), which in differential form can be written as \[{{\mathrm d}}\rho = \bigg(\frac{\partial \rho}{\partial \mu} \bigg)_{{\mathrm{P}},T} {{\mathrm d}}\mu + \bigg(\frac{\partial \rho}{\partial T} \bigg)_{{\mathrm{P}},\mu} {{\mathrm d}}T + \bigg(\frac{\partial \rho}{\partial {\mathrm{P}}} \bigg)_{\mu,T} {{\mathrm d}}{\mathrm{P}},\] , \[{{\mathrm d}}\rho = \frac{\rho}{\mu}\bigg(\frac{\partial \ln \rho}{\partial \ln \mu} \bigg)_{{\mathrm{P}},T} {{\mathrm d}}\mu + \frac{\rho}{T} \bigg(\frac{\partial \ln \rho}{\partial \ln T} \bigg)_{{\mathrm{P}},\mu} {{\mathrm d}}T + \frac{\rho}{{\mathrm{P}}} \bigg(\frac{\partial \ln \rho}{\partial \ln {\mathrm{P}}} \bigg)_{\mu,T} {{\mathrm d}}{\mathrm{P}}\] and finally \[\frac{{{\mathrm d}}\rho}{\rho} = \alpha \frac{{{\mathrm d}}{\mathrm{P}}}{{\mathrm{P}}} - \delta \frac{{{\mathrm d}}T}{T} + \phi \frac{{{\mathrm d}}\mu}{\mu}\] Where \( \phi \equiv \big(\frac{\partial \ln \rho}{\partial \ln \mu} \big)_{{\mathrm{P}},T}, \; \delta \equiv - \big(\frac{\partial \ln \rho}{\partial \ln T} \big)_{{\mathrm{P}},\mu}, \; \alpha \equiv \big(\frac{\partial \ln \rho}{\partial \ln {\mathrm{P}}} \big)_{\mu,T}\) and for an ideal gas \(\alpha = \delta = \phi = 1\), so that using pressure equilibrium (\({{\mathrm d}}{\mathrm{P}}= 0\)) we obtain \(\Delta \ln \rho = \Delta \ln T\). In the case of a mixture of perfect gas and radiation we can write \[\label{eq:deltarho} \Delta \ln \rho = -Q \Delta \ln T\] where \[Q = \frac{4-3 \beta}{\beta} - \bigg(\frac{\partial \ln \mu}{\partial \ln T} \bigg)_{{\mathrm{P}}},\] with \(\beta = {\mathrm{P}}_{Gas}/{\mathrm{P}}\). Plugging Eq. \ref{eq:deltarho} in \ref{eq:velocity2} one obtains \[\bar{v}^2 = \frac{1}{8} g Q (\Delta T (\lambda) / T) \lambda .\] and substituting the temperature excess as calculated in \ref{eq:deltat} \[\bar{v}^2 = \frac{1}{8} g Q \frac{\lambda^2}{{H_{\mathrm{P}}}}\, (\nabla - \nabla'),\] so that the average convective velocity can be calculated \[\label{eq:velocity} \bar{v} = \frac{1}{2 \sqrt{2}} \lambda \, \bigg(\frac{g\, Q}{{H_{\mathrm{P}}}}\bigg)^{1/2} \, (\nabla - \nabla')^{1/2}\\ = \frac{1}{2 \sqrt{2}} g \, \lambda \bigg(\frac{\rho\, Q }{{\mathrm{P}}}\bigg)^{1/2} \, (\nabla - \nabla')^{1/2}.\] Mixing length convective velocities depend linearly on the choice of the mixing length parameter \(\lambda\). Also they depend on \(Q^{1/2}\), which diverges for \(\beta \rightarrow 0\), that is in radiation dominated regions. In these regions MLT often results in supersonic convective velocities. Note that the above calculations hold only in the case of pressure equilibrium, which is only justified if \(\bar{v} < c_s\), i.e. for subsonic convection.
We can now rewrite the convective flux using the derived form of the convective velocities and Eq. \ref{eq:convectiveflux} \[\label{eq:convectiveflux2} F_{c} = \frac{1}{2} \rho \bar{v}\, {c_{\mathrm{P}}}\, T \frac{\lambda}{{H_{\mathrm{P}}}} ( \nabla - \nabla') \\ = \frac{1}{4 \sqrt{2}} g^2 \rho^{5/2} {c_{\mathrm{P}}}Q^{1/2} T {\mathrm{P}}^{-3/2} \lambda^2 \, (\nabla - \nabla')^{3/2}.\]