Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |