wrf-model / WRF

The official repository for the Weather Research and Forecasting (WRF) model
Other
1.18k stars 658 forks source link

Incorrect eddy diffusivity coefficient for the TKE diffusion term associate with namelist options km_opt = 2 and km_opt = 5 #2026

Closed brankokosovic closed 2 months ago

brankokosovic commented 3 months ago

Describe the bug When prognostic TKE turbulence closure is used and horizontal and vertical diffusions of TKE are computed, the scalar (i.e. heat) eddy diffusion coefficient is used instead of the momentum eddy diffusion coefficient. This results in a factor of three greater eddy diffusion coefficient for TKE under convective conditions compared to what it should be.

To Reproduce Steps to reproduce the behavior: Currently, in module_diffusion_em.F, in the subroutine horizontal_diffusion_2 the call to horizontal_diffusion_s for TKE is:

CALL horizontal_diffusion_s  ( tke_tendf(ims,kms,jms),                 &
                               config_flags,                           &
                               tke(ims,kms,jms),                       &
                               msftx, msfty, msfux, msfuy,             &
                               msfvx, msfvy, xkhh, rdx, rdy,           &
                               fnm, fnp, cf1, cf2, cf3,                &
                               zx, zy, rdz, rdzw, dnw, dn, rho,        &
                               .true.,                                 &
                               ids, ide, jds, jde, kds, kde,           &
                               ims, ime, jms, jme, kms, kme,           &
                               its, ite, jts, jte, kts, kte           )

Here xkhh is the horizontal eddy diffusion coefficient for heat, while it should be the horizontal eddy diffusion coefficient for momentum, i.e., xkmh, as follows:

CALL horizontal_diffusion_s  ( tke_tendf(ims,kms,jms),                 &
                               config_flags,                           &
                               tke(ims,kms,jms),                       &
                               msftx, msfty, msfux, msfuy,             &
                               msfvx, msfvy, xkmh, rdx, rdy,           &
                               fnm, fnp, cf1, cf2, cf3,                &
                               zx, zy, rdz, rdzw, dnw, dn, rho,        &
                               .true.,                                 &
                               ids, ide, jds, jde, kds, kde,           &
                               ims, ime, jms, jme, kms, kme,           &
                               its, ite, jts, jte, kts, kte           )

Similarly in the call to vertical_diffusion_s for TKE the eddy diffusivity coefficient should be xkmv instead of xkhv.

Expected behavior Currently, under convective conditions, when prognostic TKE turbulence closure is used and namelist parameters is km_opt = 2 the eddy diffusivity for TKE is three times larger compared to what it should be, which means that the diffusion of TKE will be reduced by a factor of three.

Screenshots If applicable, add screenshots to help explain your problem.

Attachments If applicable, attach supporting files: namelist.input, rsl.* files, etc.

Additional context Add any other context about the problem here, such as: This bug has been in WRF at least since version 3.0. The WRF Technical Note does not explicitly distinguish between eddy diffusivities for heat and momentum.

dudhia commented 3 months ago

Thanks, it was not obvious which should be used. Is there a reference?

On Thu, Mar 21, 2024 at 1:18 PM Branko Kosovic @.***> wrote:

Describe the bug When prognostic TKE turbulence closure is used and horizontal and vertical diffusions of TKE are computed, the scalar (i.e. heat) eddy diffusion coefficient is used instead of the momentum eddy diffusion coefficient. This results in a factor of three greater eddy diffusion coefficient for TKE under convective conditions compared to what it should be.

To Reproduce Steps to reproduce the behavior: Currently, in module_diffusion_em.F, in the subroutine horizontal_diffusion_2 the call to horizontal_diffusion_s for TKE is:

CALL horizontal_diffusion_s ( tke_tendf(ims,kms,jms), & config_flags, & tke(ims,kms,jms), & msftx, msfty, msfux, msfuy, & msfvx, msfvy, xkhh, rdx, rdy, & fnm, fnp, cf1, cf2, cf3, & zx, zy, rdz, rdzw, dnw, dn, rho, & .true., & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte )

Here xkhh is the horizontal eddy diffusion coefficient for heat, while it should be the horizontal eddy diffusion coefficient for momentum, i.e., xkmh, as follows:

CALL horizontal_diffusion_s ( tke_tendf(ims,kms,jms), & config_flags, & tke(ims,kms,jms), & msftx, msfty, msfux, msfuy, & msfvx, msfvy, xkmh, rdx, rdy, & fnm, fnp, cf1, cf2, cf3, & zx, zy, rdz, rdzw, dnw, dn, rho, & .true., & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte )

Similarly in the call to vertical_diffusion_s for TKE the eddy diffusivity coefficient should be xkmv instead of xkhv.

Expected behavior Currently, under convective conditions, when prognostic TKE turbulence closure is used and namelist parameters is km_opt = 2 the eddy diffusivity for TKE is three times larger compared to what it should be, which means that the diffusion of TKE will be reduced by a factor of three.

Screenshots If applicable, add screenshots to help explain your problem.

Attachments If applicable, attach supporting files: namelist.input, rsl.* files, etc.

Additional context Add any other context about the problem here, such as: This bug has been in WRF at least since version 3.0. The WRF Technical Note does not explicitly distinguish between eddy diffusivities for heat and momentum.

— Reply to this email directly, view it on GitHub https://github.com/wrf-model/WRF/issues/2026, or unsubscribe https://github.com/notifications/unsubscribe-auth/AEIZ77FIB4FDBA4F3T4OUZTYZMXBDAVCNFSM6AAAAABFCBZTDCVHI2DSMVQWIX3LMV43ASLTON2WKOZSGIYDAOJXGUZTQMI . You are receiving this because you are subscribed to this thread.Message ID: @.***>

brankokosovic commented 3 months ago

Relevant references are Lilly (1966, http://dx.doi.org/10.5065/D62R3PMM, "The Representation of Small-scale Turbulence in Numerical Simulation Experiments") equation (4.2) and Moeng (1984, https://doi-org.cuucar.idm.oclc.org/10.1175/1520-0469(1984)041%3C2052:ALESMF%3E2.0.CO;2) equation (15). The prognostic TKE equation is derived from the subgrid momentum equation and therefore the diffusion coefficient should be the diffusion coefficient for momentum.

brankokosovic commented 3 months ago

According to Lilly (1966) the diffusion coefficient is K_m. However, according to Moeng (1984), based on Deardorff (1980, Boundary-Layer Meteorology, 18, 495–527) the eddy diffusion coefficient should be: 2 K_m. So, in this case the difference under convective conditions would be a factor of 1.5 and not 3 as previously stated.

dudhia commented 3 months ago

Thanks, I am checking what MMM people think.

On Thu, Mar 21, 2024 at 2:24 PM Branko Kosovic @.***> wrote:

According to Lilly (1966) the diffusion coefficient is K_m. However, according to Moeng (1984), based on Deardorff (1980, Boundary-Layer Meteorology, 18, 495–527) the eddy diffusion coefficient should be: 2 K_m. So, in this case the difference under convective conditions would be a factor of 1.5 and not 3 as previously stated.

— Reply to this email directly, view it on GitHub https://github.com/wrf-model/WRF/issues/2026#issuecomment-2013661090, or unsubscribe https://github.com/notifications/unsubscribe-auth/AEIZ77B5A2O2CSSUXNTQ4YTYZM6XZAVCNFSM6AAAAABFCBZTDCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAMJTGY3DCMBZGA . You are receiving this because you commented.Message ID: @.***>

dudhia commented 3 months ago

George Bryan also showed me a reference with 2 Km so I think we need to do that.

On Thu, Mar 21, 2024 at 2:31 PM Jimy Dudhia @.***> wrote:

Thanks, I am checking what MMM people think.

On Thu, Mar 21, 2024 at 2:24 PM Branko Kosovic @.***> wrote:

According to Lilly (1966) the diffusion coefficient is K_m. However, according to Moeng (1984), based on Deardorff (1980, Boundary-Layer Meteorology, 18, 495–527) the eddy diffusion coefficient should be: 2 K_m. So, in this case the difference under convective conditions would be a factor of 1.5 and not 3 as previously stated.

— Reply to this email directly, view it on GitHub https://github.com/wrf-model/WRF/issues/2026#issuecomment-2013661090, or unsubscribe https://github.com/notifications/unsubscribe-auth/AEIZ77B5A2O2CSSUXNTQ4YTYZM6XZAVCNFSM6AAAAABFCBZTDCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAMJTGY3DCMBZGA . You are receiving this because you commented.Message ID: @.***>

dudhia commented 3 months ago

Turns out we already have this factor of 2 inside the diffusion code. Someone thought of it already. IF ( doing_tke ) THEN DO j = j_start, j_end DO k = kts,ktf DO i = i_start, i_end tendency(i,k,j)=tmptendf(i,k,j)+2.* & (tendency(i,k,j)-tmptendf(i,k,j)) ENDDO ENDDO ENDDO ENDIF

On Thu, Mar 28, 2024 at 12:58 PM Jimy Dudhia @.***> wrote:

George Bryan also showed me a reference with 2 Km so I think we need to do that.

On Thu, Mar 21, 2024 at 2:31 PM Jimy Dudhia @.***> wrote:

Thanks, I am checking what MMM people think.

On Thu, Mar 21, 2024 at 2:24 PM Branko Kosovic @.***> wrote:

According to Lilly (1966) the diffusion coefficient is K_m. However, according to Moeng (1984), based on Deardorff (1980, Boundary-Layer Meteorology, 18, 495–527) the eddy diffusion coefficient should be: 2 K_m. So, in this case the difference under convective conditions would be a factor of 1.5 and not 3 as previously stated.

— Reply to this email directly, view it on GitHub https://github.com/wrf-model/WRF/issues/2026#issuecomment-2013661090, or unsubscribe https://github.com/notifications/unsubscribe-auth/AEIZ77B5A2O2CSSUXNTQ4YTYZM6XZAVCNFSM6AAAAABFCBZTDCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAMJTGY3DCMBZGA . You are receiving this because you commented.Message ID: @.***>

dudhia commented 3 months ago

OK, I think your fix should still apply because it should be 2Km, not 2Kh.