Smagorinsky model¶
(Loosly following Turbulent Flows, by S.B. Pope).
The following residual stresses require closure $$\tau_{ij}^R = \overline{U_iU_j} - \overline{U}_i\overline{U}_j.$$ Closure is required because the filtered product $\overline{U_iU_j}$ is unknown.
The deviatoric part of $\tau_{i,j}^R$ is given by $$\tau_{ij}^r \equiv \tau_{ij}^R-\frac{1}{3}\tau_{kk}^R\delta_{ij},$$ that is, $$\tau_{ij}^r= \overline{U_iU_j} - \overline{U}_i\overline{U}_j - \frac{1}{3}\left(\overline{U_kU_k} - \overline{U}_k\overline{U}_k\right)\delta_{ij}.$$
This is modeled as $$\tau_{ij}^{r,m} = -2\nu_r\overline{S}_{ij},$$ where $\overline{S}_{ij}$ is the filtered deviatoric rate of strain tensor given by $$\overline{S}_{ij} = \frac{1}{2}\left(\frac{\partial \overline{U}_i}{\partial x_j} + \frac{\partial \overline{U}_j}{\partial x_i}\right)-\frac{1}{3}\frac{\partial \overline{U}_k}{\partial x_k}\delta_{ij}.$$ The residual turbulent viscosity $\nu_r$ is given by $$\nu_r = l_S^2\overline{S},$$ where $$\overline{S} = \sqrt{2\overline{S}_{ij}\overline{S}_{ij}}$$ is the mean strain rate with units of inverse time, and $l_S^2$ is the square of the Smagorinsky lengthscale, taken as $$l_S^2 = c_S\bar{\Delta}^2.$$ These give $$\tau_{ij}^{r,m} = -2c_S\bar{\Delta}^2\overline{S}\,\overline{S}_{ij}.$$
Dynamic model¶
The dynamic model computes $c_s$ dynamically and locally, rather than relying on empiricism for its evaluation. This is done by filtering the resolved scale at a larger scale, called the test filter scale $\tilde{\bar{\Delta}}$, which is nominally twice the size of $\bar{\Delta}$. The test filter is applied to the residual stress $\tau^r$. By analogy, a residual stress $T^r$ is defined on the test filter relative to the base filter. The test filter filters the resolved scale, so test filtered products of based filtered quantities, like $\widetilde{\overline{U}_{ij}\,\overline{U}_{ij}}$, can be computed using the resolved scale quantities $\overline{U}_{ij}$. $T^r$ and $\tau^r$ are modeled similarly, both using the same $c_s$, which is assumed to apply to both the base and test filter scales. Importantly, both $\widetilde{\tau^r}$ and $T^r$ contain $\widetilde{\overline{U}_{ij}\,\overline{U}_{ij}}$, so subtracting $\widetilde{\tau^r}$ from $T^r$ results in cancellation of the $\widetilde{\overline{U}_{ij}\,\overline{U}_{ij}}$ terms. The result, called $\mathcal{L}$ can be computed and equated to the difference between the model for $T^r$ and the test-filtered model for $\tau^r$ to solve for $c_s$. The details follow.
Define the double filtered residual stress to be $$T_{ij}^R = \widetilde{\overline{U_iU_j}} - \widetilde{\overline{U}}_i\widetilde{\overline{U}}_j.$$ This is the residual stress of the resolved scale at the test filter scale. The deviatoric part is $$T_{ij}^r = \widetilde{\overline{U_iU_j}} - \widetilde{\overline{U}}_i\widetilde{\overline{U}}_j - \frac{1}{3}\left(\widetilde{\overline{U_kU_k}}-\widetilde{\overline{U}}_k\widetilde{\overline{U}}_k\right)\delta_{ij}.$$ This can be modeled analogously to $\tau^{r,m}_{ij}$ as $$T_{ij}^{r,m} = -2c_S\tilde{\bar{\Delta}}^2\widetilde{\overline{S}}\,\widetilde{\overline{S}}_{ij}$$
Now, apply the test filter to $\tau_{ij}^r$ and subtract it from $T_{ij}^r$ and call the result $\mathcal{L}_{ij}^d$:
$$\mathcal{L}_{ij}^d = T_{ij}^r-\widetilde{\tau}_{ij}^r = \begin{align} &\left({\widetilde{\overline{U_iU_j}}} - \widetilde{\overline{U}}_i\widetilde{\overline{U}}_j - \frac{1}{3}\left({\widetilde{\overline{U_kU_k}}}-\widetilde{\overline{U}}_k\widetilde{\overline{U}}_k\right)\delta_{ij}\right) \\ - &\left({\widetilde{\overline{U_iU_j}}} - \widetilde{\overline{U}_i\,\overline{U}_j} - \frac{1}{3}\left({\widetilde{\overline{U_kU_k}}} - \widetilde{\overline{U}_k\overline{U}_k}\right)\delta_{ij}\right) \end{align}$$
This filtering and subtraction results in cancellation of the unknown filtered products $\widetilde{\overline{U_iU_j}}$ and $\widetilde{\overline{U_kU_k}}$, giving $$\mathcal{L}_{ij}^d = \underbrace{\widetilde{\overline{U}_i\overline{U}_j} -\widetilde{\overline{U}}_i\widetilde{\overline{U}}_j}_{\mathcal{L}_{ij}} - \frac{1}{3}\left(\widetilde{\overline{U}_k\overline{U}_k} -\widetilde{\overline{U}}_k\widetilde{\overline{U}}_k\right)\delta_{ik}.$$ Similarly, for the modeled terms, let $$\mathcal{L}_{ij}^{d,m} = T_{ij}^{r,m}-\widetilde{\tau}_{ij}^{r,m},$$ which gives $$\mathcal{L}_{ij}^{d,m} = c_S\left(-2\tilde{\bar{\Delta}}^2\widetilde{\overline{S}}\,\widetilde{\overline{S}}_{ij} + 2\bar{\Delta}^2\widetilde{\overline{S}\,\overline{S}}_{ij}\right).$$
Finally, we equate $\mathcal{L}_{ij}^{d,m}$ and $\mathcal{L}_{ij}^d$, $$\mathcal{L}_{ij}^{d} = c_S\underbrace{\left(-2\tilde{\bar{\Delta}}^2\widetilde{\overline{S}}\,\widetilde{\overline{S}}_{ij} + 2\bar{\Delta}^2\widetilde{\overline{S}\,\overline{S}}_{ij}\right)}_{M_{ij}}.$$ or $$\mathcal{L}_{ij}^{d} = c_SM_{ij}$$ This gives an equation for $c_S$, but $\mathcal{L}_{ij}^d$ and $M_{ij}$ have five independent components since the tensors are symmetric (giving six components), and deviatoric, which constrains the trace to equal zero, reducing the independent components from six to five. $c_S$ is chosen to minimize the square error.
The sum square error between $\mathcal{L}_{ij}^{d,m}$ and $\mathcal{L}_{ij}^{d}$ is \begin{align} \epsilon_{err} &= \left(\mathcal{L}_{ij}^{d,m} - \mathcal{L}_{ij}^{d}\right)\left(\mathcal{L}_{ij}^{d,m} - \mathcal{L}_{ij}^{d}\right), \\ &=\left(c_SM_{ij}-\mathcal{L}_{ij}^d\right)\left(c_SM_{ij}-\mathcal{L}_{ij}^d\right), \\ &=c_S^2M_{ij}M_{ij}-2c_SM_{ij}\mathcal{L}_{ij}^d + \mathcal{L}_{ij}^d\mathcal{L}_{ij}^d. \end{align} Then $$\frac{d\epsilon_{rr}}{dc_S} = 2c_SM_{ij}M_{ij}-2M_{ij}\mathcal{L}_{ij}^d.$$ Setting this to zero and solving for $c_S$ gives $$c_S = \frac{M_{ij}\mathcal{L}_{ij}^d}{M_{kl}M_{kl}}.$$ This minimizes the error since $$\frac{d^2\epsilon_{rr}}{dc_S^2}=2M_{ij}M_{ij} > 0.$$
Note, $\mathcal{L}_{ij}^d$ can be replaced with $\mathcal{L}_{ij}$ in the above expression for $c_S$. This follows from $$M_{ij}\mathcal{L}_{ij} = M_{ij}\left(\mathcal{L}_{ij}^d+\frac{1}{3}\mathcal{L}_{kk}\delta_{ij}\right) = M_{ii}\mathcal{L}_{ij}^d + \frac{1}{3}\mathcal{L}_{kk}M_{ij}\delta_{ij} = M_{ii}\mathcal{L}_{ij}^d,$$ where $M_{ij}\delta_{ij} = \mathrm{trace}(M_{ij})=0$ since $M_{ij}$ is deviatoric, and $M_{ij}$ is deviatoric since $S_{ij}$ is deviatoric.
Summary¶
$$\tau_{ij}^{r,m} = -2c_S\bar{\Delta}^2\overline{S}\,\overline{S}_{ij}.$$
$$\overline{S}_{ij} = \frac{1}{2}\left(\frac{\partial \overline{U}_i}{\partial x_j} + \frac{\partial \overline{U}_j}{\partial x_i}\right)-\frac{1}{3}\frac{\partial \overline{U}_k}{\partial x_k}\delta_{ij}.$$
$$\overline{S} = \sqrt{2\overline{S}_{ij}\overline{S}_{ij}}$$
$$c_S = \frac{M_{ij}\mathcal{L}_{ij}}{M_{kl}M_{kl}}.$$
$$M_{ij} = -2\tilde{\bar{\Delta}}^2\widetilde{\overline{S}}\,\widetilde{\overline{S}}_{ij} + 2\bar{\Delta}^2\widetilde{\overline{S}\,\overline{S}}_{ij}$$
$$\mathcal{L}_{ij}^d = \widetilde{\overline{U}_i\overline{U}_j} -\widetilde{\overline{U}}_i\widetilde{\overline{U}}_j - \frac{1}{3}\left(\widetilde{\overline{U}_k\overline{U}_k} -\widetilde{\overline{U}}_k\widetilde{\overline{U}}_k\right)\delta_{ik}.$$
Spatial variations¶
The above formulation assumes that $c_S$ is spatially uniform. There are also issues if $M_{kl}M_{kl}$ is zero (the denominator of the $c_S$ equation. Pope gives a discussion and provides a few references. In particular, see