Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |