To take into account the effect of turbulent fluctuations on a flow field, the Navier-Stokes equations are modified to include such effects. The obtained equations are called Reynolds-averaged Navier-Stokes (RANS) equations.
As an example, the steady incompressible momentum equation can be written as follows (in Einstein tensor notation):
$$ \rho \frac{\partial u_i}{\partial t} + \rho u_i \frac{\partial u_i}{\partial x_j} = -\frac{\partial p}{\partial x_i} + \frac{\partial t_{ij}}{\partial xj}$$
Where $t_{ij}$ is the viscous stress tensor defined as:
$$ t_{ij} = 2\mu s_{ij}$$
where $\mu$ is viscosity and $s_{ij}$ is the strain-rate tensor:
$$ s_{ij} = \frac{1}{2}(\frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i})$$
After time-averaging, one can obtain the following equation:
$$ \rho \frac{\partial U_i}{\partial t} + \rho U_i \frac{\partial U_i}{\partial x_j} = -\frac{\partial P}{\partial x_i} + \frac{\partial }{\partial xj}(2\mu s_{ij} - \overline{\rho {u_j}'{u_i}'})$$
The resulted term $- \overline{\rho {u_j}`{u_i}`}$ is called the Reynolds-stress tensor (or turbulent shear stress), which is the basis of all turbulence models.
The Boussinesq approximation is merely assuming that turbulent shear stress is analogous to viscous shear stress (by introducing a new term called eddy viscosity ${\mu}_t$). Thus, one can write:
$$ - \overline{\rho {u_j}'{u_i}'} = {\mu}_t \frac{\partial U_i}{\partial x_j}$$
The approximation didn't solve the closure problem of turbulence, but it was used to come up with turbulent models that can model that eddy viscosity. By definition it's not a complete model, so there is no direct option for this in Fluent.
