I think I have stumbled across evidence of a bug in the GLS vertical mixing scheme in ROMS, at least for the k-epsilon scheme. It remains possible that I am just doing something stupid.
I find that under the k-epsilon scheme vertical stratification diffuses away much more quickly than under MY25 or LMD schemes, and faster than I would anticipate. I have also compared ROMS k-epsilon to a k-epsilon solution in GOTM, and also find that ROMS mixes away density gradients more quickly than GOTM with the same mixing scheme.
To test this, I have a run with linear equation of state and a linear innitial temperature stratification equivalent to a rho_z=-0.01 kg/m^4 in a periodic channel with no horizontal gradients in density or wind forcing; the only thing that should happen is diffusion of stratification.
I am eager to know if anyone knows of any way I could have misconfigured k-epsilon to get this excess mixing? And if not, where would one start to debug this code? Who is the expert on the GLS schemes in ROMS? This is a common scheme in coastal applications, and it would be nice to know that it is correctly implemented.
I have attached all files needed to replicate the run below; however, briefly, my three mixing schemes are activated with the following options:
Code: Select all
#define GLS_MIXING #if defined GLS_MIXING # undef KANTHA_CLAYSON # define CANUTO_A # define N2S2_HORAVG # define RI_SPLINES #endif #undef MY25_MIXING #if defined MY25_MIXING # define KANTHA_CLAYSON # define N2S2_HORAVG # define RI_SPLINES #endif #undef LMD_MIXING #ifdef LMD_MIXING # define ANA_SRFLUX # define LMD_CONVEC # define LMD_DDMIX # define LMD_SKPP # define LMD_BKPP # define LMD_NONLOCAL # define LMD_RIMIX #endif
Code: Select all
! Turbulent closure parameters. AKK_BAK == 5.0d-6 ! m2/s AKP_BAK == 5.0d-6 ! m2/s TKENU2 == 0.0d0 ! m2/s TKENU4 == 0.0d0 ! m4/s ! Generic length-scale turbulence closure parameters. GLS_P == 3.0d0 ! K-epsilon GLS_M == 1.5d0 GLS_N == -1.0d0 GLS_Kmin == 7.6d-6 GLS_Pmin == 1.0d-12 GLS_CMU0 == 0.5477d0 GLS_C1 == 1.44d0 GLS_C2 == 1.92d0 GLS_C3M == -0.4d0 GLS_C3P == 1.0d0 GLS_SIGK == 1.0d0 GLS_SIGP == 1.30d0
By day 45 you can see that the temperature stratification is almost gone in ROMS k-epsilon, but largely intact in all the other mixing schemes, including k-epsilon in GOTM.
An animation of the evolution can be found at http://oxbow.sr.unh.edu/data/bump02_compareGOTM30.mp4
As I said above, I am eager to know if anyone knows of any way I could have misconfigured k-epsilon to get this excess mixing. And if not, where would one start to debug this code? Who is the expert on the GLS schemes in ROMS?