Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |