Stratified plane Couette flow test

[ali-ramadhan]

We agreed that stratified plane Couette flow could be a test case against which we can verify the implementation of AMD in PR #309, although the published results [Vreugdenhil & Taylor (2018)] use a slightly modified version of AMD.

They report a number of LES runs, and since we don't have a vertically stretched grid, it might be easiest to try and reproduce the Pr = 0.7 case which they do with (Nx, Ny, Nz) = (64, 49, 64) grid points [for Ri=0 or 0.01, for Ri = 0.04 it's (64, 65, 64)]. I might suggest focusing on the Ri = 0 case as they also run a resolved DNS with (Nx, Ny, Nz) = (256, 129, 256) against which we might be able to compare (might not be possible, depending on vertical grid stretching).

Unfortunately, I don't think we can reproduce their results without a vertically stretched grid... If I understand the paper correctly, the grid is stretched according to

y_j = h*tanh(Sf*(2*(j-1)/(Ny-1))) / tanh(Sf)

so for h=1 you get grid spacings of Δy ~ 0.125 away from the wall and Δy ~ 0.0000035 adjacent to the wall. So a faithful reproduction would need 285,000+ vertical levels lol.

I don't fully understand how they got their values for the vertical grid cell size adjacent to the wall ∆y_w^+, but it seems like if h = 100,000 then our values agree (they get a spacing of Δy ~ 0.35 adjacent to the wall). Either way, the ratio between the thickest and thinnest spacings is ~35,000.

Not sure if there's still a way we could compare results in this case without a vertically stretched grid...

Absolute worst case scenario, I've been working on an implementation of a vertically stretched grid (see PRs #283 and #306) but it's a work in progress and might not be the best use of our efforts right now.

Reference: Catherine A. Vreugdenhil and John R. Taylor, Large-eddy simulations of stratified plane Couette flow using the anisotropic minimum-dissipation model, Physics of Fluids 30, 085104 (2018).

cc @glwagner @rafferrari

[rafferrari]

Hi all,

Greg, do you know if Cat Vreugdenhi is visiting GFD this summer? It may be easy to ask her to run an example without stretched grid. That said, do we have a good sense of how important a stretched grid is for the accuracy of our LES simulations?

Raffaele

On Jul 10, 2019, at 9:32 AM, Ali Ramadhan notifications@github.com wrote:

We agreed this could be a test case against which we can verify the implementation of AMD in PR #309, although the published results use a slightly modified version of AMD.

They report a number of LES runs, and since we don't have a vertically stretched grid, it might be easiest to try and reproduce the Pr = 0.7 case which they do with (Nx, Ny, Nz) = (64, 49, 64) grid points [for Ri=0 or 0.01, for Ri = 0.04 it's (64, 65, 64)]. I might suggest focusing on the Ri = 0 case as they also run a resolved DNS with (Nx, Ny, Nz) = (256, 129, 256) against which we might be able to compare (might not be possible, depending on vertical grid stretching).

Unfortunately, I don't think we can reproduce their results without a vertically stretched grid... If I understand the paper correctly, the grid is stretched according to

y_j = htanh(Sf(2*(j-1)/(Ny-1))) / tanh(Sf) so for h=1 you get grid spacings of Δy ~ 0.125 away from the wall and Δy ~ 0.0000035 adjacent to the wall. So a faithful reproduction would need 35000+ vertical levels lol.

I don't fully understand how they got their values for the vertical grid cell size adjacent to the wall ∆y_w^+, but it seems like if h = 100,000 then our values agree (they get a spacing of Δy ~ 0.35 adjacent to the wall). Either way, the ratio between the thickest and thinnest spacings is ~35,000.

Not sure if there's still a way we could compare results in this case without a vertically stretched grid...

Absolute worst case scenario, I've been working on an implementation of a vertically stretched grid (see PRs #283 and #306) but it's a work in progress and might not be the best use of our efforts right now.

Reference: Catherine A. Vreugdenhil and John R. Taylor, Large-eddy simulations of stratified plane Couette flow using the anisotropic minimum-dissipation model, Physics of Fluids 30, 085104 (2018).

cc @glwagner @rafferrari

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[ali-ramadhan]

Alternatively, we can compare with the original databases of stratified plane Couette flow DNS runs by Deusebio et al. (2015) who varied Re and Ri, and Zhou et al. (2017) who varied Pr and Ri.

That way we can design our own experiment, but maybe we still can't do anything rigorous without high-resolution near the walls.

Maybe we could keep increasing the vertical resolution and see if our solutions (the statistics) converge towards the resolved DNS solution? We could go from Nz = 50 to Nz = 2000 levels or something.

References:

• Enrico Deusebio, C. P. Caulfield, and J. R. Taylor, The intermittency boundary in stratified plane Couette flow, Journal of Fluid Mechanics 781, pp. 298-329 (2015).
• Qi Zhou, John R. Taylor, and C. P. Caulfield, Self-similar mixing in stratified plane Couette flow for varying Prandtl number, Journal of Fluid Mechanics 820, pp. 86-120 (2017).

[johncmarshall54]

Ali, i thought the stretching was sorted now, after the hackathon. Without stretching the code will have limited use. John

On Wed, Jul 10, 2019, 3:32 PM Ali Ramadhan notifications@github.com wrote:

We agreed this could be a test case against which we can verify the implementation of AMD in PR #309 https://github.com/climate-machine/Oceananigans.jl/pull/309, although the published results use a slightly modified version of AMD.

They report a number of LES runs, and since we don't have a vertically stretched grid, it might be easiest to try and reproduce the Pr = 0.7 case which they do with (Nx, Ny, Nz) = (64, 49, 64) grid points [for Ri=0 or 0.01, for Ri = 0.04 it's (64, 65, 64)]. I might suggest focusing on the Ri = 0 case as they also run a resolved DNS with (Nx, Ny, Nz) = (256, 129, 256) against which we might be able to compare (might not be possible, depending on vertical grid stretching).

Unfortunately, I don't think we can reproduce their results without a vertically stretched grid... If I understand the paper correctly, the grid is stretched according to

y_j = htanh(Sf(2*(j-1)/(Ny-1))) / tanh(Sf)

so for h=1 you get grid spacings of Δy ~ 0.125 away from the wall and Δy ~ 0.0000035 adjacent to the wall. So a faithful reproduction would need 35000+ vertical levels lol.

I don't fully understand how they got their values for the vertical grid cell size adjacent to the wall ∆y_w^+, but it seems like if h = 100,000 then our values agree (they get a spacing of Δy ~ 0.35 adjacent to the wall). Either way, the ratio between the thickest and thinnest spacings is ~35,000.

Not sure if there's still a way we could compare results in this case without a vertically stretched grid...

Absolute worst case scenario, I've been working on an implementation of a vertically stretched grid (see PRs #283 https://github.com/climate-machine/Oceananigans.jl/pull/283 and #306 https://github.com/climate-machine/Oceananigans.jl/pull/306) but it's a work in progress and might not be the best use of our efforts right now.

Reference: Catherine A. Vreugdenhil and John R. Taylor, Large-eddy simulations of stratified plane Couette flow using the anisotropic minimum-dissipation model https://aip.scitation.org/doi/abs/10.1063/1.5037039, Physics of Fluids 30, 085104 (2018).

cc @glwagner https://github.com/glwagner @rafferrari https://github.com/rafferrari

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/climate-machine/Oceananigans.jl/issues/310?email_source=notifications&email_token=AKXUEQWVAM5WA445JVSWDZ3P6XQFBA5CNFSM4H7PMRLKYY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4G6LLIOA, or mute the thread https://github.com/notifications/unsubscribe-auth/AKXUEQSSGRFFX26YFGPJV43P6XQFBANCNFSM4H7PMRLA .

[glwagner]

@ali-ramadhan thanks for looking into that!

I think a convergence test with coarse horizontal resolution and varying vertical resolution is probably a good idea, perhaps producing a plot like Fig 11 in Vreugdenhil and Taylor.

In my reading, ∆y_w^+ is the vertical spacing in 'wall units', which is the physical spacing divided by nu/u_\tau, where nu is viscosity and u_\tau is the friction velocity defined in eq 16.

In Table II, values for both ∆y_w^+ and ∆y_c^+ (grid spacing at wall and at domain center) are given. It looks like many of the runs may be feasible. Here's some code:

 y(j, Ny, Sf, h) = h / tanh(Sf) * tanh(Sf*(2*(j-1)/Ny - 1)) 
Δy(j, Ny, Sf, h) = y(j+1, Ny, Sf, h) - y(j, Ny, Sf, h)

N_unstretched(Ny_VT, Sf) = ceil(Int, 2 / Δy(2, Ny_VT, Sf, 1)) 

@show N_unstretched(49, 3)
@show N_unstretched(65, 2.5)
@show N_unstretched(97, 3)
@show N_unstretched(129, 2)

It looks like some of their runs are feasible, since

julia> N_unstretched(49, 3) # runs 2-4
1141

julia> N_unstretched(65, 2.5) # run 5
768

julia> N_unstretched(97, 3) # runs 7-11
2710

julia> N_unstretched(129, 2) # run 12
805

Runs B-G (which are included to demonstrate resolutions that are apparently inadequate) may also be feasible.

One concern about the suitability of stratified Couette flow is mentioned at the end of the paper...

Stratified plane Couette flow is a challenging test case for the LES model because it has a linearly stable laminar state which introduces requirement 1. The results of Abkar and Moin13 suggest that AMD LES performs even better in other stratified wall-bounded flows. Balancing all these concerns is key to using the AMD model in LES of wall-bounded stratified flow. Nevertheless the AMD model is able to capture turbulent intermittency and mean and turbulent flowproperties in stratified plane Couette flow.

@johncmarshall54 we are planning to implement stretched vertical grid functionality after the LES is implemented.

[ali-ramadhan]

Ali, i thought the stretching was sorted now, after the hackathon. Without stretching the code will have limited use. John

