Setting PRT as a function of stratification
 

Hi!

I would like to set the turbulent Prandtl number as a function of the stratification in the k-e model. One way of doing this would be to set the following in Q1.

1. To calculate ENUT in GROUND from KE and EP with my own formulas for PRT: ENUT=GRND

2. To set the buoyancy term in the equations for KE and EP in GROUND using my own values for PRT: PATCH(patchname,PHASEM,1,NX,1,NY,1,NZ,1,LSTEP) COVAL(patchname,KE,FIXFLU,GRND1) COVAL(patchname,EP,FIXLFU,GRND2)

Assuming the above is correct, my remaining problem is the following. Is the section in GROUND where I set my own formula for ENUT called for each variable, or only once per slab/timestep/whatever? Can I test on INDVAR so I get different PRT-formulas for the viscosity and the diffusivity?

I run 3.2 on SGI IRIX.

Thanks, Olof

Re: Setting PRT as a function of stratification
 

Hi OlOf,

I may be wrong, but if my recollection is correct, I don't think you can do it this way.

- for a start, I am doubting that ENUT can be set indvidually and differently for a given variable. ENUT is ENUT, in the same way as RHO is RHO.

- there is no provision for a user's coding of PRT (turbulent) in ground (but please check).

Therefore, the only way left to alter the PRT value in ground seems to be to set PRT(indvar in question) to a large number in order to cancel out the built-in turbulent diffusive coefficient.

Then, you can rewrite your own coefficient via the PRNDTL(indvar in question), for which access is definitively provided in ground. This coding should hence reflect both the laminar and turbulent contribution.

Well, at least, sthing like that.

Good luck,

Hervé.

Re: Setting PRT as a function of stratification
 

Thanks, Hervé!

I got an answer saying more or less the same thing from Michael Malin at CHAM, though a bit more elaborated. I'll try it right away.

Olof

Re: Setting PRT as a function of stratification
 

This is to record the detailed advice for the benefit of other interested users. The recommended approach for introducing stability-dependent turbulent viscosity and scalar diffusivity is as follows:

STORE(CMUM,CMUS) in the Q1 file, so as create 3d stores for the parameters CMU and CMU', as these will be a function of local stability parameters. This storage will facilitate GROUND computations and permit printout of these parameters to the RESULT and PHI files.

In Group 9 Section 5 of GROUND compute ENUT for the current IZ slab from CMUM*(KE**0.5)*EL1 (note that ENUT is NOT computed for each INDVAR, but rather once for use in all transport equations). Note that the user will have to create the do loop and arrange F-array access to ENUT, CMUM, EL1, etc.

For those scalars that require variable turbulent Prandtl number, set PRT(PHI)=1.E10 so as to eliminate the turbulent contribution to the effective exchange coefficient for scalar diffusive transport. Then set PRNDTL(PHI)=-GRND so as to inform EARTH that you will be setting the laminar diffusivity for that PHI in Group 9 Section 7 of GROUND. However, in practice you will be setting the effective diffusivity (laminar+turbulent) because we have already suppressed the turbulent contribution by setting PRT(PHI)=1.E10. This procedure is necessary because there is no facility which permits the user to set a variable PRT directly in GROUND.

In Group 9 Section 7 of GROUND compute for the current IZ slab the effective diffusivity LAMPR for the diffusive transport of PHI, i.e. (CMUS*(KE**0.5)*EL1 + ENUL/laminar Prandtl number for PHI).Note that again the user will have to create the do loop and arrange F-array access to CMUS, CMUM, EL1, ENUL etc.

The use of variable CMU in the k-e model means that you will have to be very careful what source-term linearisation (KELIN) you employ for the built-in k-e model. This is because CMU appears explicity in the k-e source terms for some KELIN options (see attached documentation of KELIN). I recommend that the usen of KELIN=3 because this will avoid the problem.

The last area to consider is the buoyant production/destruction of k and e. You can modify do loop 30 as required in subroutine gxkegb of the grex source file gxgenk. This is the built-in k-e buoyancy term option, which uses a linearisation procedure that prevents negative values of k and e in the event of strong stable stratification. If you use your own coding with explicit source terms, there is a possibility of numerical instability under these conditions.