Home

Context Navigation

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

Last change on this file since 4010 was 3997, checked in by eckhard, 6 years ago
inifor: Read origin_z from static driver if given; alert user to warnings
Property svn:keywords set to `Id`
File size: 27.5 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-2019 Leibniz Universitaet Hannover
18	! Copyright 2017-2019 Deutscher Wetterdienst Offenbach
19	!------------------------------------------------------------------------------!
20	!
21	! Current revisions:
22	! -----------------
23	!
24	!
25	! Former revisions:
26	! -----------------
27	! $Id: inifor.f90 3997 2019-05-23 12:35:57Z hellstea $
28	! Alert user to warnings if any
29	!
30	!
31	! 3866 2019-04-05 14:25:01Z eckhard
32	! Use PALM's working precision
33	! Show error message if compiled without netCDF support
34	! Renamed run_control -> log_runtime
35	! Improved coding style added comments
36	!
37	!
38	! 3785 2019-03-06 10:41:14Z eckhard
39	! Average geostrophic wind components on coarse COSMO levels instead of fine PALM levels
40	! Remove --debug netCDF output of internal pressure profiles
41	!
42	! 3680 2019-01-18 14:54:12Z knoop
43	! Prefixed all INIFOR modules with inifor_
44	!
45	!
46	! 3615 2018-12-10 07:21:03Z raasch
47	! bugfix: abort replaced by inifor_abort
48	!
49	! 3613 2018-12-07 18:20:37Z eckhard
50	! Moved version output to setup_parameters()
51	!
52	! 3557 2018-11-22 16:01:22Z eckhard
53	! Updated documentation
54	!
55	! 3537 2018-11-20 10:53:14Z eckhard
56	! Print version number on program start
57	!
58	!
59	! 3456 2018-10-30 14:29:54Z eckhard
60	! NetCDf output of internal arrays only with --debug option
61	!
62	!
63	! 3401 2018-10-23 09:33:31Z eckhard
64	! Re-compute geostrophic winds every time step
65	!
66	! 3395 2018-10-22 17:32:49Z eckhard
67	! Added main loop support for computation of geostrophic winds and surface
68	! pressure
69	! Cleaned up terminal output, show some meVssages only with --debug option
70	!
71	! 3183 2018-07-27 14:25:55Z suehring
72	! Introduced new PALM grid stretching
73	! Renamend initial-condition mode variable 'mode' to 'ic_mode'
74	! Improved log messages
75	!
76	!
77	! 3182 2018-07-27 13:36:03Z suehring
78	! Initial revision
79	!
80	!
81	!
82	! Authors:
83	! --------
84	!> @author Eckhard Kadasch, (Deutscher Wetterdienst, Offenbach)
85	!
86	! Description:
87	! ------------
88	!> INIFOR is an interpolation tool for generating meteorological initialization
89	!> and forcing data for the urban climate model PALM-4U. The required
90	!> meteorological fields are interpolated from output data of the mesoscale
91	!> model COSMO. This is the main program file.
92	!------------------------------------------------------------------------------!
93	PROGRAM inifor
94
95	#if defined ( __netcdf )
96
97	USE inifor_control
98	USE inifor_defs
99	USE inifor_grid, &
100	ONLY: averaging_width_ns, &
101	averaging_width_ew, &
102	cfg, &
103	cosmo_grid, &
104	f3, &
105	fini_grids, &
106	fini_io_groups, &
107	fini_variables, &
108	fini_file_lists, &
109	io_group_list, &
110	lam_centre, &
111	lambda_n, &
112	nx, ny, nz, &
113	origin_lat, &
114	origin_lon, &
115	output_file, &
116	output_var_table, &
117	p0, &
118	phi_centre, &
119	phi_n, &
120	preprocess, &
121	palm_grid, &
122	setup_grids, &
123	setup_parameters, &
124	setup_variable_tables, &
125	setup_io_groups
126	USE inifor_io
127	USE inifor_transform, &
128	ONLY: average_pressure_perturbation, &
129	average_profile, &
130	extrapolate_density, &
131	extrapolate_pressure, &
132	geostrophic_winds, &
133	get_surface_pressure, &
134	interp_average_profile, &
135	interpolate_1d, &
136	interpolate_1d_arr, &
137	interpolate_2d, &
138	interpolate_3d
139	USE inifor_types
140
141	IMPLICIT NONE
142
143	INTEGER :: igroup !< loop index for IO groups
144	INTEGER :: ivar !< loop index for output variables
145	INTEGER :: iter !< loop index for time steps
146
147	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: output_arr !< array buffer for interpolated quantities
148	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_centre !< density profile of the centre averaging domain
149	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: ug_cosmo !< profile of the geostrophic wind in x direction on COSMO levels
150	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: vg_cosmo !< profile of the geostrophic wind in y direction on COSMO levels
151	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: ug_palm !< profile of the geostrophic wind in x direction interpolated onto PALM levels
152	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: vg_palm !< profile of the geostrophic wind in y direction interpolated onto PALM levels
153	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_north !< density profile of the northern averaging domain
154	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_south !< density profile of the southern averaging domain
155	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_east !< density profile of the eastern averaging domain
156	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_west !< density profile of the western averaging domain
157	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: p_north !< pressure profile of the northern averaging domain
158	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: p_south !< pressure profile of the southern averaging domain
159	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: p_east !< pressure profile of the eastern averaging domain
160	REAL(wp), ALLOCATABLE, DIMENSION(:), TARGET :: p_west !< pressure profile of the western averaging domain
161
162	REAL(wp), POINTER, DIMENSION(:) :: internal_arr !< pointer to the currently processed internal array (density, pressure)
163	REAL(wp), POINTER, DIMENSION(:) :: ug_vg_palm !< pointer to the currently processed geostrophic wind component
164
165	TYPE(nc_var), POINTER :: output_var !< pointer to the currently processed output variable
166	TYPE(io_group), POINTER :: group !< pointer to the currently processed IO group
167	TYPE(container), ALLOCATABLE :: input_buffer(:) !< buffer of the current IO group's input arrays
168
169	LOGICAL, SAVE :: ug_vg_have_been_computed = .FALSE. !< flag for managing geostrophic wind allocation and computation
170	!LOGICAL, SAVE :: debugging_output = .TRUE. !< flag controllging output of internal variables
171
172	!> \mainpage About INIFOR
173	!> ...
174	!
175	!------------------------------------------------------------------------------
176	!- Section 1: Initialization
177	!------------------------------------------------------------------------------
178	CALL log_runtime( 'init', 'void' )
179
180	!
181	!-- Initialize INIFOR's parameters from command-line interface and namelists
182	CALL setup_parameters
183
184	!
185	!-- Initialize all grids, including interpolation neighbours and weights
186	CALL setup_grids
187	CALL log_runtime( 'time', 'init' )
188
189	!
190	!-- Initialize the netCDF output file and define dimensions
191	CALL setup_netcdf_dimensions( output_file, palm_grid, cfg%start_date, &
192	origin_lon, origin_lat )
193	CALL log_runtime( 'time', 'write' )
194
195	!
196	!-- Set up the tables containing the input and output variables and set
197	!-- the corresponding netCDF dimensions for each output variable
198	CALL setup_variable_tables( cfg%ic_mode )
199	CALL log_runtime( 'time', 'write' )
200
201	!
202	!-- Add the output variables to the netCDF output file
203	CALL setup_netcdf_variables( output_file%name, output_var_table )
204
205	CALL setup_io_groups
206	CALL log_runtime( 'time', 'init' )
207
208	!------------------------------------------------------------------------------
209	!-- Section 2: Main loop
210	!------------------------------------------------------------------------------
211	!
212	!-- Input and output variables are organized into IO groups. For instance, the
213	!-- 'thermodynamics' IO group bundles the input variaebls T, P, QV and the
214	!-- output variables p, theta, rho, and qv.
215	!-- An IO group bunldes variables that are physically dependent on each other.
216	!-- In case of the 'thermodynamics' group, theta = f(P,T), rho = f(P,T,QV).
217	DO igroup = 1, SIZE( io_group_list )
218
219	group => io_group_list(igroup)
220	IF ( group%to_be_processed ) THEN
221
222	!-- Loop over all output time steps for the current group.
223	DO iter = 1, group%nt
224
225	!------------------------------------------------------------------------------
226	!-- Section 2.1: Read and preprocess input data
227	!------------------------------------------------------------------------------
228	CALL read_input_variables( group, iter, input_buffer )
229	CALL log_runtime( 'time', 'read' )
230
231	!-- Carry out all required physical conversion of the input variables
232	!-- of the current IO group on the input (COSMO) grid. For instance,
233	!-- horizontal velocities are rotated to the PALM coordinate system and
234	!-- potential temperature is computed from the absolute temperature and
235	!-- pressure.
236	CALL preprocess( group, input_buffer, cosmo_grid, iter )
237	CALL log_runtime( 'time', 'comp' )
238
239	!TODO: move this assertion into 'preprocess'.
240	IF ( .NOT. ALL(input_buffer(:)%is_preprocessed .AND. .TRUE.) ) THEN
241	message = "Input buffers for group '" // TRIM( group%kind ) // &
242	"' could not be preprocessed sucessfully."
243	CALL inifor_abort( 'main loop', message )
244	ENDIF
245
246	!------------------------------------------------------------------------------
247	!-- Section 2.2: Interpolate each output variable of the group
248	!------------------------------------------------------------------------------
249	DO ivar = 1, group%nv
250
251	output_var => group%out_vars(ivar)
252
253	IF ( output_var%to_be_processed .AND. &
254	iter .LE. output_var%nt ) THEN
255
256	message = "Processing '" // TRIM( output_var%name ) // &
257	"' (" // TRIM( output_var%kind ) // &
258	"), iteration " // TRIM( str( iter ) ) //" of " //&
259	TRIM( str( output_var%nt ) )
260	CALL report( 'main loop', message )
261
262	SELECT CASE( TRIM( output_var%task ) )
263
264	!-- 2D horizontal interpolation
265	CASE( 'interpolate_2d' )
266
267	SELECT CASE( TRIM( output_var%kind ) )
268
269	CASE( 'init soil' )
270
271	ALLOCATE( output_arr(0:output_var%grid%nx, &
272	0:output_var%grid%ny, &
273	SIZE( output_var%grid%depths )) )
274
275	CASE ( 'surface forcing' )
276
277	ALLOCATE( output_arr(0:output_var%grid%nx, &
278	0:output_var%grid%ny, 1) )
279
280	CASE DEFAULT
281
282	message = "'" // TRIM( output_var%kind ) // "' is not a soil variable"
283	CALL inifor_abort( "main loop", message )
284
285	END SELECT
286	CALL log_runtime( 'time', 'alloc' )
287
288	CALL interpolate_2d( input_buffer(output_var%input_id)%array(:,:,:), &
289	output_arr(:,:,:), output_var%intermediate_grid, output_var )
290	CALL log_runtime( 'time', 'comp' )
291
292
293	!-- Interpolation in 3D, used for atmospheric initial and
294	!-- boundary conditions.
295	CASE ( 'interpolate_3d' )
296
297	ALLOCATE( output_arr(0:output_var%grid % nx, &
298	0:output_var%grid % ny, &
299	1:output_var%grid % nz) )
300
301	CALL log_runtime( 'time', 'alloc' )
302	CALL interpolate_3d( &
303	input_buffer(output_var%input_id)%array(:,:,:), &
304	output_arr(:,:,:), &
305	output_var%intermediate_grid, &
306	output_var%grid)
307	CALL log_runtime( 'time', 'comp' )
308
309	!-- Compute initial avaerage profiles (if --init-mode profile
310	!-- is used)
311	CASE ( 'average profile' )
312
313	ALLOCATE( output_arr(0:output_var%grid%nx, &
314	0:output_var%grid%ny, &
315	1:output_var%grid%nz) )
316	CALL log_runtime( 'time', 'alloc' )
317
318	CALL interp_average_profile( &
319	input_buffer(output_var%input_id)%array(:,:,:), &
320	output_arr(0,0,:), &
321	output_var%averaging_grid )
322
323	IF ( TRIM( output_var%name ) == &
324	'surface_forcing_surface_pressure' ) THEN
325
326	IF ( cfg%p0_is_set ) THEN
327	output_arr(0,0,1) = p0
328	ELSE
329	CALL get_surface_pressure( &
330	output_arr(0,0,:), rho_centre, &
331	output_var%averaging_grid )
332	ENDIF
333
334	ENDIF
335	CALL log_runtime( 'time', 'comp' )
336
337	!-- Compute internal profiles, required for differentiation of
338	!-- geostrophic wind
339	CASE ( 'internal profile' )
340
341	message = "Averaging of internal profile for variable '" //&
342	TRIM( output_var%name ) // "' is not supported."
343
344	SELECT CASE ( TRIM( output_var%name ) )
345
346	CASE( 'internal_density_centre' )
347	ALLOCATE( rho_centre(1:cosmo_grid%nz) )
348	internal_arr => rho_centre
349
350	CASE( 'internal_density_north' )
351	ALLOCATE( rho_north(1:cosmo_grid%nz) )
352	internal_arr => rho_north
353
354	CASE( 'internal_density_south' )
355	ALLOCATE( rho_south(1:cosmo_grid%nz) )
356	internal_arr => rho_south
357
358	CASE( 'internal_density_east' )
359	ALLOCATE( rho_east(1:cosmo_grid%nz) )
360	internal_arr => rho_east
361
362	CASE( 'internal_density_west' )
363	ALLOCATE( rho_west(1:cosmo_grid%nz) )
364	internal_arr => rho_west
365
366	CASE( 'internal_pressure_north' )
367	ALLOCATE( p_north(1:cosmo_grid%nz) )
368	internal_arr => p_north
369
370	CASE( 'internal_pressure_south' )
371	ALLOCATE( p_south(1:cosmo_grid%nz) )
372	internal_arr => p_south
373
374	CASE( 'internal_pressure_east' )
375	ALLOCATE( p_east(1:cosmo_grid%nz) )
376	internal_arr => p_east
377
378	CASE( 'internal_pressure_west' )
379	ALLOCATE( p_west(1:cosmo_grid%nz) )
380	internal_arr => p_west
381
382	CASE DEFAULT
383	CALL inifor_abort( 'main loop', message )
384
385	END SELECT
386	CALL log_runtime( 'time', 'alloc' )
387
388
389	SELECT CASE( TRIM( output_var%name ) )
390
391	CASE( 'internal_pressure_north', &
392	'internal_pressure_south', &
393	'internal_pressure_east', &
394	'internal_pressure_west' )
395
396	CALL average_pressure_perturbation( &
397	input_buffer(output_var%input_id) % array(:,:,:),&
398	internal_arr(:), &
399	cosmo_grid, output_var%averaging_grid &
400	)
401
402	CASE DEFAULT
403
404	CALL average_profile( &
405	input_buffer(output_var%input_id) % array(:,:,:),&
406	internal_arr(:), &
407	output_var%averaging_grid &
408	)
409
410	END SELECT
411
412
413	!
414	!-- Output of geostrophic pressure profiles (with --debug
415	!-- option) is currently deactivated, since they are now
416	!-- defined on averaged COSMO levels instead of PALM levels
417	!-- (requires definiton of COSMO levels in netCDF output.)
418	!IF (.TRUE.) THEN
419	! ALLOCATE( output_arr(1,1,1:output_var%grid % nz) )
420	! output_arr(1,1,:) = internal_arr(:)
421	!ENDIF
422	CALL log_runtime( 'time', 'comp' )
423
424	!
425	!-- This case gets called twice, the first time for ug, the
426	!-- second time for vg. We compute ug and vg at the first call
427	!-- and keep both of them around for the second call.
428	CASE ( 'geostrophic winds' )
429
430	IF (.NOT. ug_vg_have_been_computed ) THEN
431	ALLOCATE( ug_palm(output_var%grid%nz) )
432	ALLOCATE( vg_palm(output_var%grid%nz) )
433	ALLOCATE( ug_cosmo(cosmo_grid%nz) )
434	ALLOCATE( vg_cosmo(cosmo_grid%nz) )
435
436	IF ( cfg%ug_defined_by_user ) THEN
437	ug_palm = cfg%ug
438	vg_palm = cfg%vg
439	ELSE
440	CALL geostrophic_winds( p_north, p_south, p_east, &
441	p_west, rho_centre, f3, &
442	averaging_width_ew, &
443	averaging_width_ns, &
444	phi_n, lambda_n, &
445	phi_centre, lam_centre, &
446	ug_cosmo, vg_cosmo )
447
448	CALL interpolate_1d( ug_cosmo, ug_palm, &
449	output_var%grid )
450
451	CALL interpolate_1d( vg_cosmo, vg_palm, &
452	output_var%grid )
453	ENDIF
454
455	ug_vg_have_been_computed = .TRUE.
456
457	ENDIF
458
459	!
460	!-- Select output array of current geostrophic wind component
461	SELECT CASE( TRIM( output_var%name ) )
462	CASE ( 'ls_forcing_ug' )
463	ug_vg_palm => ug_palm
464	CASE ( 'ls_forcing_vg' )
465	ug_vg_palm => vg_palm
466	END SELECT
467
468	ALLOCATE( output_arr(1,1,output_var%grid%nz) )
469	output_arr(1,1,:) = ug_vg_palm(:)
470
471	!-- User defined constant profiles
472	CASE ( 'set profile' )
473
474	ALLOCATE( output_arr(1,1,1:nz) )
475	CALL log_runtime( 'time', 'alloc' )
476
477	SELECT CASE ( TRIM( output_var%name ) )
478
479	CASE ( 'nudging_tau' )
480	output_arr(1, 1, :) = NUDGING_TAU
481
482	CASE DEFAULT
483	message = "'" // TRIM( output_var%name ) // &
484	"' is not a valid '" // TRIM( output_var%kind ) //&
485	"' variable kind."
486	CALL inifor_abort( 'main loop', message )
487	END SELECT
488	CALL log_runtime( 'time', 'comp' )
489
490	CASE DEFAULT
491	message = "Processing task '" // TRIM( output_var%task ) //&
492	"' not recognized."
493	CALL inifor_abort( '', message )
494
495	END SELECT
496	CALL log_runtime( 'time', 'comp' )
497
498	!------------------------------------------------------------------------------
499	!- Section 2.3: Write current time step of current variable
500	!------------------------------------------------------------------------------
501	!
502	!-- Output of geostrophic pressure profiles (with --debug
503	!-- option) is currently deactivated, since they are now
504	!-- defined on averaged COSMO levels instead of PALM levels
505	!-- (requires definiton of COSMO levels in netCDF output.)
506	!IF (.NOT. output_var%is_internal .OR. debugging_output) THEN
507
508	IF ( .NOT. output_var%is_internal ) THEN
509	message = "Writing variable '" // TRIM( output_var%name ) // "'."
510	CALL report( 'main loop', message )
511	CALL update_output( output_var, output_arr, iter, &
512	output_file, cfg )
513	CALL log_runtime( 'time', 'write' )
514	ENDIF
515
516	IF ( ALLOCATED( output_arr ) ) DEALLOCATE( output_arr )
517	CALL log_runtime( 'time', 'alloc' )
518
519	ENDIF
520
521	!
522	!-- output variable loop
523	ENDDO
524
525	ug_vg_have_been_computed = .FALSE.
526	IF ( group%kind == 'thermodynamics' ) THEN
527	DEALLOCATE( rho_centre )
528	DEALLOCATE( ug_palm )
529	DEALLOCATE( vg_palm )
530	DEALLOCATE( ug_cosmo )
531	DEALLOCATE( vg_cosmo )
532	IF ( .NOT. cfg%ug_defined_by_user ) THEN
533	DEALLOCATE( rho_north )
534	DEALLOCATE( rho_south )
535	DEALLOCATE( rho_east )
536	DEALLOCATE( rho_west )
537	DEALLOCATE( p_north )
538	DEALLOCATE( p_south )
539	DEALLOCATE( p_east )
540	DEALLOCATE( p_west )
541	ENDIF
542	ENDIF
543
544	!
545	!-- Keep input buffer around for averaged (radiation) and
546	!-- accumulated COSMO quantities (precipitation).
547	IF ( group%kind == 'running average' .OR. &
548	group%kind == 'accumulated' ) THEN
549	ELSE
550	CALL report( 'main loop', 'Deallocating input buffer', cfg%debug )
551	DEALLOCATE( input_buffer )
552	ENDIF
553	CALL log_runtime( 'time', 'alloc' )
554
555	!
556	!-- time steps / input files loop
557	ENDDO
558
559	IF ( ALLOCATED( input_buffer ) ) THEN
560	CALL report( 'main loop', 'Deallocating input buffer', cfg%debug )
561	DEALLOCATE( input_buffer )
562	ENDIF
563	CALL log_runtime( 'time', 'alloc' )
564
565	ELSE
566
567	message = "Skipping IO group " // TRIM( str( igroup ) ) // " '" // TRIM( group%kind ) // "'"
568	IF ( ALLOCATED( group%in_var_list ) ) THEN
569	message = TRIM( message ) // " with input variable '" // &
570	TRIM( group%in_var_list(1)%name ) // "'."
571	ENDIF
572
573	CALL report( 'main loop', message, cfg%debug )
574
575	!
576	!-- IO group%to_be_processed conditional
577	ENDIF
578
579	!
580	!-- IO groups loop
581	ENDDO
582
583	!------------------------------------------------------------------------------
584	!- Section 3: Clean up.
585	!------------------------------------------------------------------------------
586	CALL fini_file_lists
587	CALL fini_io_groups
588	CALL fini_variables
589	!CALL fini_grids
590	CALL log_runtime( 'time', 'alloc' )
591	CALL log_runtime( 'report', 'void' )
592
593	message = "Finished writing dynamic driver '" // TRIM( output_file%name )
594
595	IF ( n_wrngs > 0 ) THEN
596	message = TRIM( message ) // "' with " // TRIM( str( n_wrngs ) ) // &
597	" warning(s). Please see logfile '" // LOG_FILE_NAME // "'."
598	ELSE
599	message = TRIM( message ) // "' successfully."
600	END IF
601
602	CALL report( 'main loop', message )
603	CALL close_log
604
605	#else
606
607	USE inifor_control
608	IMPLICIT NONE
609
610	message = "INIFOR was compiled without netCDF support, which is required for it to run. " // &
611	"To use INIFOR, recompile PALM with netCDF support by adding the -D__netcdf " // &
612	"precompiler flag to your .palm.config file."
613	CALL inifor_abort( 'main loop', message )
614
615	#endif
616
617	END PROGRAM inifor

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |