MODULE nudge_mod !--------------------------------------------------------------------------------! ! This file is part of PALM. ! ! PALM is free software: you can redistribute it and/or modify it under the terms ! of the GNU General Public License as published by the Free Software Foundation, ! either version 3 of the License, or (at your option) any later version. ! ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. ! ! You should have received a copy of the GNU General Public License along with ! PALM. If not, see . ! ! Copyright 1997-2014 Leibniz Universitaet Hannover !--------------------------------------------------------------------------------! ! ! Current revisions: ! ------------------ ! ! ! Former revisions: ! ----------------- ! $Id: nudging.f90 1356 2014-04-10 10:24:12Z witha $ ! ! 1355 2014-04-10 10:21:29Z heinze ! Error message specified. ! ! 1353 2014-04-08 15:21:23Z heinze ! REAL constants provided with KIND-attribute ! ! 1320 2014-03-20 08:40:49Z raasch ! ONLY-attribute added to USE-statements, ! kind-parameters added to all INTEGER and REAL declaration statements, ! kinds are defined in new module kinds, ! old module precision_kind is removed, ! revision history before 2012 removed, ! comment fields (!:) to be used for variable explanations added to ! all variable declaration statements ! ! 1318 2014-03-17 13:35:16Z raasch ! module interfaces removed ! ! 1268 2013-12-12 09:47:53Z heinze ! bugfix: argument of calc_mean_profile corrected ! ! 1251 2013-11-07 08:14:30Z heinze ! bugfix: calculate dtm and dtp also in vector version ! ! 1249 2013-11-06 10:45:47Z heinze ! remove call of user module ! reformatting ! ! 1241 2013-10-30 11:36:58Z heinze ! Initial revision ! ! Description: ! ------------ ! Nudges u, v, pt and q to given profiles on a relaxation timescale tnudge. ! Profiles are read in from NUDGIN_DATA. Code is based on Neggers et al. (2012) ! and also part of DALES and UCLA-LES. !--------------------------------------------------------------------------------! PRIVATE PUBLIC init_nudge, nudge SAVE INTERFACE nudge MODULE PROCEDURE nudge MODULE PROCEDURE nudge_ij END INTERFACE nudge CONTAINS SUBROUTINE init_nudge USE arrays_3d, & ONLY: ptnudge, qnudge, timenudge, tnudge, unudge, vnudge, wnudge, & zu USE control_parameters, & ONLY: dt_3d, lptnudge, lqnudge, lunudge, lvnudge, lwnudge, & message_string, ntnudge USE indices, & ONLY: nzb, nzt USE kinds IMPLICIT NONE INTEGER(iwp) :: finput = 90 !: INTEGER(iwp) :: ierrn !: INTEGER(iwp) :: k !: INTEGER(iwp) :: t !: CHARACTER(1) :: hash !: REAL(wp) :: highheight !: REAL(wp) :: highqnudge !: REAL(wp) :: highptnudge !: REAL(wp) :: highunudge !: REAL(wp) :: highvnudge !: REAL(wp) :: highwnudge !: REAL(wp) :: hightnudge !: REAL(wp) :: lowheight !: REAL(wp) :: lowqnudge !: REAL(wp) :: lowptnudge !: REAL(wp) :: lowunudge !: REAL(wp) :: lowvnudge !: REAL(wp) :: lowwnudge !: REAL(wp) :: lowtnudge !: REAL(wp) :: fac !: ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), & tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge), & vnudge(nzb:nzt+1,1:ntnudge), wnudge(nzb:nzt+1,1:ntnudge) ) ALLOCATE( timenudge(0:ntnudge) ) ptnudge = 0.0_wp; qnudge = 0.0_wp; tnudge = 0.0_wp; unudge = 0.0_wp vnudge = 0.0_wp; wnudge = 0.0_wp; timenudge = 0.0_wp t = 0 OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', & FORM='FORMATTED', IOSTAT=ierrn ) IF ( ierrn /= 0 ) THEN message_string = 'file NUDGING_DATA does not exist' CALL message( 'nudging', 'PA0365', 1, 2, 0, 6, 0 ) ENDIF ierrn = 0 rloop:DO t = t + 1 hash = "#" ierrn = 1 ! not zero ! !-- Search for the next line consisting of "# time", !-- from there onwards the profiles will be read DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(t) IF ( ierrn < 0 ) EXIT rloop ENDDO ierrn = 0 READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, & lowvnudge, lowwnudge , lowptnudge, & lowqnudge IF ( ierrn /= 0 ) THEN message_string = 'errors in file NUDGING_DATA' CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 ) ENDIF ierrn = 0 READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, & highvnudge, highwnudge , highptnudge, & highqnudge IF ( ierrn /= 0 ) THEN message_string = 'errors in file NUDGING_DATA' CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 ) ENDIF DO k = nzb, nzt+1 IF ( highheight < zu(k) ) THEN lowheight = highheight lowtnudge = hightnudge lowunudge = highunudge lowvnudge = highvnudge lowwnudge = highwnudge lowptnudge = highptnudge lowqnudge = highqnudge ierrn = 0 READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , & highunudge , highvnudge , & highwnudge , highptnudge, & highqnudge IF (ierrn /= 0 ) THEN WRITE( message_string, * ) 'zu(nzt+1) = ', zu(nzt+1), 'm is ',& 'higher than the maximum height in NUDING_DATA which ', & 'is ', lowheight, 'm. Interpolation on PALM ', & 'grid is not possible.' CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 ) ENDIF ENDIF ! !-- Interpolation of prescribed profiles in space fac = ( highheight - zu(k) ) / ( highheight - lowheight ) tnudge(k,t) = fac * lowtnudge + ( 1.0_wp - fac ) * hightnudge unudge(k,t) = fac * lowunudge + ( 1.0_wp - fac ) * highunudge vnudge(k,t) = fac * lowvnudge + ( 1.0_wp - fac ) * highvnudge wnudge(k,t) = fac * lowwnudge + ( 1.0_wp - fac ) * highwnudge ptnudge(k,t) = fac * lowptnudge + ( 1.0_wp - fac ) * highptnudge qnudge(k,t) = fac * lowqnudge + ( 1.0_wp - fac ) * highqnudge ENDDO ENDDO rloop CLOSE ( finput ) ! !-- Prevent nudging if nudging profiles exhibt too small values !-- not used so far lptnudge = ANY( ABS( ptnudge ) > 1.0e-8_wp ) lqnudge = ANY( ABS( qnudge ) > 1.0e-8_wp ) lunudge = ANY( ABS( unudge ) > 1.0e-8_wp ) lvnudge = ANY( ABS( vnudge ) > 1.0e-8_wp ) lwnudge = ANY( ABS( wnudge ) > 1.0e-8_wp ) END SUBROUTINE init_nudge !------------------------------------------------------------------------------! ! Call for all grid points !------------------------------------------------------------------------------! SUBROUTINE nudge ( time, prog_var ) USE arrays_3d, & ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, & v, vnudge USE buoyancy_mod, & ONLY: calc_mean_profile USE control_parameters, & ONLY: dt_3d, message_string USE indices, & ONLY: nxl, nxr, nys, nyn, nzb, nzb_u_inner, nzt USE kinds, & ONLY: iwp, wp USE statistics, & ONLY: hom IMPLICIT NONE CHARACTER (LEN=*) :: prog_var !: REAL(wp) :: currtnudge !: REAL(wp) :: dtm !: REAL(wp) :: dtp !: REAL(wp) :: time !: INTEGER(iwp) :: i !: INTEGER(iwp) :: j !: INTEGER(iwp) :: k !: INTEGER(iwp) :: t !: t = 1 DO WHILE ( time > timenudge(t) ) t = t+1 ENDDO IF ( time /= timenudge(1) ) THEN t = t-1 ENDIF dtm = ( time - timenudge(t) ) / ( timenudge(t+1) - timenudge(t) ) dtp = ( timenudge(t+1) - time ) / ( timenudge(t+1) - timenudge(t) ) SELECT CASE ( prog_var ) CASE ( 'u' ) CALL calc_mean_profile( u, 1, 'nudging' ) DO i = nxl, nxr DO j = nys, nyn DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + & tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,1,0) - & ( unudge(k,t) * dtp + & unudge(k,t+1) * dtm ) ) / currtnudge ENDDO ENDDO ENDDO CASE ( 'v' ) CALL calc_mean_profile( v, 2, 'nudging' ) DO i = nxl, nxr DO j = nys, nyn DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + & tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,2,0) - & ( vnudge(k,t) * dtp + & vnudge(k,t+1) * dtm ) ) / currtnudge ENDDO ENDDO ENDDO CASE ( 'pt' ) CALL calc_mean_profile( pt, 4, 'nudging' ) DO i = nxl, nxr DO j = nys, nyn DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + & tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,4,0) - & ( ptnudge(k,t) * dtp + & ptnudge(k,t+1) * dtm ) ) / currtnudge ENDDO ENDDO ENDDO CASE ( 'q' ) CALL calc_mean_profile( q, 41, 'nudging' ) DO i = nxl, nxr DO j = nys, nyn DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + & tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,41,0) - & ( qnudge(k,t) * dtp + & qnudge(k,t+1) * dtm ) ) / currtnudge ENDDO ENDDO ENDDO CASE DEFAULT message_string = 'unknown prognostic variable "' // prog_var // '"' CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 ) END SELECT END SUBROUTINE nudge !------------------------------------------------------------------------------! ! Call for grid point i,j !------------------------------------------------------------------------------! SUBROUTINE nudge_ij( i, j, time, prog_var ) USE arrays_3d, & ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, & v, vnudge USE buoyancy_mod, & ONLY: calc_mean_profile USE control_parameters, & ONLY: dt_3d, message_string USE indices, & ONLY: nxl, nxr, nys, nyn, nzb, nzb_u_inner, nzt USE kinds, & ONLY: iwp, wp USE statistics, & ONLY: hom IMPLICIT NONE CHARACTER (LEN=*) :: prog_var !: REAL(wp) :: currtnudge !: REAL(wp) :: dtm !: REAL(wp) :: dtp !: REAL(wp) :: time !: INTEGER(iwp) :: i !: INTEGER(iwp) :: j !: INTEGER(iwp) :: k !: INTEGER(iwp) :: t !: t = 1 DO WHILE ( time > timenudge(t) ) t = t+1 ENDDO IF ( time /= timenudge(1) ) THEN t = t-1 ENDIF dtm = ( time - timenudge(t) ) / ( timenudge(t+1) - timenudge(t) ) dtp = ( timenudge(t+1) - time ) / ( timenudge(t+1) - timenudge(t) ) SELECT CASE ( prog_var ) CASE ( 'u' ) CALL calc_mean_profile( u, 1, 'nudging' ) DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,1,0) - & ( unudge(k,t) * dtp + & unudge(k,t+1) * dtm ) ) / currtnudge ENDDO CASE ( 'v' ) CALL calc_mean_profile( v, 2, 'nudging' ) DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,2,0) - & ( vnudge(k,t) * dtp + & vnudge(k,t+1) * dtm ) ) / currtnudge ENDDO CASE ( 'pt' ) CALL calc_mean_profile( pt, 4, 'nudging' ) DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,4,0) - & ( ptnudge(k,t) * dtp + & ptnudge(k,t+1) * dtm ) ) / currtnudge ENDDO CASE ( 'q' ) CALL calc_mean_profile( q, 41, 'nudging' ) DO k = nzb_u_inner(j,i)+1, nzt currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm ) tend(k,j,i) = tend(k,j,i) - ( hom(k,1,41,0) - & ( qnudge(k,t) * dtp + & qnudge(k,t+1) * dtm ) ) / currtnudge ENDDO CASE DEFAULT message_string = 'unknown prognostic variable "' // prog_var // '"' CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 ) END SELECT END SUBROUTINE nudge_ij END MODULE nudge_mod