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-2012 Leibniz University Hannover !--------------------------------------------------------------------------------! ! ! Current revisions: ! ------------------ ! ! ! Former revisions: ! ----------------- ! $Id: nudging.f90 1252 2013-11-07 08:15:25Z heinze $ ! ! 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 USE control_parameters USE cpulog USE indices USE interfaces USE pegrid IMPLICIT NONE INTEGER :: finput = 90, ierrn, k, t CHARACTER(1) :: hash REAL :: highheight, highqnudge, highptnudge, highunudge, highvnudge, & highwnudge, hightnudge REAL :: lowheight, lowqnudge, lowptnudge, lowunudge, lowvnudge, & lowwnudge, lowtnudge REAL :: 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; qnudge = 0.0; tnudge = 0.0; unudge = 0.0 vnudge = 0.0; wnudge = 0.0; timenudge = 0.0 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 message_string = 'errors in file NUDGING_DATA' CALL message( 'nudging', 'PA0366', 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 - fac ) * hightnudge unudge(k,t) = fac * lowunudge + ( 1 - fac ) * highunudge vnudge(k,t) = fac * lowvnudge + ( 1 - fac ) * highvnudge wnudge(k,t) = fac * lowwnudge + ( 1 - fac ) * highwnudge ptnudge(k,t) = fac * lowptnudge + ( 1 - fac ) * highptnudge qnudge(k,t) = fac * lowqnudge + ( 1 - fac ) * highqnudge ENDDO ENDDO rloop CLOSE ( finput ) ! !-- Prevent nudging if nudging profiles exhibt too small values !-- not used so far lptnudge = ANY( ABS( ptnudge ) > 1e-8 ) lqnudge = ANY( ABS( qnudge ) > 1e-8 ) lunudge = ANY( ABS( unudge ) > 1e-8 ) lvnudge = ANY( ABS( vnudge ) > 1e-8 ) lwnudge = ANY( ABS( wnudge ) > 1e-8 ) END SUBROUTINE init_nudge !------------------------------------------------------------------------------! ! Call for all grid points !------------------------------------------------------------------------------! SUBROUTINE nudge ( time, prog_var ) USE arrays_3d USE buoyancy_mod USE control_parameters USE cpulog USE indices USE interfaces USE pegrid USE statistics IMPLICIT NONE CHARACTER (LEN=*) :: prog_var REAL :: currtnudge, dtm, dtp, time INTEGER :: i, j, k, 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( v, 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( v, 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 USE buoyancy_mod USE control_parameters USE cpulog USE indices USE interfaces USE pegrid USE statistics IMPLICIT NONE CHARACTER (LEN=*) :: prog_var REAL :: currtnudge, dtm, dtp, time INTEGER :: i, j, k, 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