Explore different SGS models

[ashwinvis]

• Dynamic Smagorinsky Model: present in turbChannel example in NekExamples. See also: Turb Channel must go Nek5000/NekExamples#105. More computationally expensive, may not be needed. Wall damping function is not hardcoded, and but automatic / dynamic

• Algebraic model for computing unresolved K.E + Deodroff.

• Vreman model?: no wall damping required.

• Piomelli Guertz: See https://doi.org/10.1017/jfm.2015.29 . Grid scale independent SGS model, perhaps useful for spectral element methods, because non-uniform grid.

• Shear improved Smagorinsky model: See https://doi.org/10.1017/S0022112006003429. Simple to implement, take away the mean shear component and only fluctuation shear is retained.

• Algebraic Minimum Dissipation (AMD) model: Rozema et al., “Minimum-Dissipation Models for Large-Eddy Simulation.”

[timofeymukha]

@ashwinvis I'm wondering if the way the squard length-scale is computed is slightly incorrect.
Looking at this code

         di = (xm1(ip,j,k,e)-xm1(im,j,k,e))**2
 $            + (ym1(ip,j,k,e)-ym1(im,j,k,e))**2
 $            + (zm1(ip,j,k,e)-zm1(im,j,k,e))**2

           dj = (xm1(i,jp,k,e)-xm1(i,jm,k,e))**2
 $            + (ym1(i,jp,k,e)-ym1(i,jm,k,e))**2
 $            + (zm1(i,jp,k,e)-zm1(i,jm,k,e))**2

           dk = (xm1(i,j,kp,e)-xm1(i,j,km,e))**2
 $            + (ym1(i,j,kp,e)-ym1(i,j,km,e))**2
 $            + (zm1(i,j,kp,e)-zm1(i,j,km,e))**2

           di = di/(ip-im)
           dj = dj/(jp-jm)
           dk = dk/(kp-km)
           dg2(i,j,k,e) = (di*dj*dk)**ndim_inv

Consider a grid aligned with x-y-z. Then di is Δx^2/2, dj is Δy^2/2, dk is Δz^2/2. I would say that there is an additional 1/2 factor missing here, since what we want e.g. Δx/2 to represent the scale along x, and afterwards we square it. The resulting difference is a factor of 2, your code giving x2 of what I think would be correct.

[timofeymukha]

@ashwinvis Pretty much convinced regarding the above now. I measured some distances in a particular channel flow case in Paraview, and also cross-checked with the length scale computation in the code for the minimal dissipation model I got from J. Bae. If we change the code to

di = sqrt(di)/(ip-im)
dj = sqrt(dj)/(jp-jm)
dk = sqrt(dk)/(kp-km)
dg2(i,j,k,e) = (di*dj*dk)**(2*ndim_inv)

results are consistent.

[ashwinvis]

Now your argument is much more clear to me and it makes sense. Thanks!