diff --git a/src/ZSOLVER/tke_eps_sources.f90 b/src/ZSOLVER/tke_eps_sources.f90 deleted file mode 100644 index dd01a9a6dde07b6f18f937adde92f0baa5bd3d01..0000000000000000000000000000000000000000 --- a/src/ZSOLVER/tke_eps_sources.f90 +++ /dev/null @@ -1,691 +0,0 @@ -!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier -!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence -!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt -!MNH_LIC for details. version 1. -!----------------------------------------------------------------- -! ########################### - MODULE MODI_TKE_EPS_SOURCES -! ########################### -INTERFACE -! - SUBROUTINE TKE_EPS_SOURCES(KKA,KKU,KKL,KMI,PTKEM,PLM,PLEPS,PDP,PTRH, & - PRHODJ,PDZZ,PDXX,PDYY,PDZX,PDZY,PZZ, & - PTSTEP,PIMPL,PEXPL, & - HTURBLEN,HTURBDIM, & - TPFILE,OTURB_DIAG, & - PTP,PRTKES,PRTKESM, PRTHLS,PCOEF_DISS,PTR,PDISS ) -! -USE MODD_IO, ONLY: TFILEDATA -! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO -INTEGER, INTENT(IN) :: KMI ! model index number -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTKEM ! TKE at t-deltat -REAL, DIMENSION(:,:,:), INTENT(IN) :: PLM ! mixing length -REAL, DIMENSION(:,:,:), INTENT(IN) :: PLEPS ! dissipative length -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! density * grid volume -REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ,PDZX,PDZY - ! metric coefficients -REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! physical height w-pt -REAL, INTENT(IN) :: PTSTEP ! Time step -REAL, INTENT(IN) :: PEXPL, PIMPL ! Coef. temporal. disc. -CHARACTER(len=4), INTENT(IN) :: HTURBDIM ! dimensionality of the - ! turbulence scheme -CHARACTER(len=4), INTENT(IN) :: HTURBLEN ! kind of mixing length -TYPE(TFILEDATA), INTENT(INOUT):: TPFILE ! Output file -LOGICAL, INTENT(IN) :: OTURB_DIAG ! switch to write some - ! diagnostic fields in the syncronous FM-file -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PDP ! Dyn. prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTRH -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTP ! Ther. prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRTKES ! RHOD * Jacobian * - ! TKE at t+deltat -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRTKESM ! Advection source -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRTHLS ! Source of Theta_l -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCOEF_DISS ! 1/(Cph*Exner) -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PTR ! Transport prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDISS ! Dissipati prod. of TKE -! -! -! -END SUBROUTINE TKE_EPS_SOURCES -! -END INTERFACE -! -END MODULE MODI_TKE_EPS_SOURCES -! -! ################################################################## - SUBROUTINE TKE_EPS_SOURCES(KKA,KKU,KKL,KMI,PTKEM,PLM,PLEPS,PDP, & - PTRH,PRHODJ,PDZZ,PDXX,PDYY,PDZX,PDZY,PZZ, & - PTSTEP,PIMPL,PEXPL, & - HTURBLEN,HTURBDIM, & - TPFILE,OTURB_DIAG, & - PTP,PRTKES,PRTKESM, PRTHLS,PCOEF_DISS,PTR,PDISS ) -! ################################################################## -! -! -!!**** *TKE_EPS_SOURCES* - routine to compute the sources of the turbulent -!! evolutive variables: TKE and its dissipation when it is taken into -!! account. The contribution to the heating of tke dissipation is computed. -!! -!! PURPOSE -!! ------- -! The purpose of this routine is to compute the sources necessary for -! the evolution of the turbulent kinetic energy and its dissipation -! if necessary. -! -!!** METHOD -!! ------ -!! The vertical turbulent flux is computed in an off-centered -!! implicit scheme (a Crank-Nicholson type with coefficients different -!! than 0.5), which allows to vary the degree of implicitness of the -!! formulation. -!! In high resolution, the horizontal transport terms are also -!! calculated, but explicitly. -!! The evolution of the dissipation as a variable is made if -!! the parameter HTURBLEN is set equal to KEPS. The same reasoning -!! made for TKE applies. -!! -!! EXTERNAL -!! -------- -!! GX_U_M,GY_V_M,GZ_W_M -!! GX_M_U,GY_M_V : Cartesian vertical gradient operators -!! -!! MXF,MXM.MYF,MYM,MZF,MZM: Shuman functions (mean operators) -!! DZF : Shuman functions (difference operators) -!! -!! -!! IMPLICIT ARGUMENTS -!! ------------------ -!! Module MODD_CST : contains physical constants -!! -!! XG : gravity constant -!! -!! Module MODD_CTURB: contains the set of constants for -!! the turbulence scheme -!! -!! XCET,XCED : transport and dissipation cts. for the TKE -!! XCDP,XCDD,XCDT: constants from the parameterization of -!! the K-epsilon equation -!! XTKEMIN,XEPSMIN : minimum values for the TKE and its -!! dissipation -!! -!! Module MODD_PARAMETERS: -!! -!! JPVEXT -!! Module MODD_BUDGET: -!! NBUMOD : model in which budget is calculated -!! CBUTYPE : type of desired budget -!! 'CART' for cartesian box configuration -!! 'MASK' for budget zone defined by a mask -!! 'NONE' ' for no budget -!! LBU_RTKE : logical for budget of RTKE (turbulent kinetic energy) -!! .TRUE. = budget of RTKE -!! .FALSE. = no budget of RTKE -!! -!! -!! REFERENCE -!! --------- -!! Book 2 of documentation (routine TKE_EPS_SOURCES) -!! Book 1 of documentation (Chapter: Turbulence) -!! -!! AUTHOR -!! ------ -!! Joan Cuxart * INM and Meteo-France * -!! -!! MODIFICATIONS -!! ------------- -!! Original August 23, 1994 -!! Modifications: Feb 14, 1995 (J.Cuxart and J.Stein) -!! Doctorization and Optimization -!! June 29, 1995 (J.Stein) TKE budget -!! June 28, 1995 (J.Cuxart) Add LES tools -!! Modifications: February 29, 1996 (J. Stein) optimization -!! Modifications: May 6, 1996 (N. Wood) Extend some loops over -!! the outer points -!! Modifications: August 30, 1996 (P. Jabouille) calcul ZFLX at the -!! IKU level -!! October 10, 1996 (J.Stein) set Keff at t-deltat -!! Oct 8, 1996 (Cuxart,Sanchez) Var.LES: XETR_TF,XDISS_TF -!! December 20, 1996 (J.-P. Pinty) update the CALL BUDGET -!! November 24, 1997 (V. Masson) bug in <v'e> -!! removes the DO loops -!! Augu. 9, 1999 (J.Stein) TKE budget correction -!! Mar 07 2001 (V. Masson and J. Stein) remove the horizontal -!! turbulent transports of Tke computation -!! Nov 06, 2002 (V. Masson) LES budgets -!! July 20, 2003 (J.-P. Pinty P Jabouille) add the dissipative heating -!! May 2006 Remove KEPS -!! October 2009 (G. Tanguy) add ILENCH=LEN(YCOMMENT) after -!! change of YCOMMENT -!! 2012-02 Y. Seity, add possibility to run with reversed -!! vertical levels -!! 2015-01 (J. Escobar) missing get_halo(ZRES) for JPHEXT<> 1 -!! J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 -!! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O -! P. Wautelet 20/05/2019: add name argument to ADDnFIELD_ll + new ADD4DFIELD_ll subroutine -! P. Wautelet 02/2020: use the new data structures and subroutines for budgets -! -------------------------------------------------------------------------- -! -!* 0. DECLARATIONS -! ------------ -! -USE MODD_ARGSLIST_ll, ONLY: LIST_ll -use modd_budget, only: lbudget_tke, lbudget_th, NBUDGET_TKE, NBUDGET_TH, tbudgets -USE MODD_CONF -USE MODD_CST -USE MODD_CTURB -USE MODD_DIAG_IN_RUN, ONLY: LDIAG_IN_RUN, XCURRENT_TKE_DISS -use modd_field, only: tfielddata, TYPEREAL -USE MODD_IO, ONLY: TFILEDATA -USE MODD_LES -USE MODD_PARAMETERS -! -USE MODE_ARGSLIST_ll, ONLY: ADD3DFIELD_ll, CLEANLIST_ll -use mode_budget, only: Budget_store_add, Budget_store_end, Budget_store_init -USE MODE_EXCHANGE_ll, ONLY: UPDATE_HALO_ll -USE MODE_IO_FIELD_WRITE, only: IO_Field_write -#ifdef MNH_OPENACC -USE MODE_MNH_ZWORK, ONLY: MNH_MEM_GET, MNH_MEM_POSITION_PIN, MNH_MEM_RELEASE -#endif -use mode_mppdb -! -#if defined(MNH_BITREP) || defined(MNH_BITREP_OMP) -USE MODI_BITREP -#endif -USE MODI_GET_HALO -USE MODI_GRADIENT_M -USE MODI_GRADIENT_U -USE MODI_GRADIENT_V -USE MODI_GRADIENT_W -USE MODI_LES_MEAN_SUBGRID -#ifndef MNH_OPENACC -USE MODI_SHUMAN -#else -USE MODI_SHUMAN_DEVICE -#endif -USE MODI_TRIDIAG_TKE -! -IMPLICIT NONE -! -! -!* 0.1 declarations of arguments -! -! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO - -INTEGER, INTENT(IN) :: KMI ! model index number -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTKEM ! TKE at t-deltat -REAL, DIMENSION(:,:,:), INTENT(IN) :: PLM ! mixing length -REAL, DIMENSION(:,:,:), INTENT(IN) :: PLEPS ! dissipative length -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! density * grid volume -REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ,PDZX,PDZY - ! metric coefficients -REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! physical height w-pt -REAL, INTENT(IN) :: PTSTEP ! Time step -REAL, INTENT(IN) :: PEXPL, PIMPL ! Coef. temporal. disc. -CHARACTER(len=4), INTENT(IN) :: HTURBDIM ! dimensionality of the - ! turbulence scheme -CHARACTER(len=4), INTENT(IN) :: HTURBLEN ! kind of mixing length -TYPE(TFILEDATA), INTENT(INOUT):: TPFILE ! Output file -LOGICAL, INTENT(IN) :: OTURB_DIAG ! switch to write some - ! diagnostic fields in the syncronous FM-file -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PDP ! Dyn. prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTRH -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTP ! Ther. prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRTKES ! RHOD * Jacobian * - ! TKE at t+deltat -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRTHLS ! Source of Theta_l -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCOEF_DISS ! 1/(Cph*Exner) -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRTKESM ! Advection source -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PTR ! Transport prod. of TKE -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDISS ! Dissipati prod. of TKE -! -!* 0.2 declaration of local variables -! -REAL, DIMENSION(:,:,:), pointer , contiguous :: & - ZA, & ! under diagonal elements of the tri-diagonal matrix involved - ! in the temporal implicit scheme - ZRES, & ! treated variable at t+ deltat when the turbu- - ! lence is the only source of evolution added to the ones - ! considered in ZSOURCE. This variable is also used to - ! temporarily store some diagnostics stored in FM file - ZFLX, & ! horizontal or vertical flux of the treated variable - ZSOURCE, & ! source of evolution for the treated variable - ZKEFF ! effectif diffusion coeff = LT * SQRT( TKE ) -INTEGER :: IIB,IIE,IJB,IJE,IKB,IKE - ! Index values for the Beginning and End - ! mass points of the domain -INTEGER :: IIU,IJU,IKU ! array size in the 3 dimensions -! -TYPE(LIST_ll), POINTER :: TZFIELDDISS_ll ! list of fields to exchange -INTEGER :: IINFO_ll ! return code of parallel routine -TYPE(TFIELDDATA) :: TZFIELD -! -#ifdef MNH_OPENACC -REAL, DIMENSION(:,:,:), pointer , contiguous :: ZTMP1_DEVICE, ZTMP2_DEVICE, ZTMP3_DEVICE, ZTMP4_DEVICE -#endif -! -INTEGER :: JIU,JJU,JKU -INTEGER :: JI,JJ,JK -LOGICAL :: GTURBDIM_3DIM -!---------------------------------------------------------------------------- -! -GTURBDIM_3DIM = HTURBDIM=='3DIM' -! -!$acc data present( PTKEM, PLM, PLEPS, PDP, PTRH, & -!$acc & PRHODJ, PDZZ, PDXX, PDYY, PDZX, PDZY, PZZ, & -!$acc & PTP, PRTKES, PRTKESM, PRTHLS, PCOEF_DISS, & -!$acc & PTR, PDISS ) - -if ( mppdb_initialized ) then - !Check all in arrays - call Mppdb_check( ptkem, "Tke_eps_sources beg:ptkem" ) - call Mppdb_check( plm, "Tke_eps_sources beg:plm" ) - call Mppdb_check( pleps, "Tke_eps_sources beg:pleps" ) - call Mppdb_check( prhodj, "Tke_eps_sources beg:prhodj" ) - call Mppdb_check( pdxx, "Tke_eps_sources beg:pdxx" ) - call Mppdb_check( pdyy, "Tke_eps_sources beg:pdyy" ) - call Mppdb_check( pdzz, "Tke_eps_sources beg:pdzz" ) - call Mppdb_check( pdzx, "Tke_eps_sources beg:pdzx" ) - call Mppdb_check( pdzy, "Tke_eps_sources beg:pdzy" ) - call Mppdb_check( pzz, "Tke_eps_sources beg:pzz" ) - call Mppdb_check( ptrh, "Tke_eps_sources beg:ptrh" ) - call Mppdb_check( ptp, "Tke_eps_sources beg:ptp" ) - call Mppdb_check( pcoef_diss, "Tke_eps_sources beg:pcoef_diss" ) - call Mppdb_check( prtkesm, "Tke_eps_sources beg:prtkesm" ) - !check all inout arrays - call Mppdb_check( pdp, "Tke_eps_sources beg:pdp" ) - call Mppdb_check( prtkes, "Tke_eps_sources beg:prtkes" ) - call Mppdb_check( prthls, "Tke_eps_sources beg:prthls" ) -end if - -JIU = size(ptkem, 1 ) -JJU = size(ptkem, 2 ) -JKU = size(ptkem, 3 ) - -#ifndef MNH_OPENACC -allocate( za (JIU,JJU,JKU ) ) -allocate( zres (JIU,JJU,JKU ) ) -allocate( zflx (JIU,JJU,JKU ) ) -allocate( zsource(JIU,JJU,JKU ) ) -allocate( zkeff (JIU,JJU,JKU ) ) -#else -!Pin positions in the pools of MNH memory -CALL MNH_MEM_POSITION_PIN() - -CALL MNH_MEM_GET( za ,JIU,JJU,JKU) -CALL MNH_MEM_GET( zres ,JIU,JJU,JKU) -CALL MNH_MEM_GET( zflx ,JIU,JJU,JKU) -CALL MNH_MEM_GET( zsource,JIU,JJU,JKU) -CALL MNH_MEM_GET( zkeff ,JIU,JJU,JKU) - -CALL MNH_MEM_GET( ztmp1_device,JIU,JJU,JKU) -CALL MNH_MEM_GET( ztmp2_device,JIU,JJU,JKU) -CALL MNH_MEM_GET( ztmp3_device,JIU,JJU,JKU) -CALL MNH_MEM_GET( ztmp4_device,JIU,JJU,JKU) -#endif - -!$acc data present( ZA, ZRES, ZFLX, ZSOURCE, ZKEFF, ZTMP1_DEVICE, ZTMP2_DEVICE, ZTMP3_DEVICE, ZTMP4_DEVICE ) - -NULLIFY(TZFIELDDISS_ll) -! -!* 1. PRELIMINARY COMPUTATIONS -! ------------------------ -! -CALL GET_INDICE_ll (IIB,IJB,IIE,IJE) -!$acc kernels -IIU=SIZE(PTKEM,1) -IJU=SIZE(PTKEM,2) -IKB=KKA+JPVEXT_TURB*KKL -IKE=KKU-JPVEXT_TURB*KKL -! -! compute the effective diffusion coefficient at the mass point -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) -ZKEFF(:,:,:) = PLM(:,:,:) * SQRT(PTKEM(:,:,:)) -#else -#ifdef MNH_COMPILER_NVHPC -!$acc loop independent collapse(3) -#endif -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - ZKEFF(JI,JJ,JK) = PLM(JI,JJ,JK) * BR_POW(PTKEM(JI,JJ,JK),0.5) -END DO -#endif -! -!$acc end kernels -if (lbudget_th) call Budget_store_init( tbudgets(NBUDGET_TH), 'DISSH', prthls(:, :, :) ) -!$acc kernels -!---------------------------------------------------------------------------- -! -!* 2. TKE EQUATION -! ------------ -! -!* 2.1 Horizontal turbulent explicit transport -! -! -! Complete the sources of TKE with the horizontal turbulent explicit transport -! -IF (GTURBDIM_3DIM) THEN - PTR(:,:,:) = PTRH(:,:,:) -ELSE - PTR(:,:,:) = 0. -END IF -! -! -! -!* 2.2 Explicit TKE sources except horizontal turbulent transport -! -! -! extrapolate the dynamic production with a 1/Z law from its value at the -! W(IKB+1) value stored in PDP(IKB) to the mass localization tke(IKB) -PDP(:,:,IKB) = PDP(:,:,IKB) * (1. + PDZZ(:,:,IKB+KKL)/PDZZ(:,:,IKB)) -! -! Compute the source terms for TKE: ( ADVECtion + NUMerical DIFFusion + ..) -! + (Dynamical Production) + (Thermal Production) - (dissipation) -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) -ZFLX(:,:,:) = XCED * SQRT(PTKEM(:,:,:)) / PLEPS(:,:,:) -#else -#ifdef MNH_COMPILER_NVHPC -!$acc loop independent collapse(3) -#endif -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - ZFLX(JI,JJ,JK) = XCED * BR_POW(PTKEM(JI,JJ,JK),0.5) / PLEPS(JI,JJ,JK) -END DO -#endif -#ifdef MNH_COMPILER_NVHPC -!$acc loop independent collapse(3) -#endif -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - ZSOURCE(JI,JJ,JK) = ( PRTKES(JI,JJ,JK) + PRTKESM(JI,JJ,JK) ) / PRHODJ(JI,JJ,JK) & - - PTKEM(JI,JJ,JK) / PTSTEP & - + PDP(JI,JJ,JK) + PTP(JI,JJ,JK) + PTR(JI,JJ,JK) - PEXPL * ZFLX(JI,JJ,JK) * PTKEM(JI,JJ,JK) -END DO !CONCURRENT -!$acc end kernels -! -!* 2.2 implicit vertical TKE transport -! -! -! Compute the vector giving the elements just under the diagonal for the -! matrix inverted in TRIDIAG -! -#ifndef MNH_OPENACC -ZA(:,:,:) = - PTSTEP * XCET * & -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) - MZM(ZKEFF) * MZM(PRHODJ) / PDZZ(:,:,:)**2 -#else - MZM(ZKEFF) * MZM(PRHODJ) / BR_P2(PDZZ(:,:,:)) -#endif -#else -CALL MZM_DEVICE(ZKEFF, ZTMP1_DEVICE) !Warning: re-used later -CALL MZM_DEVICE(PRHODJ,ZTMP2_DEVICE) !Warning: re-used later -!$acc kernels ! present(ZA) -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) -ZA(:,:,:) = - PTSTEP * XCET * ZTMP1_DEVICE(:,:,:) * ZTMP2_DEVICE(:,:,:) / PDZZ(:,:,:)**2 -#else -!$acc_nv loop independent collapse(3) -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - ZA(JI,JJ,JK) = - PTSTEP * XCET * ZTMP1_DEVICE(JI,JJ,JK) * ZTMP2_DEVICE(JI,JJ,JK) / BR_P2(PDZZ(JI,JJ,JK)) -END DO !CONCURRENT -#endif -!$acc end kernels -#endif -! -! Compute TKE at time t+deltat: ( stored in ZRES ) -! -#ifndef MNH_OPENACC -CALL TRIDIAG_TKE(KKA,KKU,KKL,PTKEM,ZA,PTSTEP,PEXPL,PIMPL,PRHODJ,& - & ZSOURCE,PTSTEP*ZFLX,ZRES) -CALL GET_HALO(ZRES) -#else -!$acc kernels -#ifdef MNH_COMPILER_NVHPC -!$acc loop independent collapse(3) -#endif -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - ZTMP3_DEVICE(JI,JJ,JK) = PTSTEP*ZFLX(JI,JJ,JK) -END DO !CONCURRENT -!$acc end kernels -CALL TRIDIAG_TKE(KKA,KKU,KKL,PTKEM,ZA,PTSTEP,PEXPL,PIMPL,PRHODJ,& - & ZSOURCE,ZTMP3_DEVICE,ZRES) -! acc update self(ZRES) -CALL GET_HALO_D(ZRES) -! acc update device(ZRES) -#endif -! -!* diagnose the dissipation -! -IF (LDIAG_IN_RUN) THEN -!$acc data copyout(XCURRENT_TKE_DISS) - XCURRENT_TKE_DISS = ZFLX(:,:,:) * PTKEM(:,:,:) & - *(PEXPL*PTKEM(:,:,:) + PIMPL*ZRES(:,:,:)) -!$acc end data - CALL ADD3DFIELD_ll( TZFIELDDISS_ll, XCURRENT_TKE_DISS, 'TKE_EPS_SOURCES::XCURRENT_TKE_DISS' ) - CALL UPDATE_HALO_ll(TZFIELDDISS_ll,IINFO_ll) - CALL CLEANLIST_ll(TZFIELDDISS_ll) -ENDIF -! -! TKE must be greater than its minimum value -! -! CL : Now done at the end of the time step in ADVECTION_METSV -!GTKENEG = ZRES <= XTKEMIN -!WHERE ( GTKENEG ) -! ZRES = XTKEMIN -!END WHERE -! -IF ( LLES_CALL .OR. & - (OTURB_DIAG .AND. tpfile%lopened) ) THEN -! -! Compute the cartesian vertical flux of TKE in ZFLX -! -#ifndef MNH_OPENACC - ZFLX(:,:,:) = - XCET * MZM(ZKEFF) * & - DZM(PIMPL * ZRES + PEXPL * PTKEM ) / PDZZ(:,:,:) -#else -!$acc kernels - ZTMP3_DEVICE(:,:,:) = PIMPL * ZRES(:,:,:) + PEXPL * PTKEM(:,:,:) -!$acc end kernels - CALL DZM_DEVICE( ZTMP3_DEVICE, ZTMP4_DEVICE ) -!$acc kernels - ZFLX(:,:,:) = - XCET * ZTMP1_DEVICE(:,:,:) * ZTMP4_DEVICE(:,:,:) / PDZZ(:,:,:) !Re-use of ZTMP1_DEVICE -#endif -! - ZFLX(:,:,IKB) = 0. - ZFLX(:,:,KKA) = 0. -! -! Compute the whole turbulent TRansport of TKE: -! -#ifndef MNH_OPENACC - PTR(:,:,:)= PTR - DZF( MZM(PRHODJ) * ZFLX / PDZZ ) /PRHODJ -#else - ZTMP1_DEVICE(:,:,:) = ZTMP2_DEVICE(:,:,:) * ZFLX(:,:,:) / PDZZ(:,:,:) !Re-use of ZTMP2_DEVICE -!$acc end kernels - CALL DZF_DEVICE( ZTMP1_DEVICE, ZTMP2_DEVICE ) -!$acc kernels - PTR(:,:,:)= PTR(:,:,:) - ZTMP2_DEVICE(:,:,:) / PRHODJ(:,:,:) -!$acc end kernels -#endif -! -! Storage in the LES configuration -! - IF (LLES_CALL) THEN -#ifndef MNH_OPENACC - CALL LES_MEAN_SUBGRID( MZF(ZFLX), X_LES_SUBGRID_WTke ) - CALL LES_MEAN_SUBGRID( -PTR, X_LES_SUBGRID_ddz_WTke ) -#else -!$acc data copy(X_LES_SUBGRID_WTke,X_LES_SUBGRID_ddz_WTke) - CALL MZF_DEVICE( ZFLX, ZTMP1_DEVICE ) - CALL LES_MEAN_SUBGRID( ZTMP1_DEVICE, X_LES_SUBGRID_WTke ) -!$acc kernels - ZTMP1_DEVICE(:,:,:) = -PTR(:,:,:) -!$acc end kernels - CALL LES_MEAN_SUBGRID( ZTMP1_DEVICE, X_LES_SUBGRID_ddz_WTke ) -!$acc end data -#endif - END IF -! -END IF -! -!* 2.4 stores the explicit sources for budget purposes -! -if (lbudget_tke) then - ! Dynamical production - call Budget_store_add( tbudgets(NBUDGET_TKE), 'DP', pdp(:, :, :) * prhodj(:, :, :) ) - ! Thermal production - call Budget_store_add( tbudgets(NBUDGET_TKE), 'TP', ptp(:, :, :) * prhodj(:, :, :) ) - ! Dissipation - call Budget_store_add( tbudgets(NBUDGET_TKE), 'DISS', -xced * sqrt( ptkem(:, :, :) ) / pleps(:, :, :) & - * ( pexpl * ptkem(:, :, :) + pimpl * zres(:, :, :) ) * prhodj(:, :, :) ) -end if -! -!* 2.5 computes the final RTKE and stores the whole turbulent transport -! with the removal of the advection part - -if (lbudget_tke) then - !Store the previous source terms in prtkes before initializing the next one - !$acc kernels present(ZRES) - PRTKES(:,:,:) = PRTKES(:,:,:) + PRHODJ(:,:,:) * & - ( PDP(:,:,:) + PTP(:,:,:) & - - XCED * SQRT(PTKEM(:,:,:)) / PLEPS(:,:,:) * ( PEXPL*PTKEM(:,:,:) + PIMPL*ZRES(:,:,:) ) ) - !$acc end kernels - call Budget_store_init( tbudgets(NBUDGET_TKE), 'TR', prtkes(:, :, :) ) -end if - -!$acc kernels -DO CONCURRENT ( JI=1:JIU,JJ=1:JJU,JK=1:JKU) - PRTKES(JI,JJ,JK) = ZRES(JI,JJ,JK) * PRHODJ(JI,JJ,JK) / PTSTEP - PRTKESM(JI,JJ,JK) -ENDDO -!$acc end kernels -! -! stores the whole turbulent transport -! -if (lbudget_tke) call Budget_store_end( tbudgets(NBUDGET_TKE), 'TR', prtkes(:, :, :) ) - -!---------------------------------------------------------------------------- -! -!* 3. COMPUTE THE DISSIPATIVE HEATING -! ------------------------------- -! -!$acc kernels -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) -!dir$ concurrent -PRTHLS(:,:,:) = PRTHLS(:,:,:) + XCED * SQRT(PTKEM(:,:,:)) / PLEPS(:,:,:) * & - (PEXPL*PTKEM(:,:,:) + PIMPL*ZRES(:,:,:)) * PRHODJ(:,:,:) * PCOEF_DISS(:,:,:) -#else -DO CONCURRENT (JI=1:JIU, JJ=1:JJU ,JK=1:JKU) - PRTHLS(JI,JJ,JK) = PRTHLS(JI,JJ,JK) + XCED * BR_POW(PTKEM(JI,JJ,JK),0.5) / PLEPS(JI,JJ,JK) * & - (PEXPL*PTKEM(JI,JJ,JK) + PIMPL*ZRES(JI,JJ,JK)) * PRHODJ(JI,JJ,JK) * PCOEF_DISS(JI,JJ,JK) -ENDDO -#endif -!$acc end kernels - -if (lbudget_th) call Budget_store_end( tbudgets(NBUDGET_TH), 'DISSH', prthls(:, :, :) ) - -!---------------------------------------------------------------------------- -! -!* 4. STORES SOME DIAGNOSTICS -! ----------------------- -! -!$acc kernels -#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP) -PDISS(:,:,:) = -XCED * (PTKEM(:,:,:)**1.5) / PLEPS(:,:,:) -#else -DO CONCURRENT (JI=1:JIU, JJ=1:JJU, JK=1:JKU) - PDISS(JI,JJ,JK) = -XCED * BR_POW(PTKEM(JI,JJ,JK),1.5) / PLEPS(JI,JJ,JK) -ENDDO -#endif -!$acc end kernels -! -IF ( OTURB_DIAG .AND. tpfile%lopened ) THEN -!$acc update self(PDP,PTP,PTR,PDISS) -! -! stores the dynamic production -! - TZFIELD%CMNHNAME = 'DP' - TZFIELD%CSTDNAME = '' - TZFIELD%CLONGNAME = 'DP' - TZFIELD%CUNITS = 'm2 s-3' - TZFIELD%CDIR = 'XY' - TZFIELD%CCOMMENT = 'X_Y_Z_DP' - TZFIELD%NGRID = 1 - TZFIELD%NTYPE = TYPEREAL - TZFIELD%NDIMS = 3 - TZFIELD%LTIMEDEP = .TRUE. - CALL IO_Field_write(TPFILE,TZFIELD,PDP) -! -! stores the thermal production -! - TZFIELD%CMNHNAME = 'TP' - TZFIELD%CSTDNAME = '' - TZFIELD%CLONGNAME = 'TP' - TZFIELD%CUNITS = 'm2 s-3' - TZFIELD%CDIR = 'XY' - TZFIELD%CCOMMENT = 'X_Y_Z_TP' - TZFIELD%NGRID = 1 - TZFIELD%NTYPE = TYPEREAL - TZFIELD%NDIMS = 3 - TZFIELD%LTIMEDEP = .TRUE. - CALL IO_Field_write(TPFILE,TZFIELD,PTP) -! -! stores the whole turbulent transport -! - TZFIELD%CMNHNAME = 'TR' - TZFIELD%CSTDNAME = '' - TZFIELD%CLONGNAME = 'TR' - TZFIELD%CUNITS = 'm2 s-3' - TZFIELD%CDIR = 'XY' - TZFIELD%CCOMMENT = 'X_Y_Z_TR' - TZFIELD%NGRID = 1 - TZFIELD%NTYPE = TYPEREAL - TZFIELD%NDIMS = 3 - TZFIELD%LTIMEDEP = .TRUE. - CALL IO_Field_write(TPFILE,TZFIELD,PTR) -! -! stores the dissipation of TKE -! - TZFIELD%CMNHNAME = 'DISS' - TZFIELD%CSTDNAME = '' - TZFIELD%CLONGNAME = 'DISS' - TZFIELD%CUNITS = 'm2 s-3' - TZFIELD%CDIR = 'XY' - TZFIELD%CCOMMENT = 'X_Y_Z_DISS' - TZFIELD%NGRID = 1 - TZFIELD%NTYPE = TYPEREAL - TZFIELD%NDIMS = 3 - TZFIELD%LTIMEDEP = .TRUE. - CALL IO_Field_write(TPFILE,TZFIELD,PDISS) -END IF -! -! Storage in the LES configuration of the Dynamic Production of TKE and -! the dissipation of TKE -! -IF (LLES_CALL ) THEN -!$acc data copy(X_LES_SUBGRID_DISS_Tke) - CALL LES_MEAN_SUBGRID( PDISS, X_LES_SUBGRID_DISS_Tke ) -!$acc end data -END IF - -if ( mppdb_initialized ) then - !check all inout arrays - call Mppdb_check( pdp, "Tke_eps_sources end:pdp" ) - call Mppdb_check( prtkes, "Tke_eps_sources end:prtkes" ) - call Mppdb_check( prthls, "Tke_eps_sources end:prthls" ) - !check all out arrays - call Mppdb_check( ptr, "Tke_eps_sources end:ptr" ) - call Mppdb_check( pdiss, "Tke_eps_sources end:pdiss" ) -end if - -!$acc end data - -#ifndef MNH_OPENACC -deallocate( ZA, ZRES, ZFLX, ZSOURCE, ZKEFF ) -#else -!Release all memory allocated with MNH_MEM_GET calls since last call to MNH_MEM_POSITION_PIN -CALL MNH_MEM_RELEASE() -#endif - -!$acc end data - -END SUBROUTINE TKE_EPS_SOURCES