Yes well, we know how to do a very fast Poisson solve in a vertically stretched grid but still working on integrating it (PR #306). The fast solver (which relies on cyclic reduction) may need some massaging to work with Neumann boundary conditions for pressure... The plan is to merge the vertically stretched grid with the finite volume operators.

[ali-ramadhan]

@glwagner Never mind, I am an idiot. I divided by Ny-1 by mistake. The matching factor of 0.35 might have been pure coincidence.

We need closer to O(1000) vertical levels, which should be pretty doable.

[ali-ramadhan]

Temperature and velocity profiles look qualitatively correct which is promising but the sub-grid-scale viscosity and diffusivity don't match. At least the order of magnitude on each is roughly correct.

Most notably their version of AMD allows the SGS ν and κ to go below the background values, which seems essential for ν, although not sure what that physically means. Seems that the VerstappenAnisotropicMinimumDissipation closure in the new-closures branch allows for this so I'll try rerunning with it.

The SGS profiles seem kind of noisy but maybe I just need to bump up the resolution.

The SGS/turbulent Prandtl number is still too low but this should get fixed if the SGS ν and κ profiles are right.

The wall velocity Uw=1 so the x-axis is U/Uw.

The wall temperature is Tw = 0.01 so this plot should match the profile from Vreugdenhil & Taylor.

[glwagner]

I won’t get the chance to look closely at these profiles until later, but here are some thoughts in the intervening time:

The nonlinear, subgrid scale diffusivities model the effect of unresolved, subgrid scale turbulence on the resolved flow. This is their physical meaning. When the SGS diffusicities are 0, this physically means that there is no subgrid scale turbulence. This occurs when turbulence is completely resolved or if there is no turbulence, as in the case of laminar flow. In the stratified couette example, the near-wall motions are “resolved”, so that near the wall the SGS diffusivities should vanish. This is because the near wall motions are nearly laminar (very strongly sheared) with only small turbulence fluctuations. The interaction of the flow with the wall is not part of the AMD model, which means that this part of the flow must either be resolved or modeled with a “wall model”, which somehow accounts for the effect of wall effects on the nonlinear subgrid scale diffusivities.

When the nonlinear diffusivities are calculated in Oceananigans, the constant molecular contribution is added so that the resulting diffusivity field is the sum of the nonlinear and molecular components. Therefore to calculate the nonlinear part correctly from a saved diffusivity field, the molecular value must subtracted. This behavior is identical between RozemaAnisotropicMinimumDissipation and VerstappenAnisotropicMinimimDissipation.

Because of the fact that this experiment requires near wall motions to be resolved to be correct, the results should depend strongly on resolution.

On Sun, Aug 18, 2019 at 7:11 AM Ali Ramadhan notifications@github.com wrote:

Temperature and velocity profiles look qualitatively correct which is promising but the sub-grid-scale viscosity and diffusivity don't match. At least the order of magnitude on each is roughly correct.

Most notably their version of AMD allows the SGS ν and κ to go below the background values, which seems essential for ν, although not sure what that physically means. Seems that the VerstappenAnisotropicMinimumDissipation closure in the new-closures branch allows for this so I'll try rerunning with it.

The SGS profiles seem kind of noisy but maybe I just need to bump up the resolution.

The SGS/turbulent Prandtl number is still too low but this should get fixed if the SGS ν and κ profiles are right.

The wall velocity Uw=1 so the x-axis is U/Uw. [image: plots_stratified_couette_flow_u_profiles] https://user-images.githubusercontent.com/20099589/63219768-2343de80-c147-11e9-8591-7d52fe4a27c5.png

The wall temperature is Tw = 0.01 so this plot should match the profile from Vreugdenhil & Taylor. [image: plots_stratified_couette_flow_T_profiles] https://user-images.githubusercontent.com/20099589/63219770-26d76580-c147-11e9-83e1-f8a076dc7810.png

[image: correct_profiles] https://user-images.githubusercontent.com/20099589/63219833-4fac2a80-c148-11e9-8a0e-55d8b31b541a.png

[image: plots_stratified_couette_flow_nu_profiles] https://user-images.githubusercontent.com/20099589/63219836-55097500-c148-11e9-811c-e520e8dd1c4b.png

[image: plots_stratified_couette_flow_kappaT_profiles] https://user-images.githubusercontent.com/20099589/63219837-589cfc00-c148-11e9-8ff9-e9afc5f69472.png

[image: correct_SGS_profiles] https://user-images.githubusercontent.com/20099589/63219838-5b97ec80-c148-11e9-800e-7b3a9923511e.png

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/climate-machine/Oceananigans.jl/issues/310?email_source=notifications&email_token=ADUQQBTLE6XCGRLSC2QHD7LQFEU5DA5CNFSM4H7PMRLKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD4Q5U2Y#issuecomment-522312299, or mute the thread https://github.com/notifications/unsubscribe-auth/ADUQQBWKU7TT45JZYQO5CJTQFEU5DANCNFSM4H7PMRLA .

[ali-ramadhan]

Thanks for the explanation @glwagner.

I didn't realize that the background diffusivity had to be subtracted from the nonlinear diffusivity. In that case the ν profile isn't super far off although both nonlinear diffusivity profiles (especially the κ profile) suggest that perhaps the boundary layer is not being resolved.

I'm rerunning the simulation but with the vertical resolution tripled. Hopefully this will improve the comparison.

[glwagner]

@ali-ramadhan a few more comments:

• I don't understand why the temperature profiles you've plotted have spikes in them, especially since they are horizontal averages. What's going on? The sub grid scale diffusivities also seem to have far too much variance.

• Can you plot multiple equilibrated horizontally averaged solutions for differing model constant C? I think we may need to determine the correct value of this constant for our numerical method in order for our results to match. It looks like the SGS diffusivity values are too small, so this might help...

• We can't be sure what resolution an Oceananigans.jl model needs to be to match Diablo results, since Diablo uses a spectral discretization in the horizontal. That's something to keep in mind. Note also that V&T take the spectral discretization into account when defining their 'filter widths' in their AMD implementation.

• If our results matched Diablo without resolving the viscous boundary layer near the walls, I think this would not constitute a verification of Oceananigans (instead, it would indicate a problem with Oceananigans, or the presence of some kind of 'cancellation of errors'). We must resolve the viscous boundary layers near the walls in this experiment for the results to be valid.

For reference, the AMD viscosity is calculated here:

https://github.com/climate-machine/Oceananigans.jl/blob/352c018b6cb8801500bae61ac858ddc146bf3f48/src/turbulence_closures/anisotropic_minimum_dissipation.jl#L50

the diffusivities are similar.

[ali-ramadhan]

Thanks for the comments! I've rerun the simulation with Ri = 0.04 so that the viscous boundary layer is resolved. Horizontal resolution may still be an issue but I can roughly double it. We may also have to figure out which C to use, but might be more useful to do when we know we have a run that looks reasonable?

Just posting vertical profiles from that run. They look better but non-linear diffusivity profiles are still off. The shape of the profiles looks different from the figure I posted above. The magnitudes also seem very low (they're barely above the background values). Magnitude could be a plotting error (if it is I can't see it...).

Discrepency could be because I'm doing the Pr = 0.7 case while the profiles posted in V&T figure 7 are for the Pr = 7 case. I assume things will be different, but I don't have the intuition to tell how big of a difference to expect.

Note: nu and kappaT profiles are nu_SGS/nu and kappa_SGS/kappa.

[ali-ramadhan]

Just reporting some dimensionless numbers for that run (it's number 5 from table 1 of V&T).

The friction Reynolds number is a little different between the two walls but is generally between 230-250. This does not match run 5's value of 183 but seems to match the Ri=0 values. If we're not resolving something, maybe we're effectively running at a lower Ri number? Not sure if it should be fluctuating this much?

The friction Nusselt number is quite different between the top and bottom walls, but generally around 5 and 7. Bottom wall matches table 1's value of 7.1 pretty closely. Again, not sure if it should fluctuate this much, and if so how to time average properly.

---

[ali-ramadhan]

Just to document what's been done:

So I ran the Ri = 0.01, Pr = 0.7 case that Cat ran at three different horizontal resolutions and unfortunately the results are quite different from hers:

1. U and T don't fully go to zero at the top and bottom. This could be because we only have 64 vertical levels, not enough to resolve the boundary layers, but otherwise the profiles look qualitatively correct.
2. Our νSGS/ν and κSGS/κ are too low by a factor of ~2. They decrease as you increase the resolution which makes sense. κ profile shape looks similar to Cat's. ν profile shape should be flatter.
3. Re is roughly a factor of 2 too large. Nu is too large by ~60%. (edited) I might try changing the Poincare "constant" to see if it makes a difference here. Currently C = 1/12, so I'll try 1/24 and 1/6. Might also try 128 vertical levels, but that's no longer a fair comparison.

[ali-ramadhan]

To make sure it's not just low resolution being bad, I reran everything at 256x256x128 (to keep aspect ratio the same). SGS diffusivity profiles look okay at C = 1/6 (although shapes still don't fully match Cat's) so it's also possible that I'm calculating Re and Nu incorrectly...

[ali-ramadhan]

So going up to 256x256x128 didn't help. Re and Nu don't actually change much. Re and Nu are still too big by a factor of ~2...

[ali-ramadhan]

Apparently I'm just an idiot.

h = model.grid.Lz
Reτ = h * √uτ²⁺ / ν
Nu = (qw⁺ * h)/(κ * Θw)

The height of the domain is actually 2h based on the domain geometry described in Cat's Table 1. So I should have used h = model.grid.Lz / 2.

There's the factor of 2.

[glwagner]

🤷‍♂️🧠🏃‍♂️🙉🦉

[glwagner]

How does the flow field look?

[ali-ramadhan]

Rerunning Ri = 0.01 and Ri = 0.04 to make sure both come out matching Cat's simulations.

Will also plot some flow fields from the 256x256x128 runs to see what they look like.

[rafferrari]

I am not sure I understand the definitions of Re and Nu. Re should be a velocity scale times a length scale divided my nu. Nu should be the ratio of the adjective heat flux w theta over the diffusion heat flux kappa grad theta.

Raffaele

On Sep 6, 2019, at 7:31 AM, Ali Ramadhan notifications@github.com wrote:

Apparently I'm just an idiot.

h = model.grid.Lz Reτ = h √uτ²⁺ / ν Nu = (qw⁺ h)/(κ * Θw) The height of the domain is actually 2h based on the domain geometry described in Cat's Table 1. So I should have used h = model.grid.Lz / 2.

There's the factor of 2.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[glwagner]

@rafferrari

√uτ²⁺ is a velocity scale (it's written oddly I think...), h is a length scale, and ν is ν.

@ali-ramadhan should correct me, but it looks like κ * Θw / h is a diffusive heat flux scaling across the domain, whereas qw⁺ is the wall temperature flux (not sure how its calculated).

[glwagner]

The Nusselt number here is the ratio between the measured turbulent temperature flux and the temperature flux we expect in the laminar solution?

[rafferrari]

I thought one of the advantages of Julia is that one can use easy to read notation. The notation for Re and Nu is anything but intuitive. Should we discuss how to improve on that?

Raffaele

On Sep 6, 2019, at 8:08 AM, Gregory L. Wagner notifications@github.com wrote:

@rafferrari

√uτ²⁺ is a velocity scale (it's written oddly I think...), h is a length scale, and ν is ν.

@ali-ramadhan should correct me, but it looks like κ * Θw / h is a diffusive heat flux scaling across the domain, whereas qw⁺ is the wall temperature flux (not sure how its calculated).

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[ali-ramadhan]

@rafferrari we're using the frictional Re and Nu as defined by equation (20) of Vreugdenhil & Taylor (2018).

frictional Re = (friction velocity) * (domain height / 2) / ν

Re is still a velocity scale U times a length scale L divided by ν just with specific choices for U and L.

advective heat transfer = (wall heat flux qw) / (wall temperature Θw)
conductive heat transfer = κ / (domain height / 2)
frictional Nu = (advective heat transfer) / (conductive heat transfer)

So Nu is still the ratio of advective to conductive heat transfer, just with specific definitions again. Well, my "heat transfers" have units of m/s I guess so maybe they're not exactly heat transfer rates...

@glwagner Yeah it's written a little oddly as we can easily calculate the friction velocity squared (which I then take square roots of), but perhaps the script will be clearer if I just use uτ.

The superscript + and - indicate the top and bottom wall values, but as @glwagner has pointed out the superscript + is reserved for wall units so perhaps I need to improve my notation here as well.

[ali-ramadhan]

I thought one of the advantages of Julia is that one can use easy to read notation. The notation for Re and Nu is anything but intuitive. Should we discuss how to improve on that?

Yes there's definitely room for improvement here.

[glwagner]

@ali-ramadhan how is qw calculated?

Are you sure your definition of Nu is correct? I would interpret it this way:

Laminar solution heat flux = Θw κ / h

(Because the wall-to-wall temperature difference is 2Θw, and the wall spacing is 2h, so the laminar temperature gradient would be Θw/h.)

Turbulent heat flux = qw.

Then

Nusselt number = Turbulent heat flux / laminar heat flux = qw * h / (Θw * κ).

I'll review the notation in the pull request.

[ali-ramadhan]

@glwagner Quite possibly.

qw is calculated from equation (16) of Vreugdenhil & Taylor (2018)

qw = (friction velocity) * (friction temperature)

[glwagner]

Thanks @ali-ramadhan.

For the record, this means that the wall heat flux is define in a rational, sane manner, ie:

qw = κ | ∂y Θ |_wall |

aka, the (absolute value of the) heat flux out of the wall.

(In V & T, y is the wall-normal coordinate.)

@ali-ramadhan technically we should be using the total diffusivity (LES + molecular) to calculate qw (because these new coarse simulations do not resolve the wall boundary layer). Since the SGS diffusivity is quite small next to the wall (some fraction of the molecular diffusivity), this correction may not be huge.

[rafferrari]

Terrible notation in my humble opinion. w stands for vertical velocity. Why is the advective heat flux qw instead of w theta?

Raffaele

On Sep 6, 2019, at 9:58 AM, Gregory L. Wagner notifications@github.com wrote:

Thanks @ali-ramadhan https://github.com/ali-ramadhan.

For the record, this means that its

qw = κ ∂y Θ |_wall

aka, the heat flux out of the wall.

(In V & T, y is the wall-normal coordinate.)

@ali-ramadhan https://github.com/ali-ramadhan technically we should be using the total diffusivity (LES + molecular) to calculate qw (because these new coarse simulations do not resolve the wall boundary layer). Since the SGS diffusivity is quite small (some fraction of the molecular diffusivity), this correction may not be huge.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/climate-machine/Oceananigans.jl/issues/310?email_source=notifications&email_token=AK24ROM5GSG5ZSWBZ5KICRTQIJOXJA5CNFSM4H7PMRLKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD6C5R3Q#issuecomment-528865518, or mute the thread https://github.com/notifications/unsubscribe-auth/AK24ROM3CIZRZH5K2POAN73QIJOXJANCNFSM4H7PMRLA.

[glwagner]

@rafferrari I agree its confusing for us --- it comes from engineering, where the vertical velocity is "u₂"...

These simulations run to steady-state, in which case the horizontally averaged total vertical temperature flux is constant at every level. Thus it is sufficient to calculate the heat flux at the wall.

Technically speaking, the Nusselt number is the ratio between turbulent heat flux and laminar heat flux. In situations with strong turbulence of course the heat flux is dominated by the advective component... ? In LES we need to be somewhat delicate, because the total heat flux includes "resolved-advective" and "unresolved-advective" (subgrid scale) components. Luckily in equilibrated problems these distinctions are moot.

[ali-ramadhan]

Terrible notation in my humble opinion. w stands for vertical velocity. Why is the advective heat flux qw instead of w theta?

We should improve it before merging the test and polishing it off. @glwagner also suggested agreeing on some common notation for our papers.

Following the V&T paper it's q_w (subscript w for wall) but there is no subscript w in Unicode (or if there is, it's not available in Julia) so I ended up writing qw but maybe q_w, qʷ, or even q_wall would be better notation.

[glwagner]

+1 for q_wall for the purposes of this verification example.

[ali-ramadhan]

So I reran the simulations with Ri = 0.01 and Ri = 0.04 to compare with Cat's coarse uniform grid DIABLO simulations.

Some comments:

1. Horizontally averaged velocity and temperature profiles look qualitatively correct. Oceananigans Ri = 0.01 velocity curve could be smoother but perhaps we just need higher resolution.
2. ν and κ profiles are different. Oceananigans' magnitude is lower and shape is a bit different. Perhaps this is expected due to differences in the numerical method. Perhaps the Poincaré "constant" needs to be set to a suitable value for our finite volume numerics. Maybe we don't expect exact agreement here.
3. Most importantly, the frictional Reynolds (Re) and Nusselt (Nu) numbers match more closely now. Re values are within ~10% of each other. Oceananigans Nu values are ~20% smaller. Probably higher resolution is required to match more closely with Diablo and DNS?
4. Oceananigans uses z for the vertical direction while DIABLO uses y.