Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |