Home

Context Navigation

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

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

Line
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-2014 Leibniz Universitaet Hannover
18	!--------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! ------------------
22	! Error message specified.
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: nudging.f90 1355 2014-04-10 10:21:29Z heinze $
27	!
28	! 1353 2014-04-08 15:21:23Z heinze
29	! REAL constants provided with KIND-attribute
30	!
31	! 1320 2014-03-20 08:40:49Z raasch
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
39	!
40	! 1318 2014-03-17 13:35:16Z raasch
41	! module interfaces removed
42	!
43	! 1268 2013-12-12 09:47:53Z heinze
44	! bugfix: argument of calc_mean_profile corrected
45	!
46	! 1251 2013-11-07 08:14:30Z heinze
47	! bugfix: calculate dtm and dtp also in vector version
48	!
49	! 1249 2013-11-06 10:45:47Z heinze
50	! remove call of user module
51	! reformatting
52	!
53	! 1241 2013-10-30 11:36:58Z heinze
54	! Initial revision
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
76	USE arrays_3d, &
77	ONLY: ptnudge, qnudge, timenudge, tnudge, unudge, vnudge, wnudge, &
78	zu
79
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
89	IMPLICIT NONE
90
91
92	INTEGER(iwp) :: finput = 90 !:
93	INTEGER(iwp) :: ierrn !:
94	INTEGER(iwp) :: k !:
95	INTEGER(iwp) :: t !:
96
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
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
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
126
127	t = 0
128	OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', &
129	FORM='FORMATTED', IOSTAT=ierrn )
130
131	IF ( ierrn /= 0 ) THEN
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 = "#"
141	ierrn = 1 ! not zero
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
147	READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(t)
148	IF ( ierrn < 0 ) EXIT rloop
149
150	ENDDO
151
152	ierrn = 0
153	READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, &
154	lowvnudge, lowwnudge , lowptnudge, &
155	lowqnudge
156
157	IF ( ierrn /= 0 ) THEN
158	message_string = 'errors in file NUDGING_DATA'
159	CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
160	ENDIF
161
162	ierrn = 0
163	READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, &
164	highvnudge, highwnudge , highptnudge, &
165	highqnudge
166
167	IF ( ierrn /= 0 ) THEN
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
173	IF ( highheight < zu(k) ) THEN
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
183	READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , &
184	highunudge , highvnudge , &
185	highwnudge , highptnudge, &
186	highqnudge
187	IF (ierrn /= 0 ) THEN
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 )
193	ENDIF
194	ENDIF
195
196	!
197	!-- Interpolation of prescribed profiles in space
198
199	fac = ( highheight - zu(k) ) / ( highheight - lowheight )
200
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
207	ENDDO
208
209	ENDDO rloop
210
211	CLOSE ( finput )
212
213	!
214	!-- Prevent nudging if nudging profiles exhibt too small values
215	!-- not used so far
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 )
221
222	END SUBROUTINE init_nudge
223
224	!------------------------------------------------------------------------------!
225	! Call for all grid points
226	!------------------------------------------------------------------------------!
227	SUBROUTINE nudge ( time, prog_var )
228
229	USE arrays_3d, &
230	ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, &
231	v, vnudge
232
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
248	IMPLICIT NONE
249
250	CHARACTER (LEN=*) :: prog_var !:
251
252	REAL(wp) :: currtnudge !:
253	REAL(wp) :: dtm !:
254	REAL(wp) :: dtp !:
255	REAL(wp) :: time !:
256
257	INTEGER(iwp) :: i !:
258	INTEGER(iwp) :: j !:
259	INTEGER(iwp) :: k !:
260	INTEGER(iwp) :: t !:
261
262
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
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
316	CALL calc_mean_profile( pt, 4, 'nudging' )
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
335	CALL calc_mean_profile( q, 41, 'nudging' )
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
367	USE arrays_3d, &
368	ONLY: pt, ptnudge, q, qnudge, tend, timenudge, tnudge, u, unudge, &
369	v, vnudge
370
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
386	IMPLICIT NONE
387
388
389	CHARACTER (LEN=*) :: prog_var !:
390
391	REAL(wp) :: currtnudge !:
392	REAL(wp) :: dtm !:
393	REAL(wp) :: dtp !:
394	REAL(wp) :: time !:
395
396	INTEGER(iwp) :: i !:
397	INTEGER(iwp) :: j !:
398	INTEGER(iwp) :: k !:
399	INTEGER(iwp) :: t !:
400
401
402	t = 1
403	DO WHILE ( time > timenudge(t) )
404	t = t+1
405	ENDDO
406	IF ( time /= timenudge(1) ) THEN
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 |