Home

Context Navigation

source: palm/trunk/SOURCE/nudging.f90 @ 1355

Last change on this file since 1355 was 1355, checked in by heinze, 11 years ago
Error message specified.
Property svn:keywords set to `Id`
File size: 15.3 KB

Rev	Line
[1239]	1	MODULE nudge_mod
	2
	3	!--------------------------------------------------------------------------------!
	4	! This file is part of PALM.
	5	!
	6	! PALM is free software: you can redistribute it and/or modify it under the terms
	7	! of the GNU General Public License as published by the Free Software Foundation,
	8	! either version 3 of the License, or (at your option) any later version.
	9	!
	10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
	11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
	12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
	13	!
	14	! You should have received a copy of the GNU General Public License along with
	15	! PALM. If not, see <http://www.gnu.org/licenses/>.
	16	!
[1310]	17	! Copyright 1997-2014 Leibniz Universitaet Hannover
[1239]	18	!--------------------------------------------------------------------------------!
	19	!
	20	! Current revisions:
	21	! ------------------
[1355]	22	! Error message specified.
[1354]	23	!
[1321]	24	! Former revisions:
	25	! -----------------
	26	! $Id: nudging.f90 1355 2014-04-10 10:21:29Z heinze $
	27	!
[1354]	28	! 1353 2014-04-08 15:21:23Z heinze
	29	! REAL constants provided with KIND-attribute
	30	!
[1321]	31	! 1320 2014-03-20 08:40:49Z raasch
[1320]	32	! ONLY-attribute added to USE-statements,
	33	! kind-parameters added to all INTEGER and REAL declaration statements,
	34	! kinds are defined in new module kinds,
	35	! old module precision_kind is removed,
	36	! revision history before 2012 removed,
	37	! comment fields (!:) to be used for variable explanations added to
	38	! all variable declaration statements
[1239]	39	!
[1319]	40	! 1318 2014-03-17 13:35:16Z raasch
	41	! module interfaces removed
	42	!
[1269]	43	! 1268 2013-12-12 09:47:53Z heinze
	44	! bugfix: argument of calc_mean_profile corrected
	45	!
[1252]	46	! 1251 2013-11-07 08:14:30Z heinze
	47	! bugfix: calculate dtm and dtp also in vector version
	48	!
[1250]	49	! 1249 2013-11-06 10:45:47Z heinze
	50	! remove call of user module
	51	! reformatting
	52	!
[1242]	53	! 1241 2013-10-30 11:36:58Z heinze
	54	! Initial revision
[1239]	55	!
	56	! Description:
	57	! ------------
	58	! Nudges u, v, pt and q to given profiles on a relaxation timescale tnudge.
	59	! Profiles are read in from NUDGIN_DATA. Code is based on Neggers et al. (2012)
	60	! and also part of DALES and UCLA-LES.
	61	!--------------------------------------------------------------------------------!
	62
	63	PRIVATE
	64	PUBLIC init_nudge, nudge
	65	SAVE
	66
	67	INTERFACE nudge
	68	MODULE PROCEDURE nudge
	69	MODULE PROCEDURE nudge_ij
	70	END INTERFACE nudge
	71
	72	CONTAINS
	73
	74	SUBROUTINE init_nudge
	75
[1320]	76	USE arrays_3d, &
	77	ONLY: ptnudge, qnudge, timenudge, tnudge, unudge, vnudge, wnudge, &
	78	zu
[1239]	79
[1320]	80	USE control_parameters, &
	81	ONLY: dt_3d, lptnudge, lqnudge, lunudge, lvnudge, lwnudge, &
	82	message_string, ntnudge
	83
	84	USE indices, &
	85	ONLY: nzb, nzt
	86
	87	USE kinds
	88
[1239]	89	IMPLICIT NONE
	90
	91
[1320]	92	INTEGER(iwp) :: finput = 90 !:
	93	INTEGER(iwp) :: ierrn !:
	94	INTEGER(iwp) :: k !:
	95	INTEGER(iwp) :: t !:
[1239]	96
[1320]	97	CHARACTER(1) :: hash !:
	98
	99	REAL(wp) :: highheight !:
	100	REAL(wp) :: highqnudge !:
	101	REAL(wp) :: highptnudge !:
	102	REAL(wp) :: highunudge !:
	103	REAL(wp) :: highvnudge !:
	104	REAL(wp) :: highwnudge !:
	105	REAL(wp) :: hightnudge !:
	106
	107	REAL(wp) :: lowheight !:
	108	REAL(wp) :: lowqnudge !:
	109	REAL(wp) :: lowptnudge !:
	110	REAL(wp) :: lowunudge !:
	111	REAL(wp) :: lowvnudge !:
	112	REAL(wp) :: lowwnudge !:
	113	REAL(wp) :: lowtnudge !:
	114
	115	REAL(wp) :: fac !:
	116
