Unclear how heat interacts with Navier Stokes I am playing around with an Navier stokes solver and I'm having trouble introducing heat.
Am I right in thinking this would be introduced in the ${\bf f}$ term of ${\partial{\bf u}\over\partial t} = -({\bf u}\cdot\nabla){\bf u}+v\nabla^2{\bf u}+{\bf f}$?
I find the lack of mass term disturbing. Please reassure me!
Also, how do I calculate the force vector, given that I have computed the scalar $dQ\over dt$ from Newton's law of cooling? Or is that the wrong approach?
 A: The Navier-Stokes equation to which you refer is more generally the first moment of velocity of the Boltzmann equation.  In order to get a proper connection to heating, you need a second-velocity-moment Navier-Stokes equation.  The Boltzmann equation keeps track of distributions of particles.  This changes the question from "What is the density and flow of a fluid at a point $\mathbf{x}$ at a time $t$?" to "What is the probability of finding a particle between $\mathbf{x}$ and $\mathbf{x}+\mathrm{d}\mathbf{x}$, with a velocity between $\mathbf{v}$ and $\mathbf{v}+\mathrm{d}\mathbf{v}$, at time $t$?"  A nice transition between the two formalisms is discussed in these notes (although they include gravity as an external force the way plasma physicists include the Lorentz force... just imagine the equations without those terms for a plain fluid, keeping only the collisonal term).
It's worth noting that each moment depends on a term from the next higher moment ($d\rho/d t$ depends on $\mathbf{u}$, $d\mathbf{u}/d t$ depends on $\overleftrightarrow{P}$,  $d E/d t$ (which is the same moment as $\overleftrightarrow{P}$) depends on the heat conduction, which is a 3rd order moment... Any equation that cuts off has to assume some kind of closure method.  For example, to close at first order, you might assume that the pressure is isotropic.  Or to close at second order, you might assume that the conduction is infinite (compared to the timescales of interest).
To answer your specific question, volumetric heating can result in a change in pressure, but you need an equation of state linking pressure, temperature and density.  (Heating steam will have a very different response from heating water.)  The modified pressure term can in turn couple to the first velocity moment (the Navier-Stokes equation you have written).
A: To allow for heat effects in a fluid, you need to couple the Navier-Stokes equations, (momentum conservation) which BTW contain the continuity equation for mass conservation too, to the energy (or temperature) equation (energy conservation).
Momentum dissipation in the momentum equation
$$
\frac{\partial v}{\partial t} + (\vec{v}\cdot\nabla)\vec{v} 
  = -\frac{\nabla p}{\rho} + \frac{1}{\rho}\nabla S
$$
is more correctly described by the divergence of the symmetric stress tensor $S$
$$
S = \rho \,\nu (\nabla \circ\vec{v} + (\nabla \circ\vec{v})^{T}) + \rho \,\eta\, I (\nabla \cdot \vec{v})
$$
The coefficients $\mu$ and $\nu$ denote the dynamic and kinematic viscosity respectively, $\circ$ is the tensor (outer) product, and $I$ is the unity tensor.
In the energy equation, the momentum dissipation leads to a corresponding positive definite dissipation (frictional heating) $\epsilon$
$$
\epsilon = \frac{1}{\rho}(S \, \nabla) \cdot \vec{v}
$$
The energy equation can, dependent on the system considered, have different sources of diabatic heating apart from the radiative heating (of which Newtonian cooling es a spacial case), such as latent heating due to phase transitions for example.
A: In fluid dynamics, especially in modeling, there are different flavors for including heat.
First of all, heat can be a passive tracer which does not influence the flow, this basically means, that you solve the scalar transport-diffusion equation for the temperature (or energy if you which), which reads, in your notation
$${\partial{T}\over\partial t} +({\bf u}\cdot\nabla){T}=\alpha\nabla^2{T}+Q$$
where $\alpha$ is thermal diffusion coefficient and $Q$ are heat sources (which may be dependent on $T$ for radiation).
The density of fluids is temperature-dependent. You should incorporate this in the Navier-Stokes equations. If the density differences are only relevant for the buoyant terms (e.g. gravity), it is enough to include the effect in $\bf f$, by the Boussinesq approximation. This basically means that you linearize the $\rho(T)$ curve and include the force term as $\beta g (T-T_{ref})$ with $\beta$ a thermal expansion coefficient.
The next step in complexity, arises when $\rho$ changes such, that it takes over control for more terms. The Navier-Stokes equation you gave, assumes constant density, which therefor drops out of the equations. This is generally not true, and you have to take into account the variable density in all terms.
