Home

Context Navigation

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

Last change on this file since 1313 was 1310, checked in by raasch, 11 years ago
update of GPL copyright
Property svn:keywords set to `Id`
File size: 12.1 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	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: nudging.f90 1310 2014-03-14 08:01:56Z fricke $
27	!
28	! 1268 2013-12-12 09:47:53Z heinze
29	! bugfix: argument of calc_mean_profile corrected
30	!
31	! 1251 2013-11-07 08:14:30Z heinze
32	! bugfix: calculate dtm and dtp also in vector version
33	!
34	! 1249 2013-11-06 10:45:47Z heinze
35	! remove call of user module
36	! reformatting
37	!
38	! 1241 2013-10-30 11:36:58Z heinze
39	! Initial revision
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
90	OPEN ( finput, FILE='NUDGING_DATA', STATUS='OLD', &
91	FORM='FORMATTED', IOSTAT=ierrn )
92
93	IF ( ierrn /= 0 ) THEN
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
109	READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(t)
110	IF ( ierrn < 0 ) EXIT rloop
111
112	ENDDO
113
114	ierrn = 0
115	READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, &
116	lowvnudge, lowwnudge , lowptnudge, &
117	lowqnudge
118
119	IF ( ierrn /= 0 ) THEN
120	message_string = 'errors in file NUDGING_DATA'
121	CALL message( 'nudging', 'PA0366', 1, 2, 0, 6, 0 )
122	ENDIF
123
124	ierrn = 0
125	READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, &
126	highvnudge, highwnudge , highptnudge, &
127	highqnudge
128
129	IF ( ierrn /= 0 ) THEN
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
135	IF ( highheight < zu(k) ) THEN
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
145	READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , &
146	highunudge , highvnudge , &
147	highwnudge , highptnudge, &
148	highqnudge
149	IF (ierrn /= 0 ) THEN
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
158	fac = ( highheight - zu(k) ) / ( highheight - lowheight )
159
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
166	ENDDO
167
168	ENDDO rloop
169
170	CLOSE ( finput )
171
172	!
173	!-- Prevent nudging if nudging profiles exhibt too small values
174	!-- not used so far
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
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
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
259	CALL calc_mean_profile( pt, 4, 'nudging' )
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
278	CALL calc_mean_profile( q, 41, '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,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
332	IF ( time /= timenudge(1) ) THEN
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 |