[1239]	117	ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), &
	118	tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge), &
	119	vnudge(nzb:nzt+1,1:ntnudge), wnudge(nzb:nzt+1,1:ntnudge) )
	120
	121	ALLOCATE( timenudge(0:ntnudge) )
	122
	123
[1353]	124	ptnudge = 0.0_wp; qnudge = 0.0_wp; tnudge = 0.0_wp; unudge = 0.0_wp
	125	vnudge = 0.0_wp; wnudge = 0.0_wp; timenudge = 0.0_wp
[1239]	126
	127	t = 0
[1249]	128	OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', &
	129	FORM='FORMATTED', IOSTAT=ierrn )
[1239]	130
[1249]	131	IF ( ierrn /= 0 ) THEN
[1239]	132	message_string = 'file NUDGING_DATA does not exist'
	133	CALL message( 'nudging', 'PA0365', 1, 2, 0, 6, 0 )
	134	ENDIF
	135
	136	ierrn = 0
	137
	138	rloop:DO
	139	t = t + 1
	140	hash = "#"
[1320]	141	ierrn = 1 ! not zero
[1239]	142	!
	143	!-- Search for the next line consisting of "# time",
	144	!-- from there onwards the profiles will be read
	145	DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) )
	146
[1249]	147	READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(t)
	148	IF ( ierrn < 0 ) EXIT rloop
[1239]	149
	150	ENDDO
	151
	152	ierrn = 0
[1249]	153	READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, &
	154	lowvnudge, lowwnudge , lowptnudge, &
	155	lowqnudge
[1239]	156
[1249]	157	IF ( ierrn /= 0 ) THEN
[1239]	158	message_string = 'errors in file NUDGING_DATA'
	159	CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
	160	ENDIF
	161
	162	ierrn = 0
[1249]	163	READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, &
	164	highvnudge, highwnudge , highptnudge, &
	165	highqnudge
[1239]	166
[1249]	167	IF ( ierrn /= 0 ) THEN
[1239]	168	message_string = 'errors in file NUDGING_DATA'
	169	CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
	170	ENDIF
	171
	172	DO k = nzb, nzt+1
[1249]	173	IF ( highheight < zu(k) ) THEN
[1239]	174	lowheight = highheight
	175	lowtnudge = hightnudge
	176	lowunudge = highunudge
	177	lowvnudge = highvnudge
	178	lowwnudge = highwnudge
	179	lowptnudge = highptnudge
	180	lowqnudge = highqnudge
	181
	182	ierrn = 0
[1249]	183	READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , &
	184	highunudge , highvnudge , &
	185	highwnudge , highptnudge, &
	186	highqnudge
	187	IF (ierrn /= 0 ) THEN
[1355]	188	WRITE( message_string, * ) 'zu(nzt+1) = ', zu(nzt+1), 'm is ',&
	189	'higher than the maximum height in NUDING_DATA which ', &
	190	'is ', lowheight, 'm. Interpolation on PALM ', &
	191	'grid is not possible.'
	192	CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 )
[1239]	193	ENDIF
	194	ENDIF
	195
	196	!
	197	!-- Interpolation of prescribed profiles in space
	198
[1249]	199	fac = ( highheight - zu(k) ) / ( highheight - lowheight )
[1239]	200
[1353]	201	tnudge(k,t) = fac * lowtnudge + ( 1.0_wp - fac ) * hightnudge
	202	unudge(k,t) = fac * lowunudge + ( 1.0_wp - fac ) * highunudge
	203	vnudge(k,t) = fac * lowvnudge + ( 1.0_wp - fac ) * highvnudge
	204	wnudge(k,t) = fac * lowwnudge + ( 1.0_wp - fac ) * highwnudge
	205	ptnudge(k,t) = fac * lowptnudge + ( 1.0_wp - fac ) * highptnudge
	206	qnudge(k,t) = fac * lowqnudge + ( 1.0_wp - fac ) * highqnudge
[1239]	207	ENDDO
	208
	209	ENDDO rloop
	210
[1249]	211	CLOSE ( finput )
[1239]	212
	213	!
	214	!-- Prevent nudging if nudging profiles exhibt too small values
[1241]	215	!-- not used so far
[1353]	216	lptnudge = ANY( ABS( ptnudge ) > 1.0e-8_wp )
	217	lqnudge = ANY( ABS( qnudge ) > 1.0e-8_wp )
	218	lunudge = ANY( ABS( unudge ) > 1.0e-8_wp )
	219	lvnudge = ANY( ABS( vnudge ) > 1.0e-8_wp )
	220	lwnudge = ANY( ABS( wnudge ) > 1.0e-8_wp )
[1239]	221
	222	END SUBROUTINE init_nudge
	223
	224	!------------------------------------------------------------------------------!
	225	! Call for all grid points
	226	!------------------------------------------------------------------------------!
	227	SUBROUTINE nudge ( time, prog_var )
	228
[1320]	229	USE arrays_3d, &
	230	ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, &
	231	v, vnudge
[1239]	232
[1320]	233	USE buoyancy_mod, &
	234	ONLY: calc_mean_profile
	235
	236	USE control_parameters, &
	237	ONLY: dt_3d, message_string
	238
	239	USE indices, &
	240	ONLY: nxl, nxr, nys, nyn, nzb, nzb_u_inner, nzt
	241
	242	USE kinds, &
	243	ONLY: iwp, wp
	244
	245	USE statistics, &
	246	ONLY: hom
	247
[1239]	248	IMPLICIT NONE
	249
[1320]	250	CHARACTER (LEN=*) :: prog_var !:
[1239]	251
[1320]	252	REAL(wp) :: currtnudge !:
	253	REAL(wp) :: dtm !:
	254	REAL(wp) :: dtp !:
	255	REAL(wp) :: time !:
[1239]	256
[1320]	257	INTEGER(iwp) :: i !:
	258	INTEGER(iwp) :: j !:
	259	INTEGER(iwp) :: k !:
	260	INTEGER(iwp) :: t !:
[1239]	261
	262
[1251]	263	t = 1
	264	DO WHILE ( time > timenudge(t) )
	265	t = t+1
	266	ENDDO
	267	IF ( time /= timenudge(1) ) THEN
	268	t = t-1
	269	ENDIF
	270
	271	dtm = ( time - timenudge(t) ) / ( timenudge(t+1) - timenudge(t) )
	272	dtp = ( timenudge(t+1) - time ) / ( timenudge(t+1) - timenudge(t) )
	273
[1239]	274	SELECT CASE ( prog_var )
	275
	276	CASE ( 'u' )
	277
	278	CALL calc_mean_profile( u, 1, 'nudging' )
	279
	280	DO i = nxl, nxr
	281	DO j = nys, nyn
	282	DO k = nzb_u_inner(j,i)+1, nzt
	283
	284	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + &
	285	tnudge(k,t+1) * dtm )
	286
	287
	288	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,1,0) - &
	289	( unudge(k,t) * dtp + &
	290	unudge(k,t+1) * dtm ) ) / currtnudge
	291	ENDDO
	292	ENDDO
	293	ENDDO
	294
	295	CASE ( 'v' )
	296
	297	CALL calc_mean_profile( v, 2, 'nudging' )
	298
	299	DO i = nxl, nxr
	300	DO j = nys, nyn
	301	DO k = nzb_u_inner(j,i)+1, nzt
	302
	303	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + &
	304	tnudge(k,t+1) * dtm )
	305
	306
	307	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,2,0) - &
	308	( vnudge(k,t) * dtp + &
	309	vnudge(k,t+1) * dtm ) ) / currtnudge
	310	ENDDO
	311	ENDDO
	312	ENDDO
	313
	314	CASE ( 'pt' )
	315
[1268]	316	CALL calc_mean_profile( pt, 4, 'nudging' )
[1239]	317
	318	DO i = nxl, nxr
	319	DO j = nys, nyn
	320	DO k = nzb_u_inner(j,i)+1, nzt
	321
	322	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + &
	323	tnudge(k,t+1) * dtm )
	324
	325
	326	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,4,0) - &
	327	( ptnudge(k,t) * dtp + &
	328	ptnudge(k,t+1) * dtm ) ) / currtnudge
	329	ENDDO
	330	ENDDO
	331	ENDDO
	332
	333	CASE ( 'q' )
	334
[1268]	335	CALL calc_mean_profile( q, 41, 'nudging' )
[1239]	336
	337	DO i = nxl, nxr
	338	DO j = nys, nyn
	339	DO k = nzb_u_inner(j,i)+1, nzt
	340
	341	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + &
	342	tnudge(k,t+1) * dtm )
	343
	344
	345	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,41,0) - &
	346	( qnudge(k,t) * dtp + &
	347	qnudge(k,t+1) * dtm ) ) / currtnudge
	348	ENDDO
	349	ENDDO
	350	ENDDO
	351
	352	CASE DEFAULT
	353	message_string = 'unknown prognostic variable "' // prog_var // '"'
	354	CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
	355
	356	END SELECT
	357
	358	END SUBROUTINE nudge
	359
	360
	361	!------------------------------------------------------------------------------!
	362	! Call for grid point i,j
	363	!------------------------------------------------------------------------------!
	364
	365	SUBROUTINE nudge_ij( i, j, time, prog_var )
	366
[1320]	367	USE arrays_3d, &
	368	ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, &
	369	v, vnudge
[1239]	370
[1320]	371	USE buoyancy_mod, &
	372	ONLY: calc_mean_profile
	373
	374	USE control_parameters, &
	375	ONLY: dt_3d, message_string
	376
	377	USE indices, &
	378	ONLY: nxl, nxr, nys, nyn, nzb, nzb_u_inner, nzt
	379
	380	USE kinds, &
	381	ONLY: iwp, wp
	382
	383	USE statistics, &
	384	ONLY: hom
	385
[1239]	386	IMPLICIT NONE
	387
	388
[1320]	389	CHARACTER (LEN=*) :: prog_var !:
[1239]	390
[1320]	391	REAL(wp) :: currtnudge !:
	392	REAL(wp) :: dtm !:
	393	REAL(wp) :: dtp !:
	394	REAL(wp) :: time !:
[1239]	395
[1320]	396	INTEGER(iwp) :: i !:
	397	INTEGER(iwp) :: j !:
	398	INTEGER(iwp) :: k !:
	399	INTEGER(iwp) :: t !:
[1239]	400
[1320]	401
[1239]	402	t = 1
	403	DO WHILE ( time > timenudge(t) )
	404	t = t+1
	405	ENDDO
[1249]	406	IF ( time /= timenudge(1) ) THEN
[1239]	407	t = t-1
	408	ENDIF
	409
	410	dtm = ( time - timenudge(t) ) / ( timenudge(t+1) - timenudge(t) )
	411	dtp = ( timenudge(t+1) - time ) / ( timenudge(t+1) - timenudge(t) )
	412
	413	SELECT CASE ( prog_var )
	414
	415	CASE ( 'u' )
	416
	417	CALL calc_mean_profile( u, 1, 'nudging' )
	418
	419	DO k = nzb_u_inner(j,i)+1, nzt
	420
	421	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm )
	422
	423	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,1,0) - &
	424	( unudge(k,t) * dtp + &
	425	unudge(k,t+1) * dtm ) ) / currtnudge
	426	ENDDO
	427
	428	CASE ( 'v' )
	429
	430	CALL calc_mean_profile( v, 2, 'nudging' )
	431
	432	DO k = nzb_u_inner(j,i)+1, nzt
	433
	434	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm )
	435
	436	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,2,0) - &
	437	( vnudge(k,t) * dtp + &
	438	vnudge(k,t+1) * dtm ) ) / currtnudge
	439	ENDDO
	440
	441
	442	CASE ( 'pt' )
	443
	444	CALL calc_mean_profile( pt, 4, 'nudging' )
	445
	446	DO k = nzb_u_inner(j,i)+1, nzt
	447
	448	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm )
	449
	450	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,4,0) - &
	451	( ptnudge(k,t) * dtp + &
	452	ptnudge(k,t+1) * dtm ) ) / currtnudge
	453	ENDDO
	454
	455
	456	CASE ( 'q' )
	457
	458	CALL calc_mean_profile( q, 41, 'nudging' )
	459
	460	DO k = nzb_u_inner(j,i)+1, nzt
	461
	462	currtnudge = MAX( dt_3d, tnudge(k,t) * dtp + tnudge(k,t+1) * dtm )
	463
	464	tend(k,j,i) = tend(k,j,i) - ( hom(k,1,41,0) - &
	465	( qnudge(k,t) * dtp + &
	466	qnudge(k,t+1) * dtm ) ) / currtnudge
	467	ENDDO
	468
	469	CASE DEFAULT
	470	message_string = 'unknown prognostic variable "' // prog_var // '"'
	471	CALL message( 'nudge', 'PA0367', 1, 2, 0, 6, 0 )
	472
	473	END SELECT
	474
	475
	476	END SUBROUTINE nudge_ij
	477
	478	END MODULE nudge_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |