Home

Context Navigation

source: palm/trunk/UTIL/inifor/src/inifor.f90 @ 3448

Last change on this file since 3448 was 3401, checked in by eckhard, 6 years ago
inifor: Re-compute geostrophic winds every time step
Property svn:keywords set to `Id`
File size: 20.8 KB

Line
1	!> @file src/inifor.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! 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 2017-2018 Leibniz Universitaet Hannover
18	! Copyright 2017-2018 Deutscher Wetterdienst Offenbach
19	!------------------------------------------------------------------------------!
20	!
21	! Current revisions:
22	! -----------------
23	!
24	!
25	! Former revisions:
26	! -----------------
27	! $Id: inifor.f90 3401 2018-10-23 09:33:31Z kanani $
28	! Re-compute geostrophic winds every time step
29	!
30	! 3395 2018-10-22 17:32:49Z eckhard
31	! Added main loop support for computation of geostrophic winds and surface
32	! pressure
33	! Cleaned up terminal output, show some meVssages only with --debug option
34	!
35	! 3183 2018-07-27 14:25:55Z suehring
36	! Introduced new PALM grid stretching
37	! Renamend initial-condition mode variable 'mode' to 'ic_mode'
38	! Improved log messages
39	!
40	!
41	! 3182 2018-07-27 13:36:03Z suehring
42	! Initial revision
43	!
44	!
45	!
46	! Authors:
47	! --------
48	! @author Eckhard Kadasch
49	!
50	! Description:
51	! ------------
52	!> INIFOR is an interpolation tool for generating meteorological initialization
53	!> and forcing data for the urban climate model PALM-4U. The required
54	!> meteorological fields are interpolated from output data of the mesoscale
55	!> model COSMO-DE. This is the main program file.
56	!------------------------------------------------------------------------------!
57	PROGRAM inifor
58
59	USE control
60	USE defs
61	USE grid, &
62	ONLY: setup_parameters, setup_grids, setup_variable_tables, &
63	setup_io_groups, fini_grids, fini_variables, fini_io_groups, &
64	fini_file_lists, preprocess, origin_lon, origin_lat, &
65	output_file, io_group_list, output_var_table, &
66	cosmo_grid, palm_grid, nx, ny, nz, p0, cfg, f3, &
67	averaging_width_ns, averaging_width_ew, phi_n, lambda_n, &
68	lam_centre, phi_centre
69	USE io
70	USE transform, &
71	ONLY: average_profile, interpolate_2d, interpolate_3d, &
72	geostrophic_winds, extrapolate_density, extrapolate_pressure, &
73	get_surface_pressure
74	USE types
75
76	IMPLICIT NONE
77
78	INTEGER :: igroup
79	INTEGER :: ivar
80	INTEGER :: iter
81	REAL(dp), ALLOCATABLE, DIMENSION(:,:,:) :: output_arr
82	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_centre
83	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: ug_arr
84	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: vg_arr
85	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_north
86	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_south
87	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_east
88	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_west
89	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_north
90	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_south
91	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_east
92	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_west
93	REAL(dp), POINTER, DIMENSION(:) :: internal_arr
94	REAL(dp), POINTER, DIMENSION(:) :: ug_vg_arr
95	TYPE(nc_var), POINTER :: output_var
96	TYPE(io_group), POINTER :: group
97	TYPE(container), ALLOCATABLE :: input_buffer(:)
98	LOGICAL, SAVE :: ug_vg_have_been_computed = .FALSE.
99	LOGICAL, SAVE :: debugging_output = .TRUE.
100
101	!> \mainpage About INIFOR
102	!> ...
103	!
104	!------------------------------------------------------------------------------
105	!- Section 1: Initialization
106	!------------------------------------------------------------------------------
107	CALL run_control('init', 'void')
108
109	! Initialize INIFOR's parameters from command-line interface and namelists
110	CALL setup_parameters()
111
112	! Initialize all grids, including interpolation neighbours and weights
113	CALL setup_grids()
114	CALL run_control('time', 'init')
115
116	! Initialize the netCDF output file and define dimensions
117	CALL setup_netcdf_dimensions(output_file, palm_grid, cfg % start_date, &
118	origin_lon, origin_lat)
119	CALL run_control('time', 'write')
120
121	! Set up the tables containing the input and output variables and set
122	! the corresponding netCDF dimensions for each output variable
123	CALL setup_variable_tables(cfg % ic_mode)
124	CALL run_control('time', 'write')
125
126	! Add the output variables to the netCDF output file
127	CALL setup_netcdf_variables(output_file % name, output_var_table)
128
129	CALL setup_io_groups()
130	CALL run_control('time', 'init')
131
132	!------------------------------------------------------------------------------
133	!- Section 2: Main loop
134	!------------------------------------------------------------------------------
135
136	DO igroup = 1, SIZE(io_group_list)
137
138	group => io_group_list(igroup)
139	IF ( group % to_be_processed ) THEN
140
141	DO iter = 1, group % nt
142
143	!------------------------------------------------------------------------------
144	!- Section 2.1: Read and preprocess input data
145	!------------------------------------------------------------------------------
146	CALL read_input_variables(group, iter, input_buffer)
147	CALL run_control('time', 'read')
148
149	CALL preprocess(group, input_buffer, cosmo_grid, iter)
150	CALL run_control('time', 'comp')
151
152	!TODO: move this assertion into 'preprocess'.
153	IF ( .NOT. ALL(input_buffer(:) % is_preprocessed .AND. .TRUE.) ) THEN
154	message = "Input buffers for group '" // TRIM(group % kind) // &
155	"' could not be preprocessed sucessfully."
156	CALL abort('main loop', message)
157	END IF
158
159	!------------------------------------------------------------------------------
160	!- Section 2.2: Interpolate each output variable of the group
161	!------------------------------------------------------------------------------
162	DO ivar = 1, group % nv
163
164	output_var => group % out_vars( ivar )
165
166	IF ( output_var % to_be_processed .AND. &
167	iter .LE. output_var % nt ) THEN
168
169	message = "Processing '" // TRIM(output_var % name) // &
170	"' (" // TRIM(output_var % kind) // &
171	"), iteration " // TRIM(str(iter)) //" of " // &
172	TRIM(str(output_var % nt))
173	CALL report('main loop', message)
174
175	SELECT CASE( TRIM(output_var % task) )
176
177	CASE( 'interpolate_2d' )
178
179	SELECT CASE( TRIM(output_var % kind) )
180
181	CASE( 'init soil' )
182
183	ALLOCATE( output_arr( 0:output_var % grid % nx, &
184	0:output_var % grid % ny, &
185	SIZE(output_var % grid % depths) ) )
186
187	CASE ( 'surface forcing' )
188
189	ALLOCATE( output_arr( 0:output_var % grid % nx, &
190	0:output_var % grid % ny, 1 ) )
191
192	CASE DEFAULT
193
194	message = "'" // TRIM(output_var % kind) // "' is not a soil variable"
195	CALL abort("main loop", message)
196
197	END SELECT
198	CALL run_control('time', 'alloc')
199
200	CALL interpolate_2d(input_buffer(output_var % input_id) % array(:,:,:), &
201	output_arr(:,:,:), output_var % intermediate_grid, output_var)
202	CALL run_control('time', 'comp')
203
204
205	CASE ( 'interpolate_3d' )
206
207	ALLOCATE( output_arr( 0:output_var % grid % nx, &
208	0:output_var % grid % ny, &
209	1:output_var % grid % nz ) )
210
211	CALL run_control('time', 'alloc')
212	CALL interpolate_3d( &
213	input_buffer(output_var % input_id) % array(:,:,:), &
214	output_arr(:,:,:), &
215	output_var % intermediate_grid, &
216	output_var % grid)
217	CALL run_control('time', 'comp')
218
219	CASE ( 'average profile' )
220
221	ALLOCATE( output_arr( 0:output_var % grid % nx, &
222	0:output_var % grid % ny, &
223	1:output_var % grid % nz ) )
224	CALL run_control('time', 'alloc')
225
226	CALL average_profile( &
227	input_buffer(output_var % input_id) % array(:,:,:), &
228	output_arr(0,0,:), &
229	output_var % averaging_grid)
230
231	IF ( TRIM(output_var % name) == &
232	'surface_forcing_surface_pressure' ) THEN
233
234	IF ( cfg % p0_is_set ) THEN
235	output_arr(0,0,1) = p0
236	ELSE
237	CALL get_surface_pressure( &
238	output_arr(0,0,:), rho_centre, &
239	output_var % averaging_grid)
240	END IF
241
242	END IF
243	CALL run_control('time', 'comp')
244
245	CASE ( 'internal profile' )
246
247	message = "Averaging of internal profile for variable '" //&
248	TRIM(output_var % name) // "' is not supported."
249
250	SELECT CASE (TRIM(output_var % name))
251
252	CASE('internal_density_centre')
253	ALLOCATE( rho_centre( 1:output_var % grid % nz ) )
254	internal_arr => rho_centre
255
256	CASE('internal_density_north')
257	ALLOCATE( rho_north( 1:output_var % grid % nz ) )
258	internal_arr => rho_north
259
260	CASE('internal_density_south')
261	ALLOCATE( rho_south( 1:output_var % grid % nz ) )
262	internal_arr => rho_south
263
264	CASE('internal_density_east')
265	ALLOCATE( rho_east( 1:output_var % grid % nz) )
266	internal_arr => rho_east
267
268	CASE('internal_density_west')
269	ALLOCATE( rho_west( 1:output_var % grid % nz ) )
270	internal_arr => rho_west
271
272	CASE('internal_pressure_north')
273	ALLOCATE( p_north( 1:output_var % grid % nz ) )
274	internal_arr => p_north
275
276	CASE('internal_pressure_south')
277	ALLOCATE( p_south( 1:output_var % grid % nz ) )
278	internal_arr => p_south
279
280	CASE('internal_pressure_east')
281	ALLOCATE( p_east( 1:output_var % grid % nz) )
282	internal_arr => p_east
283
284	CASE('internal_pressure_west')
285	ALLOCATE( p_west( 1:output_var % grid % nz ) )
286	internal_arr => p_west
287
288	CASE DEFAULT
289	CALL abort('main loop', message)
290
291	END SELECT
292	CALL run_control('time', 'alloc')
293
294
295	CALL average_profile( &
296	input_buffer(output_var % input_id) % array(:,:,:),&
297	internal_arr(:), &
298	output_var % averaging_grid)
299
300	SELECT CASE (TRIM(output_var % name))
301
302	CASE('internal_density_centre', &
303	'internal_density_north', &
304	'internal_density_south', &
305	'internal_density_east', &
306	'internal_density_west')
307	CALL extrapolate_density(internal_arr, &
308	output_var % averaging_grid)
309
310	CASE('internal_pressure_north')
311	CALL extrapolate_pressure(internal_arr, rho_north, &
312	output_var % averaging_grid)
313
314	CASE('internal_pressure_south')
315	CALL extrapolate_pressure(internal_arr, rho_south, &
316	output_var % averaging_grid)
317
318	CASE('internal_pressure_east')
319	CALL extrapolate_pressure(internal_arr, rho_east, &
320	output_var % averaging_grid)
321
322	CASE('internal_pressure_west')
323	CALL extrapolate_pressure(internal_arr, rho_west, &
324	output_var % averaging_grid)
325
326	CASE DEFAULT
327	CALL abort('main loop', message)
328
329	END SELECT
330
331	IF (.TRUE.) THEN
332	ALLOCATE( output_arr(1,1,1:output_var % grid % nz) )
333	output_arr(1,1,:) = internal_arr(:)
334	END IF
335	CALL run_control('time', 'comp')
336
337
338	! This case gets called twice, the first time for ug, the
339	! second time for vg. We compute ug and vg at the first call
340	! and keep vg (and ug for that matter) around for the second
341	! call.
342	CASE ( 'geostrophic winds' )
343
344	IF (.NOT. ug_vg_have_been_computed ) THEN
345	ALLOCATE( ug_arr(output_var % grid % nz) )
346	ALLOCATE( vg_arr(output_var % grid % nz) )
347
348	IF ( cfg % ug_is_set ) THEN
349	ug_arr = cfg % ug
350	vg_arr = cfg % vg
351	ELSE
352	CALL geostrophic_winds( p_north, p_south, p_east, &
353	p_west, rho_centre, f3, &
354	averaging_width_ew, &
355	averaging_width_ns, &
356	phi_n, lambda_n, phi_centre, lam_centre, &
357	ug_arr, vg_arr )
358	END IF
359
360	ug_vg_have_been_computed = .TRUE.
361
362	END IF
363
364	! Prepare output of geostrophic winds
365	SELECT CASE(TRIM(output_var % name))
366	CASE ('ls_forcing_ug')
367	ug_vg_arr => ug_arr
368	CASE ('ls_forcing_vg')
369	ug_vg_arr => vg_arr
370	END SELECT
371
372	ALLOCATE( output_arr(1,1,output_var % grid % nz) )
373	output_arr(1,1,:) = ug_vg_arr(:)
374
375	CASE ( 'average scalar' )
376
377	ALLOCATE( output_arr(1,1,1) )
378	CALL run_control('time', 'alloc')
379	output_arr(1,1,1) = p0
380	CALL run_control('time', 'comp')
381
382	CASE ( 'set profile' )
383
384	ALLOCATE( output_arr( 1, 1, 1:nz ) )
385	CALL run_control('time', 'alloc')
386
387	SELECT CASE (TRIM(output_var % name))
388
389	!CASE('ls_forcing_ug')
390	! output_arr(1, 1, :) = ug
391
392	!CASE('ls_forcing_vg')
393	! output_arr(1, 1, :) = vg
394
395	CASE('nudging_tau')
396	output_arr(1, 1, :) = NUDGING_TAU
397
398	CASE DEFAULT
399	message = "'" // TRIM(output_var % name) // &
400	"' is not a valid '" // TRIM(output_var % kind) //&
401	"' variable kind."
402	CALL abort('main loop', message)
403	END SELECT
404	CALL run_control('time', 'comp')
405
406	CASE('average large-scale profile')
407	message = "Averaging of large-scale forcing profiles " //&
408	"has not been implemented, yet."
409	CALL abort('main loop', message)
410	!ALLOCATE( output_arr( 1, 1, 1:nz ) )
411
412	CASE DEFAULT
413	message = "Processing task '" // TRIM(output_var % task) //&
414	"' not recognized."
415	CALL abort('', message)
416
417	END SELECT
418	CALL run_control('time', 'comp')
419
420	!------------------------------------------------------------------------------
421	!- Section 2.3: Write current time step of current variable
422	!------------------------------------------------------------------------------
423	IF (.NOT. output_var % is_internal .OR. debugging_output) THEN
424	message = "Writing variable '" // TRIM(output_var%name) // "'."
425	CALL report('main loop', message)
426	CALL update_output(output_var, output_arr, iter, output_file)
427	CALL run_control('time', 'write')
428	END IF
429
430	IF (ALLOCATED(output_arr)) DEALLOCATE(output_arr)
431	CALL run_control('time', 'alloc')
432
433	END IF
434
435	END DO ! ouput variables
436
437	ug_vg_have_been_computed = .FALSE.
438	IF ( group % kind == 'thermodynamics' ) THEN
439	DEALLOCATE( rho_centre )
440	DEALLOCATE( ug_arr )
441	DEALLOCATE( vg_arr )
442	IF ( .NOT. cfg % ug_is_set ) THEN
443	DEALLOCATE( rho_north )
444	DEALLOCATE( rho_south )
445	DEALLOCATE( rho_east )
446	DEALLOCATE( rho_west )
447	DEALLOCATE( p_north )
448	DEALLOCATE( p_south )
449	DEALLOCATE( p_east )
450	DEALLOCATE( p_west )
451	END IF
452	END IF
453
454	IF ( group % kind == 'running average' .OR. &
455	group % kind == 'accumulated' ) THEN
456	! Keep input buffer around for averaged (radiation) and
457	! accumulated COSMO-DE quantities (precipitation).
458	ELSE
459	CALL report('main loop', 'Deallocating input buffer', cfg % debug)
460	DEALLOCATE(input_buffer)
461	END IF
462	CALL run_control('time', 'alloc')
463
464	END DO ! time steps / input files
465
466	IF (ALLOCATED(input_buffer)) THEN
467	CALL report('main loop', 'Deallocating input buffer', cfg % debug)
468	DEALLOCATE(input_buffer)
469	END IF
470	CALL run_control('time', 'alloc')
471
472	ELSE
473
474	message = "Skipping IO group " // TRIM(str(igroup)) // " '" // TRIM(group % kind) // "'"
475	IF ( ALLOCATED(group % in_var_list) ) THEN
476	message = TRIM(message) // " with input variable '" // &
477	TRIM(group % in_var_list(1) % name) // "'."
478	END IF
479
480	CALL report('main loop', message, cfg % debug)
481
482	END IF ! IO group % to_be_processed
483
484	END DO ! IO groups
485
486	!------------------------------------------------------------------------------
487	!- Section 3: Clean up.
488	!------------------------------------------------------------------------------
489	CALL fini_file_lists()
490	CALL fini_io_groups()
491	CALL fini_variables()
492	!CALL fini_grids()
493	CALL run_control('time', 'alloc')
494	CALL run_control('report', 'void')
495
496	message = "Finished writing dynamic driver '" // TRIM(output_file % name) // &
497	"' successfully."
498	CALL report('main loop', message)
499
500
501	END PROGRAM inifor

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |