Home

Context Navigation

source: palm/trunk/SOURCE/dynamics_mod.f90 @ 4436

Last change on this file since 4436 was 4360, checked in by suehring, 5 years ago
Bugfix in output of time-averaged plant-canopy quanities; Output of plant-canopy data only where tall canopy is defined; land-surface model: fix wrong location strings; tests: update urban test case; all source code files: copyright update
Property svn:keywords set to `Id`
File size: 59.5 KB

Rev	Line
[4047]	1	!> @file dynamics_mod.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	!
[4360]	17	! Copyright 1997-2020 Leibniz Universitaet Hannover
[4047]	18	!--------------------------------------------------------------------------------------------------!
	19	!
	20	! Current revisions:
	21	! -----------------
[4097]	22	!
[4360]	23	!
[4047]	24	! Former revisions:
	25	! -----------------
	26	! $Id: dynamics_mod.f90 4360 2020-01-07 11:25:50Z oliver.maas $
[4360]	27	! Bugfix for last commit.
	28	!
	29	! 4359 2019-12-30 13:36:50Z suehring
[4358]	30	! Refine post-initialization check for realistically inital values of mixing ratio. Give an error
	31	! message for faulty initial values, but only a warning in a restart run.
	32	!
	33	! 4347 2019-12-18 13:18:33Z suehring
[4347]	34	! Implement post-initialization check for realistically inital values of mixing ratio
	35	!
	36	! 4281 2019-10-29 15:15:39Z schwenkel
[4281]	37	! Moved boundary conditions in dynamics module
	38	!
	39	! 4097 2019-07-15 11:59:11Z suehring
[4097]	40	! Avoid overlong lines - limit is 132 characters per line
	41	!
	42	! 4047 2019-06-21 18:58:09Z knoop
[4047]	43	! Initial introduction of the dynamics module with only dynamics_swap_timelevel implemented
	44	!
	45	!
	46	! Description:
	47	! ------------
	48	!> This module contains the dynamics of PALM.
	49	!--------------------------------------------------------------------------------------------------!
	50	MODULE dynamics_mod
	51
	52
	53	USE arrays_3d, &
[4281]	54	ONLY: c_u, c_u_m, c_u_m_l, c_v, c_v_m, c_v_m_l, c_w, c_w_m, c_w_m_l, &
	55	dzu, &
[4347]	56	exner, &
	57	hyp, &
[4281]	58	pt, pt_1, pt_2, pt_init, pt_p, &
[4047]	59	q, q_1, q_2, q_p, &
	60	s, s_1, s_2, s_p, &
[4281]	61	u, u_1, u_2, u_init, u_p, u_m_l, u_m_n, u_m_r, u_m_s, &
	62	v, v_1, v_2, v_p, v_init, v_m_l, v_m_n, v_m_r, v_m_s, &
	63	w, w_1, w_2, w_p, w_m_l, w_m_n, w_m_r, w_m_s
[4047]	64
[4347]	65	USE basic_constants_and_equations_mod, &
	66	ONLY: magnus, &
	67	rd_d_rv
	68
[4047]	69	USE control_parameters, &
[4281]	70	ONLY: bc_dirichlet_l, &
	71	bc_dirichlet_s, &
	72	bc_radiation_l, &
	73	bc_radiation_n, &
	74	bc_radiation_r, &
	75	bc_radiation_s, &
	76	bc_pt_t_val, &
	77	bc_q_t_val, &
	78	bc_s_t_val, &
	79	child_domain, &
	80	coupling_mode, &
	81	dt_3d, &
	82	ibc_pt_b, &
	83	ibc_pt_t, &
	84	ibc_q_b, &
	85	ibc_q_t, &
	86	ibc_s_b, &
	87	ibc_s_t, &
	88	ibc_uv_b, &
	89	ibc_uv_t, &
[4359]	90	initializing_actions, &
[4281]	91	intermediate_timestep_count, &
	92	length, &
[4347]	93	message_string, &
[4281]	94	nesting_offline, &
	95	nudging, &
[4047]	96	restart_string, &
	97	humidity, &
	98	neutral, &
[4281]	99	passive_scalar, &
	100	tsc, &
	101	use_cmax
[4047]	102
[4281]	103	USE grid_variables, &
	104	ONLY: ddx, &
	105	ddy, &
	106	dx, &
	107	dy
	108
[4047]	109	USE indices, &
	110	ONLY: nbgp, &
[4281]	111	nx, &
[4047]	112	nxl, &
[4281]	113	nxlg, &
[4047]	114	nxr, &
[4281]	115	nxrg, &
	116	ny, &
[4047]	117	nys, &
[4281]	118	nysg, &
[4047]	119	nyn, &
[4281]	120	nyng, &
[4047]	121	nzb, &
	122	nzt
	123
	124	USE kinds
	125
[4281]	126	USE pegrid
	127
	128	USE pmc_interface, &
	129	ONLY : nesting_mode
	130
	131	USE surface_mod, &
	132	ONLY : bc_h
	133
	134
[4047]	135	IMPLICIT NONE
	136
	137	LOGICAL :: dynamics_module_enabled = .FALSE. !<
	138
	139	SAVE
	140
	141	PRIVATE
	142
	143	!
	144	!-- Public functions
	145	PUBLIC &
	146	dynamics_parin, &
	147	dynamics_check_parameters, &
	148	dynamics_check_data_output_ts, &
	149	dynamics_check_data_output_pr, &
	150	dynamics_check_data_output, &
	151	dynamics_init_masks, &
	152	dynamics_define_netcdf_grid, &
	153	dynamics_init_arrays, &
	154	dynamics_init, &
	155	dynamics_init_checks, &
	156	dynamics_header, &
	157	dynamics_actions, &
	158	dynamics_non_advective_processes, &
	159	dynamics_exchange_horiz, &
	160	dynamics_prognostic_equations, &
[4281]	161	dynamics_boundary_conditions, &
[4047]	162	dynamics_swap_timelevel, &
	163	dynamics_3d_data_averaging, &
	164	dynamics_data_output_2d, &
	165	dynamics_data_output_3d, &
	166	dynamics_statistics, &
	167	dynamics_rrd_global, &
	168	dynamics_rrd_local, &
	169	dynamics_wrd_global, &
	170	dynamics_wrd_local, &
	171	dynamics_last_actions
	172
	173	!
	174	!-- Public parameters, constants and initial values
	175	PUBLIC &
	176	dynamics_module_enabled
	177
	178	INTERFACE dynamics_parin
	179	MODULE PROCEDURE dynamics_parin
	180	END INTERFACE dynamics_parin
	181
	182	INTERFACE dynamics_check_parameters
	183	MODULE PROCEDURE dynamics_check_parameters
	184	END INTERFACE dynamics_check_parameters
	185
	186	INTERFACE dynamics_check_data_output_ts
	187	MODULE PROCEDURE dynamics_check_data_output_ts
	188	END INTERFACE dynamics_check_data_output_ts
	189
	190	INTERFACE dynamics_check_data_output_pr
	191	MODULE PROCEDURE dynamics_check_data_output_pr
	192	END INTERFACE dynamics_check_data_output_pr
	193
	194	INTERFACE dynamics_check_data_output
	195	MODULE PROCEDURE dynamics_check_data_output
	196	END INTERFACE dynamics_check_data_output
	197
	198	INTERFACE dynamics_init_masks
	199	MODULE PROCEDURE dynamics_init_masks
	200	END INTERFACE dynamics_init_masks
	201
	202	INTERFACE dynamics_define_netcdf_grid
	203	MODULE PROCEDURE dynamics_define_netcdf_grid
	204	END INTERFACE dynamics_define_netcdf_grid
	205
	206	INTERFACE dynamics_init_arrays
	207	MODULE PROCEDURE dynamics_init_arrays
	208	END INTERFACE dynamics_init_arrays
	209
	210	INTERFACE dynamics_init
	211	MODULE PROCEDURE dynamics_init
	212	END INTERFACE dynamics_init
	213
	214	INTERFACE dynamics_init_checks
	215	MODULE PROCEDURE dynamics_init_checks
	216	END INTERFACE dynamics_init_checks
	217
	218	INTERFACE dynamics_header
	219	MODULE PROCEDURE dynamics_header
	220	END INTERFACE dynamics_header
	221
	222	INTERFACE dynamics_actions
	223	MODULE PROCEDURE dynamics_actions
	224	MODULE PROCEDURE dynamics_actions_ij
	225	END INTERFACE dynamics_actions
	226
	227	INTERFACE dynamics_non_advective_processes
	228	MODULE PROCEDURE dynamics_non_advective_processes
	229	MODULE PROCEDURE dynamics_non_advective_processes_ij
	230	END INTERFACE dynamics_non_advective_processes
	231
	232	INTERFACE dynamics_exchange_horiz
	233	MODULE PROCEDURE dynamics_exchange_horiz
	234	END INTERFACE dynamics_exchange_horiz
	235
	236	INTERFACE dynamics_prognostic_equations
	237	MODULE PROCEDURE dynamics_prognostic_equations
	238	MODULE PROCEDURE dynamics_prognostic_equations_ij
	239	END INTERFACE dynamics_prognostic_equations
	240
[4281]	241	INTERFACE dynamics_boundary_conditions
	242	MODULE PROCEDURE dynamics_boundary_conditions
	243	END INTERFACE dynamics_boundary_conditions
	244
[4047]	245	INTERFACE dynamics_swap_timelevel
	246	MODULE PROCEDURE dynamics_swap_timelevel
	247	END INTERFACE dynamics_swap_timelevel
	248
	249	INTERFACE dynamics_3d_data_averaging
	250	MODULE PROCEDURE dynamics_3d_data_averaging
	251	END INTERFACE dynamics_3d_data_averaging
	252
	253	INTERFACE dynamics_data_output_2d
	254	MODULE PROCEDURE dynamics_data_output_2d
	255	END INTERFACE dynamics_data_output_2d
	256
	257	INTERFACE dynamics_data_output_3d
	258	MODULE PROCEDURE dynamics_data_output_3d
	259	END INTERFACE dynamics_data_output_3d
	260
	261	INTERFACE dynamics_statistics
	262	MODULE PROCEDURE dynamics_statistics
	263	END INTERFACE dynamics_statistics
	264
	265	INTERFACE dynamics_rrd_global
	266	MODULE PROCEDURE dynamics_rrd_global
	267	END INTERFACE dynamics_rrd_global
	268
	269	INTERFACE dynamics_rrd_local
	270	MODULE PROCEDURE dynamics_rrd_local
	271	END INTERFACE dynamics_rrd_local
	272
	273	INTERFACE dynamics_wrd_global
	274	MODULE PROCEDURE dynamics_wrd_global
	275	END INTERFACE dynamics_wrd_global
	276
	277	INTERFACE dynamics_wrd_local
	278	MODULE PROCEDURE dynamics_wrd_local
	279	END INTERFACE dynamics_wrd_local
	280
	281	INTERFACE dynamics_last_actions
	282	MODULE PROCEDURE dynamics_last_actions
	283	END INTERFACE dynamics_last_actions
	284
	285
	286	CONTAINS
	287
	288
	289	!--------------------------------------------------------------------------------------------------!
	290	! Description:
	291	! ------------
	292	!> Read module-specific namelist
	293	!--------------------------------------------------------------------------------------------------!
	294	SUBROUTINE dynamics_parin
	295
	296
	297	CHARACTER (LEN=80) :: line !< dummy string that contains the current line of the parameter file
	298
	299	NAMELIST /dynamics_parameters/ &
	300	dynamics_module_enabled
	301
	302	line = ' '
	303	!
	304	!-- Try to find module-specific namelist
	305	REWIND ( 11 )
	306	line = ' '
	307	DO WHILE ( INDEX( line, '&dynamics_parameters' ) == 0 )
	308	READ ( 11, '(A)', END=12 ) line
	309	ENDDO
	310	BACKSPACE ( 11 )
	311
	312	!-- Set default module switch to true
	313	dynamics_module_enabled = .TRUE.
	314
	315	!-- Read user-defined namelist
	316	READ ( 11, dynamics_parameters, ERR = 10 )
	317
	318	GOTO 12
	319
	320	10 BACKSPACE( 11 )
	321	READ( 11 , '(A)') line
	322	CALL parin_fail_message( 'dynamics_parameters', line )
	323
	324	12 CONTINUE
	325
	326	END SUBROUTINE dynamics_parin
	327
	328
	329	!--------------------------------------------------------------------------------------------------!
	330	! Description:
	331	! ------------
	332	!> Check control parameters and deduce further quantities.
	333	!--------------------------------------------------------------------------------------------------!
	334	SUBROUTINE dynamics_check_parameters
	335
	336
	337	END SUBROUTINE dynamics_check_parameters
	338
	339
	340	!--------------------------------------------------------------------------------------------------!
	341	! Description:
	342	! ------------
	343	!> Set module-specific timeseries units and labels
	344	!--------------------------------------------------------------------------------------------------!
	345	SUBROUTINE dynamics_check_data_output_ts( dots_max, dots_num, dots_label, dots_unit )
	346
	347
	348	INTEGER(iwp), INTENT(IN) :: dots_max
	349	INTEGER(iwp), INTENT(INOUT) :: dots_num
	350	CHARACTER (LEN=*), DIMENSION(dots_max), INTENT(INOUT) :: dots_label
	351	CHARACTER (LEN=*), DIMENSION(dots_max), INTENT(INOUT) :: dots_unit
	352
	353	!
	354	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	355	IF ( dots_num == 0 .OR. dots_label(1)(1:1) == ' ' .OR. dots_unit(1)(1:1) == ' ' ) CONTINUE
	356
	357
	358	END SUBROUTINE dynamics_check_data_output_ts
	359
	360
	361	!--------------------------------------------------------------------------------------------------!
	362	! Description:
	363	! ------------
	364	!> Set the unit of module-specific profile output quantities. For those variables not recognized,
	365	!> the parameter unit is set to "illegal", which tells the calling routine that the output variable
	366	!> is not defined and leads to a program abort.
	367	!--------------------------------------------------------------------------------------------------!
	368	SUBROUTINE dynamics_check_data_output_pr( variable, var_count, unit, dopr_unit )
	369
	370
	371	CHARACTER (LEN=*) :: unit !<
	372	CHARACTER (LEN=*) :: variable !<
	373	CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit
	374
	375	INTEGER(iwp) :: var_count !<
	376
	377	!
	378	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	379	IF ( unit(1:1) == ' ' .OR. dopr_unit(1:1) == ' ' .OR. var_count == 0 ) CONTINUE
	380
	381	SELECT CASE ( TRIM( variable ) )
	382
	383	! CASE ( 'var_name' )
	384
	385	CASE DEFAULT
	386	unit = 'illegal'
	387
	388	END SELECT
	389
	390
	391	END SUBROUTINE dynamics_check_data_output_pr
	392
	393
	394	!--------------------------------------------------------------------------------------------------!
	395	! Description:
	396	! ------------
	397	!> Set the unit of module-specific output quantities. For those variables not recognized,
	398	!> the parameter unit is set to "illegal", which tells the calling routine that the output variable
	399	!< is not defined and leads to a program abort.
	400	!--------------------------------------------------------------------------------------------------!
	401	SUBROUTINE dynamics_check_data_output( variable, unit )
	402
	403
	404	CHARACTER (LEN=*) :: unit !<
	405	CHARACTER (LEN=*) :: variable !<
	406
	407	SELECT CASE ( TRIM( variable ) )
	408
	409	! CASE ( 'u2' )
	410
	411	CASE DEFAULT
	412	unit = 'illegal'
	413
	414	END SELECT
	415
	416
	417	END SUBROUTINE dynamics_check_data_output
	418
	419
	420	!------------------------------------------------------------------------------!
	421	!
	422	! Description:
	423	! ------------
	424	!> Initialize module-specific masked output
	425	!------------------------------------------------------------------------------!
	426	SUBROUTINE dynamics_init_masks( variable, unit )
	427
	428
	429	CHARACTER (LEN=*) :: unit !<
	430	CHARACTER (LEN=*) :: variable !<
	431
	432
	433	SELECT CASE ( TRIM( variable ) )
	434
	435	! CASE ( 'u2' )
	436
	437	CASE DEFAULT
	438	unit = 'illegal'
	439
	440	END SELECT
	441
	442
	443	END SUBROUTINE dynamics_init_masks
	444
	445
	446	!--------------------------------------------------------------------------------------------------!
	447	! Description:
	448	! ------------
	449	!> Initialize module-specific arrays
	450	!--------------------------------------------------------------------------------------------------!
	451	SUBROUTINE dynamics_init_arrays
	452
	453
	454	END SUBROUTINE dynamics_init_arrays
	455
	456
	457	!--------------------------------------------------------------------------------------------------!
	458	! Description:
	459	! ------------
	460	!> Execution of module-specific initializing actions
	461	!--------------------------------------------------------------------------------------------------!
	462	SUBROUTINE dynamics_init
	463
	464
	465	END SUBROUTINE dynamics_init
	466
	467
	468	!--------------------------------------------------------------------------------------------------!
	469	! Description:
	470	! ------------
	471	!> Perform module-specific post-initialization checks
	472	!--------------------------------------------------------------------------------------------------!
	473	SUBROUTINE dynamics_init_checks
	474
[4347]	475	INTEGER(iwp) :: i !< loop index in x-direction
	476	INTEGER(iwp) :: j !< loop index in y-direction
	477	INTEGER(iwp) :: k !< loop index in z-direction
[4047]	478
[4347]	479	LOGICAL :: realistic_q = .TRUE. !< flag indicating realistic mixing ratios
	480
	481	REAL(wp) :: e_s !< saturation water vapor pressure
	482	REAL(wp) :: q_s !< saturation mixing ratio
	483	REAL(wp) :: t_l !< actual temperature
	484
	485	!
	486	!-- Check for realistic initial mixing ratio. This must be in a realistic phyiscial range and must
	487	!-- not exceed the saturation mixing ratio by more than 2 percent. Please note, the check is
	488	!-- performed for each grid point (not just for a vertical profile), in order to cover also
[4358]	489	!-- three-dimensional initialization. Note, this check gives an error only for the initial run not
	490	!-- for a restart run. In case there are no cloud physics considered, the mixing ratio can exceed
	491	!-- the saturation moisture. This case a warning is given.
[4347]	492	IF ( humidity .AND. .NOT. neutral ) THEN
	493	DO i = nxl, nxr
	494	DO j = nys, nyn
	495	DO k = nzb+1, nzt
	496	!
	497	!-- Calculate actual temperature, water vapor saturation pressure, and based on this
	498	!-- the saturation mixing ratio.
	499	t_l = exner(k) * pt(k,j,i)
	500	e_s = magnus( t_l )
	501	q_s = rd_d_rv * e_s / ( hyp(k) - e_s )
	502
	503	IF ( q(k,j,i) > 1.02_wp * q_s ) realistic_q = .FALSE.
	504	ENDDO
	505	ENDDO
	506	ENDDO
	507	!
	508	!-- Since the check is performed locally, merge the logical flag from all mpi ranks,
	509	!-- in order to do not print the error message multiple times.
	510	#if defined( __parallel )
	511	CALL MPI_ALLREDUCE( MPI_IN_PLACE, realistic_q, 1, MPI_LOGICAL, MPI_LAND, comm2d, ierr)
	512	#endif
	513
[4358]	514	IF ( .NOT. realistic_q .AND. &
	515	TRIM( initializing_actions ) /= 'read_restart_data' ) THEN
[4347]	516	message_string = 'The initial mixing ratio exceeds the saturation mixing ratio.'
	517	CALL message( 'dynamic_init_checks', 'PA0697', 2, 2, 0, 6, 0 )
[4358]	518	ELSEIF ( .NOT. realistic_q .AND. &
	519	TRIM( initializing_actions ) == 'read_restart_data' ) THEN
	520	message_string = 'The mixing ratio exceeds the saturation mixing ratio.'
	521	CALL message( 'dynamic_init_checks', 'PA0697', 0, 1, 0, 6, 0 )
[4347]	522	ENDIF
	523	ENDIF
	524
[4047]	525	END SUBROUTINE dynamics_init_checks
	526
	527
	528	!--------------------------------------------------------------------------------------------------!
	529	! Description:
	530	! ------------
	531	!> Set the grids on which module-specific output quantities are defined. Allowed values for
	532	!> grid_x are "x" and "xu", for grid_y "y" and "yv", and for grid_z "zu" and "zw".
	533	!--------------------------------------------------------------------------------------------------!
	534	SUBROUTINE dynamics_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
	535
	536
	537	CHARACTER (LEN=*) :: grid_x !<
	538	CHARACTER (LEN=*) :: grid_y !<
	539	CHARACTER (LEN=*) :: grid_z !<
	540	CHARACTER (LEN=*) :: variable !<
	541
	542	LOGICAL :: found !<
	543
	544
	545	SELECT CASE ( TRIM( variable ) )
	546
	547	! CASE ( 'u2' )
	548
	549	CASE DEFAULT
	550	found = .FALSE.
	551	grid_x = 'none'
	552	grid_y = 'none'
	553	grid_z = 'none'
	554
	555	END SELECT
	556
	557
	558	END SUBROUTINE dynamics_define_netcdf_grid
	559
	560
	561	!--------------------------------------------------------------------------------------------------!
	562	! Description:
	563	! ------------
	564	!> Print a header with module-specific information.
	565	!--------------------------------------------------------------------------------------------------!
	566	SUBROUTINE dynamics_header( io )
	567
	568
	569	INTEGER(iwp) :: io !<
	570
	571	!
	572	!-- If no module-specific variables are read from the namelist-file, no information will be printed.
	573	IF ( .NOT. dynamics_module_enabled ) THEN
	574	WRITE ( io, 100 )
	575	RETURN
	576	ENDIF
	577
	578	!
	579	!-- Printing the information.
	580	WRITE ( io, 110 )
	581
	582	!
	583	!-- Format-descriptors
	584	100 FORMAT (//' *** dynamic module disabled'/)
	585	110 FORMAT (//1X,78('#') &
	586	//' User-defined variables and actions:'/ &
	587	' -----------------------------------'//)
	588
	589	END SUBROUTINE dynamics_header
	590
	591
	592	!--------------------------------------------------------------------------------------------------!
	593	! Description:
	594	! ------------
	595	!> Execute module-specific actions for all grid points
	596	!--------------------------------------------------------------------------------------------------!
	597	SUBROUTINE dynamics_actions( location )
	598
	599
	600	CHARACTER (LEN=*) :: location !<
	601
	602	! INTEGER(iwp) :: i !<
	603	! INTEGER(iwp) :: j !<
	604	! INTEGER(iwp) :: k !<
	605
	606	!
	607	!-- Here the user-defined actions follow
	608	!-- No calls for single grid points are allowed at locations before and
	609	!-- after the timestep, since these calls are not within an i,j-loop
	610	SELECT CASE ( location )
	611
	612	CASE ( 'before_timestep' )
	613
	614
	615	CASE ( 'before_prognostic_equations' )
	616
	617
	618	CASE ( 'after_integration' )
	619
	620
	621	CASE ( 'after_timestep' )
	622
	623
	624	CASE ( 'u-tendency' )
	625
	626
	627	CASE ( 'v-tendency' )
	628
	629
	630	CASE ( 'w-tendency' )
	631
	632
	633	CASE ( 'pt-tendency' )
	634
	635
	636	CASE ( 'sa-tendency' )
	637
	638
	639	CASE ( 'e-tendency' )
	640
	641
	642	CASE ( 'q-tendency' )
	643
	644
	645	CASE ( 's-tendency' )
	646
	647
	648	CASE DEFAULT
	649	CONTINUE
	650
	651	END SELECT
	652
	653	END SUBROUTINE dynamics_actions
	654
	655
	656	!--------------------------------------------------------------------------------------------------!
	657	! Description:
	658	! ------------
	659	!> Execute module-specific actions for grid point i,j
	660	!--------------------------------------------------------------------------------------------------!
	661	SUBROUTINE dynamics_actions_ij( i, j, location )
	662
	663
	664	CHARACTER (LEN=*) :: location
	665
	666	INTEGER(iwp) :: i
	667	INTEGER(iwp) :: j
	668
	669	!
	670	!-- Here the user-defined actions follow
	671	SELECT CASE ( location )
	672
	673	CASE ( 'u-tendency' )
	674
	675	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	676	IF ( i + j < 0 ) CONTINUE
	677
	678	CASE ( 'v-tendency' )
	679
	680
	681	CASE ( 'w-tendency' )
	682
	683
	684	CASE ( 'pt-tendency' )
	685
	686
	687	CASE ( 'sa-tendency' )
	688
	689
	690	CASE ( 'e-tendency' )
	691
	692
	693	CASE ( 'q-tendency' )
	694
	695
	696	CASE ( 's-tendency' )
	697
	698
	699	CASE DEFAULT
	700	CONTINUE
	701
	702	END SELECT
	703
	704	END SUBROUTINE dynamics_actions_ij
	705
	706
	707	!--------------------------------------------------------------------------------------------------!
	708	! Description:
	709	! ------------
	710	!> Compute module-specific non-advective processes for all grid points
	711	!--------------------------------------------------------------------------------------------------!
	712	SUBROUTINE dynamics_non_advective_processes
	713
	714
	715
	716	END SUBROUTINE dynamics_non_advective_processes
	717
	718
	719	!--------------------------------------------------------------------------------------------------!
	720	! Description:
	721	! ------------
	722	!> Compute module-specific non-advective processes for grid points i,j
	723	!--------------------------------------------------------------------------------------------------!
	724	SUBROUTINE dynamics_non_advective_processes_ij( i, j )
	725
	726
	727	INTEGER(iwp) :: i !<
	728	INTEGER(iwp) :: j !<
	729
	730	!
	731	!-- Next line is just to avoid compiler warnings about unused variables. You may remove it.
	732	IF ( i + j < 0 ) CONTINUE
	733
	734
	735	END SUBROUTINE dynamics_non_advective_processes_ij
	736
	737
	738	!--------------------------------------------------------------------------------------------------!
	739	! Description:
	740	! ------------
	741	!> Perform module-specific horizontal boundary exchange
	742	!--------------------------------------------------------------------------------------------------!
	743	SUBROUTINE dynamics_exchange_horiz
	744
	745
	746
	747	END SUBROUTINE dynamics_exchange_horiz
	748
	749
	750	!--------------------------------------------------------------------------------------------------!
	751	! Description:
	752	! ------------
	753	!> Compute module-specific prognostic equations for all grid points
	754	!--------------------------------------------------------------------------------------------------!
	755	SUBROUTINE dynamics_prognostic_equations
	756
	757
	758
	759	END SUBROUTINE dynamics_prognostic_equations
	760
	761
	762	!--------------------------------------------------------------------------------------------------!
	763	! Description:
	764	! ------------
	765	!> Compute module-specific prognostic equations for grid point i,j
	766	!--------------------------------------------------------------------------------------------------!
	767	SUBROUTINE dynamics_prognostic_equations_ij( i, j, i_omp_start, tn )
	768
	769
	770	INTEGER(iwp), INTENT(IN) :: i !< grid index in x-direction
	771	INTEGER(iwp), INTENT(IN) :: j !< grid index in y-direction
	772	INTEGER(iwp), INTENT(IN) :: i_omp_start !< first loop index of i-loop in prognostic_equations
	773	INTEGER(iwp), INTENT(IN) :: tn !< task number of openmp task
	774
	775	!
	776	!-- Next line is just to avoid compiler warnings about unused variables. You may remove it.
	777	IF ( i + j + i_omp_start + tn < 0 ) CONTINUE
	778
	779	END SUBROUTINE dynamics_prognostic_equations_ij
	780
	781
[4281]	782	!--------------------------------------------------------------------------------------------------!
	783	! Description:
	784	! ------------
	785	!> Compute boundary conditions of dynamics model
	786	!--------------------------------------------------------------------------------------------------!
	787	SUBROUTINE dynamics_boundary_conditions
	788
	789	IMPLICIT NONE
	790
	791	INTEGER(iwp) :: i !< grid index x direction
	792	INTEGER(iwp) :: j !< grid index y direction
	793	INTEGER(iwp) :: k !< grid index z direction
	794	INTEGER(iwp) :: l !< running index boundary type, for up- and downward-facing walls
	795	INTEGER(iwp) :: m !< running index surface elements
	796
	797	REAL(wp) :: c_max !< maximum phase velocity allowed by CFL criterion, used for outflow boundary condition
	798	REAL(wp) :: denom !< horizontal gradient of velocity component normal to the outflow boundary
	799
	800	!
	801	!-- Bottom boundary
	802	IF ( ibc_uv_b == 1 ) THEN
	803	u_p(nzb,:,:) = u_p(nzb+1,:,:)
	804	v_p(nzb,:,:) = v_p(nzb+1,:,:)
	805	ENDIF
	806	!
	807	!-- Set zero vertical velocity at topography top (l=0), or bottom (l=1) in case
	808	!-- of downward-facing surfaces.
	809	DO l = 0, 1
	810	!$OMP PARALLEL DO PRIVATE( i, j, k )
	811	!$ACC PARALLEL LOOP PRIVATE(i, j, k) &
	812	!$ACC PRESENT(bc_h, w_p)
	813	DO m = 1, bc_h(l)%ns
	814	i = bc_h(l)%i(m)
	815	j = bc_h(l)%j(m)
	816	k = bc_h(l)%k(m)
	817	w_p(k+bc_h(l)%koff,j,i) = 0.0_wp
	818	ENDDO
	819	ENDDO
	820
	821	!
	822	!-- Top boundary. A nested domain ( ibc_uv_t = 3 ) does not require settings.
	823	IF ( ibc_uv_t == 0 ) THEN
	824	!$ACC KERNELS PRESENT(u_p, v_p, u_init, v_init)
	825	u_p(nzt+1,:,:) = u_init(nzt+1)
	826	v_p(nzt+1,:,:) = v_init(nzt+1)
	827	!$ACC END KERNELS
	828	ELSEIF ( ibc_uv_t == 1 ) THEN
	829	u_p(nzt+1,:,:) = u_p(nzt,:,:)
	830	v_p(nzt+1,:,:) = v_p(nzt,:,:)
	831	ENDIF
	832
	833	!
	834	!-- Vertical nesting: Vertical velocity not zero at the top of the fine grid
	835	IF ( .NOT. child_domain .AND. .NOT. nesting_offline .AND. &
	836	TRIM(coupling_mode) /= 'vnested_fine' ) THEN
	837	!$ACC KERNELS PRESENT(w_p)
	838	w_p(nzt:nzt+1,:,:) = 0.0_wp !< nzt is not a prognostic level (but cf. pres)
	839	!$ACC END KERNELS
	840	ENDIF
	841
	842	!
	843	!-- Temperature at bottom and top boundary.
	844	!-- In case of coupled runs (ibc_pt_b = 2) the temperature is given by
	845	!-- the sea surface temperature of the coupled ocean model.
	846	!-- Dirichlet
	847	IF ( .NOT. neutral ) THEN
	848	IF ( ibc_pt_b == 0 ) THEN
	849	DO l = 0, 1
	850	!$OMP PARALLEL DO PRIVATE( i, j, k )
	851	DO m = 1, bc_h(l)%ns
	852	i = bc_h(l)%i(m)
	853	j = bc_h(l)%j(m)
	854	k = bc_h(l)%k(m)
	855	pt_p(k+bc_h(l)%koff,j,i) = pt(k+bc_h(l)%koff,j,i)
	856	ENDDO
	857	ENDDO
	858	!
	859	!-- Neumann, zero-gradient
	860	ELSEIF ( ibc_pt_b == 1 ) THEN
	861	DO l = 0, 1
	862	!$OMP PARALLEL DO PRIVATE( i, j, k )
	863	!$ACC PARALLEL LOOP PRIVATE(i, j, k) &
	864	!$ACC PRESENT(bc_h, pt_p)
	865	DO m = 1, bc_h(l)%ns
	866	i = bc_h(l)%i(m)
	867	j = bc_h(l)%j(m)
	868	k = bc_h(l)%k(m)
	869	pt_p(k+bc_h(l)%koff,j,i) = pt_p(k,j,i)
	870	ENDDO
	871	ENDDO
	872	ENDIF
	873
	874	!
	875	!-- Temperature at top boundary
	876	IF ( ibc_pt_t == 0 ) THEN
	877	pt_p(nzt+1,:,:) = pt(nzt+1,:,:)
	878	!
	879	!-- In case of nudging adjust top boundary to pt which is
	880	!-- read in from NUDGING-DATA
	881	IF ( nudging ) THEN
	882	pt_p(nzt+1,:,:) = pt_init(nzt+1)
	883	ENDIF
	884	ELSEIF ( ibc_pt_t == 1 ) THEN
	885	pt_p(nzt+1,:,:) = pt_p(nzt,:,:)
	886	ELSEIF ( ibc_pt_t == 2 ) THEN
	887	!$ACC KERNELS PRESENT(pt_p, dzu)
	888	pt_p(nzt+1,:,:) = pt_p(nzt,:,:) + bc_pt_t_val * dzu(nzt+1)
	889	!$ACC END KERNELS
	890	ENDIF
	891	ENDIF
	892	!
	893	!-- Boundary conditions for total water content,
	894	!-- bottom and top boundary (see also temperature)
	895	IF ( humidity ) THEN
	896	!
	897	!-- Surface conditions for constant_humidity_flux
	898	!-- Run loop over all non-natural and natural walls. Note, in wall-datatype
	899	!-- the k coordinate belongs to the atmospheric grid point, therefore, set
	900	!-- q_p at k-1
	901	IF ( ibc_q_b == 0 ) THEN
	902
	903	DO l = 0, 1
	904	!$OMP PARALLEL DO PRIVATE( i, j, k )
	905	DO m = 1, bc_h(l)%ns
	906	i = bc_h(l)%i(m)
	907	j = bc_h(l)%j(m)
	908	k = bc_h(l)%k(m)
	909	q_p(k+bc_h(l)%koff,j,i) = q(k+bc_h(l)%koff,j,i)
	910	ENDDO
	911	ENDDO
	912
	913	ELSE
	914
	915	DO l = 0, 1
	916	!$OMP PARALLEL DO PRIVATE( i, j, k )
	917	DO m = 1, bc_h(l)%ns
	918	i = bc_h(l)%i(m)
	919	j = bc_h(l)%j(m)
	920	k = bc_h(l)%k(m)
	921	q_p(k+bc_h(l)%koff,j,i) = q_p(k,j,i)
	922	ENDDO
	923	ENDDO
	924	ENDIF
	925	!
	926	!-- Top boundary
	927	IF ( ibc_q_t == 0 ) THEN
	928	q_p(nzt+1,:,:) = q(nzt+1,:,:)
	929	ELSEIF ( ibc_q_t == 1 ) THEN
	930	q_p(nzt+1,:,:) = q_p(nzt,:,:) + bc_q_t_val * dzu(nzt+1)
	931	ENDIF
	932	ENDIF
	933	!
	934	!-- Boundary conditions for scalar,
	935	!-- bottom and top boundary (see also temperature)
	936	IF ( passive_scalar ) THEN
	937	!
	938	!-- Surface conditions for constant_humidity_flux
	939	!-- Run loop over all non-natural and natural walls. Note, in wall-datatype
	940	!-- the k coordinate belongs to the atmospheric grid point, therefore, set
	941	!-- s_p at k-1
	942	IF ( ibc_s_b == 0 ) THEN
	943
	944	DO l = 0, 1
	945	!$OMP PARALLEL DO PRIVATE( i, j, k )
	946	DO m = 1, bc_h(l)%ns
	947	i = bc_h(l)%i(m)
	948	j = bc_h(l)%j(m)
	949	k = bc_h(l)%k(m)
	950	s_p(k+bc_h(l)%koff,j,i) = s(k+bc_h(l)%koff,j,i)
	951	ENDDO
	952	ENDDO
	953
	954	ELSE
	955
	956	DO l = 0, 1
	957	!$OMP PARALLEL DO PRIVATE( i, j, k )
	958	DO m = 1, bc_h(l)%ns
	959	i = bc_h(l)%i(m)
	960	j = bc_h(l)%j(m)
	961	k = bc_h(l)%k(m)
	962	s_p(k+bc_h(l)%koff,j,i) = s_p(k,j,i)
	963	ENDDO
	964	ENDDO
	965	ENDIF
	966	!
	967	!-- Top boundary condition
	968	IF ( ibc_s_t == 0 ) THEN
	969	s_p(nzt+1,:,:) = s(nzt+1,:,:)
	970	ELSEIF ( ibc_s_t == 1 ) THEN
	971	s_p(nzt+1,:,:) = s_p(nzt,:,:)
	972	ELSEIF ( ibc_s_t == 2 ) THEN
	973	s_p(nzt+1,:,:) = s_p(nzt,:,:) + bc_s_t_val * dzu(nzt+1)
	974	ENDIF
	975
	976	ENDIF
	977	!
	978	!-- In case of inflow or nest boundary at the south boundary the boundary for v
	979	!-- is at nys and in case of inflow or nest boundary at the left boundary the
	980	!-- boundary for u is at nxl. Since in prognostic_equations (cache optimized
	981	!-- version) these levels are handled as a prognostic level, boundary values
	982	!-- have to be restored here.
	983	IF ( bc_dirichlet_s ) THEN
	984	v_p(:,nys,:) = v_p(:,nys-1,:)
	985	ELSEIF ( bc_dirichlet_l ) THEN
	986	u_p(:,:,nxl) = u_p(:,:,nxl-1)
	987	ENDIF
	988
	989	!
	990	!-- The same restoration for u at i=nxl and v at j=nys as above must be made
	991	!-- in case of nest boundaries. This must not be done in case of vertical nesting
	992	!-- mode as in that case the lateral boundaries are actually cyclic.
	993	!-- Lateral oundary conditions for TKE and dissipation are set
	994	!-- in tcm_boundary_conds.
	995	IF ( nesting_mode /= 'vertical' .OR. nesting_offline ) THEN
	996	IF ( bc_dirichlet_s ) THEN
	997	v_p(:,nys,:) = v_p(:,nys-1,:)
	998	ENDIF
	999	IF ( bc_dirichlet_l ) THEN
	1000	u_p(:,:,nxl) = u_p(:,:,nxl-1)
	1001	ENDIF
	1002	ENDIF
	1003
	1004	!
	1005	!-- Lateral boundary conditions for scalar quantities at the outflow.
	1006	!-- Lateral oundary conditions for TKE and dissipation are set
	1007	!-- in tcm_boundary_conds.
	1008	IF ( bc_radiation_s ) THEN
	1009	pt_p(:,nys-1,:) = pt_p(:,nys,:)
	1010	IF ( humidity ) THEN
	1011	q_p(:,nys-1,:) = q_p(:,nys,:)
	1012	ENDIF
	1013	IF ( passive_scalar ) s_p(:,nys-1,:) = s_p(:,nys,:)
	1014	ELSEIF ( bc_radiation_n ) THEN
	1015	pt_p(:,nyn+1,:) = pt_p(:,nyn,:)
	1016	IF ( humidity ) THEN
	1017	q_p(:,nyn+1,:) = q_p(:,nyn,:)
	1018	ENDIF
	1019	IF ( passive_scalar ) s_p(:,nyn+1,:) = s_p(:,nyn,:)
	1020	ELSEIF ( bc_radiation_l ) THEN
	1021	pt_p(:,:,nxl-1) = pt_p(:,:,nxl)
	1022	IF ( humidity ) THEN
	1023	q_p(:,:,nxl-1) = q_p(:,:,nxl)
	1024	ENDIF
	1025	IF ( passive_scalar ) s_p(:,:,nxl-1) = s_p(:,:,nxl)
	1026	ELSEIF ( bc_radiation_r ) THEN
	1027	pt_p(:,:,nxr+1) = pt_p(:,:,nxr)
	1028	IF ( humidity ) THEN
	1029	q_p(:,:,nxr+1) = q_p(:,:,nxr)
	1030	ENDIF
	1031	IF ( passive_scalar ) s_p(:,:,nxr+1) = s_p(:,:,nxr)
	1032	ENDIF
	1033
	1034	!
	1035	!-- Radiation boundary conditions for the velocities at the respective outflow.
	1036	!-- The phase velocity is either assumed to the maximum phase velocity that
	1037	!-- ensures numerical stability (CFL-condition) or calculated after
	1038	!-- Orlanski(1976) and averaged along the outflow boundary.
	1039	IF ( bc_radiation_s ) THEN
	1040
	1041	IF ( use_cmax ) THEN
	1042	u_p(:,-1,:) = u(:,0,:)
	1043	v_p(:,0,:) = v(:,1,:)
	1044	w_p(:,-1,:) = w(:,0,:)
	1045	ELSEIF ( .NOT. use_cmax ) THEN
	1046
	1047	c_max = dy / dt_3d
	1048
	1049	c_u_m_l = 0.0_wp
	1050	c_v_m_l = 0.0_wp
	1051	c_w_m_l = 0.0_wp
	1052
	1053	c_u_m = 0.0_wp
	1054	c_v_m = 0.0_wp
	1055	c_w_m = 0.0_wp
	1056
	1057	!
	1058	!-- Calculate the phase speeds for u, v, and w, first local and then
	1059	!-- average along the outflow boundary.
	1060	DO k = nzb+1, nzt+1
	1061	DO i = nxl, nxr
	1062
	1063	denom = u_m_s(k,0,i) - u_m_s(k,1,i)
	1064
	1065	IF ( denom /= 0.0_wp ) THEN
	1066	c_u(k,i) = -c_max * ( u(k,0,i) - u_m_s(k,0,i) ) / ( denom * tsc(2) )
	1067	IF ( c_u(k,i) < 0.0_wp ) THEN
	1068	c_u(k,i) = 0.0_wp
	1069	ELSEIF ( c_u(k,i) > c_max ) THEN
	1070	c_u(k,i) = c_max
	1071	ENDIF
	1072	ELSE
	1073	c_u(k,i) = c_max
	1074	ENDIF
	1075
	1076	denom = v_m_s(k,1,i) - v_m_s(k,2,i)
	1077
	1078	IF ( denom /= 0.0_wp ) THEN
	1079	c_v(k,i) = -c_max * ( v(k,1,i) - v_m_s(k,1,i) ) / ( denom * tsc(2) )
	1080	IF ( c_v(k,i) < 0.0_wp ) THEN
	1081	c_v(k,i) = 0.0_wp
	1082	ELSEIF ( c_v(k,i) > c_max ) THEN
	1083	c_v(k,i) = c_max
	1084	ENDIF
	1085	ELSE
	1086	c_v(k,i) = c_max
	1087	ENDIF
	1088
	1089	denom = w_m_s(k,0,i) - w_m_s(k,1,i)
	1090
	1091	IF ( denom /= 0.0_wp ) THEN
	1092	c_w(k,i) = -c_max * ( w(k,0,i) - w_m_s(k,0,i) ) / ( denom * tsc(2) )
	1093	IF ( c_w(k,i) < 0.0_wp ) THEN
	1094	c_w(k,i) = 0.0_wp
	1095	ELSEIF ( c_w(k,i) > c_max ) THEN
	1096	c_w(k,i) = c_max
	1097	ENDIF
	1098	ELSE
	1099	c_w(k,i) = c_max
	1100	ENDIF
	1101
	1102	c_u_m_l(k) = c_u_m_l(k) + c_u(k,i)
	1103	c_v_m_l(k) = c_v_m_l(k) + c_v(k,i)
	1104	c_w_m_l(k) = c_w_m_l(k) + c_w(k,i)
	1105
	1106	ENDDO
	1107	ENDDO
	1108
	1109	#if defined( __parallel )
	1110	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1111	CALL MPI_ALLREDUCE( c_u_m_l(nzb+1), c_u_m(nzb+1), nzt-nzb, MPI_REAL, &
	1112	MPI_SUM, comm1dx, ierr )
	1113	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1114	CALL MPI_ALLREDUCE( c_v_m_l(nzb+1), c_v_m(nzb+1), nzt-nzb, MPI_REAL, &
	1115	MPI_SUM, comm1dx, ierr )
	1116	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1117	CALL MPI_ALLREDUCE( c_w_m_l(nzb+1), c_w_m(nzb+1), nzt-nzb, MPI_REAL, &
	1118	MPI_SUM, comm1dx, ierr )
	1119	#else
	1120	c_u_m = c_u_m_l
	1121	c_v_m = c_v_m_l
	1122	c_w_m = c_w_m_l
	1123	#endif
	1124
	1125	c_u_m = c_u_m / (nx+1)
	1126	c_v_m = c_v_m / (nx+1)
	1127	c_w_m = c_w_m / (nx+1)
	1128
	1129	!
	1130	!-- Save old timelevels for the next timestep
	1131	IF ( intermediate_timestep_count == 1 ) THEN
	1132	u_m_s(:,:,:) = u(:,0:1,:)
	1133	v_m_s(:,:,:) = v(:,1:2,:)
	1134	w_m_s(:,:,:) = w(:,0:1,:)
	1135	ENDIF
	1136
	1137	!
	1138	!-- Calculate the new velocities
	1139	DO k = nzb+1, nzt+1
	1140	DO i = nxlg, nxrg
	1141	u_p(k,-1,i) = u(k,-1,i) - dt_3d * tsc(2) * c_u_m(k) * &
	1142	( u(k,-1,i) - u(k,0,i) ) * ddy
	1143
	1144	v_p(k,0,i) = v(k,0,i) - dt_3d * tsc(2) * c_v_m(k) * &
	1145	( v(k,0,i) - v(k,1,i) ) * ddy
	1146
	1147	w_p(k,-1,i) = w(k,-1,i) - dt_3d * tsc(2) * c_w_m(k) * &
	1148	( w(k,-1,i) - w(k,0,i) ) * ddy
	1149	ENDDO
	1150	ENDDO
	1151
	1152	!
	1153	!-- Bottom boundary at the outflow
	1154	IF ( ibc_uv_b == 0 ) THEN
	1155	u_p(nzb,-1,:) = 0.0_wp
	1156	v_p(nzb,0,:) = 0.0_wp
	1157	ELSE
	1158	u_p(nzb,-1,:) = u_p(nzb+1,-1,:)
	1159	v_p(nzb,0,:) = v_p(nzb+1,0,:)
	1160	ENDIF
	1161	w_p(nzb,-1,:) = 0.0_wp
	1162
	1163	!
	1164	!-- Top boundary at the outflow
	1165	IF ( ibc_uv_t == 0 ) THEN
	1166	u_p(nzt+1,-1,:) = u_init(nzt+1)
	1167	v_p(nzt+1,0,:) = v_init(nzt+1)
	1168	ELSE
	1169	u_p(nzt+1,-1,:) = u_p(nzt,-1,:)
	1170	v_p(nzt+1,0,:) = v_p(nzt,0,:)
	1171	ENDIF
	1172	w_p(nzt:nzt+1,-1,:) = 0.0_wp
	1173
	1174	ENDIF
	1175
	1176	ENDIF
	1177
	1178	IF ( bc_radiation_n ) THEN
	1179
	1180	IF ( use_cmax ) THEN
	1181	u_p(:,ny+1,:) = u(:,ny,:)
	1182	v_p(:,ny+1,:) = v(:,ny,:)
	1183	w_p(:,ny+1,:) = w(:,ny,:)
	1184	ELSEIF ( .NOT. use_cmax ) THEN
	1185
	1186	c_max = dy / dt_3d
	1187
	1188	c_u_m_l = 0.0_wp
	1189	c_v_m_l = 0.0_wp
	1190	c_w_m_l = 0.0_wp
	1191
	1192	c_u_m = 0.0_wp
	1193	c_v_m = 0.0_wp
	1194	c_w_m = 0.0_wp
	1195
	1196	!
	1197	!-- Calculate the phase speeds for u, v, and w, first local and then
	1198	!-- average along the outflow boundary.
	1199	DO k = nzb+1, nzt+1
	1200	DO i = nxl, nxr
	1201
	1202	denom = u_m_n(k,ny,i) - u_m_n(k,ny-1,i)
	1203
	1204	IF ( denom /= 0.0_wp ) THEN
	1205	c_u(k,i) = -c_max * ( u(k,ny,i) - u_m_n(k,ny,i) ) / ( denom * tsc(2) )
	1206	IF ( c_u(k,i) < 0.0_wp ) THEN
	1207	c_u(k,i) = 0.0_wp
	1208	ELSEIF ( c_u(k,i) > c_max ) THEN
	1209	c_u(k,i) = c_max
	1210	ENDIF
	1211	ELSE
	1212	c_u(k,i) = c_max
	1213	ENDIF
	1214
	1215	denom = v_m_n(k,ny,i) - v_m_n(k,ny-1,i)
	1216
	1217	IF ( denom /= 0.0_wp ) THEN
	1218	c_v(k,i) = -c_max * ( v(k,ny,i) - v_m_n(k,ny,i) ) / ( denom * tsc(2) )
	1219	IF ( c_v(k,i) < 0.0_wp ) THEN
	1220	c_v(k,i) = 0.0_wp
	1221	ELSEIF ( c_v(k,i) > c_max ) THEN
	1222	c_v(k,i) = c_max
	1223	ENDIF
	1224	ELSE
	1225	c_v(k,i) = c_max
	1226	ENDIF
	1227
	1228	denom = w_m_n(k,ny,i) - w_m_n(k,ny-1,i)
	1229
	1230	IF ( denom /= 0.0_wp ) THEN
	1231	c_w(k,i) = -c_max * ( w(k,ny,i) - w_m_n(k,ny,i) ) / ( denom * tsc(2) )
	1232	IF ( c_w(k,i) < 0.0_wp ) THEN
	1233	c_w(k,i) = 0.0_wp
	1234	ELSEIF ( c_w(k,i) > c_max ) THEN
	1235	c_w(k,i) = c_max
	1236	ENDIF
	1237	ELSE
	1238	c_w(k,i) = c_max
	1239	ENDIF
	1240
	1241	c_u_m_l(k) = c_u_m_l(k) + c_u(k,i)
	1242	c_v_m_l(k) = c_v_m_l(k) + c_v(k,i)
	1243	c_w_m_l(k) = c_w_m_l(k) + c_w(k,i)
	1244
	1245	ENDDO
	1246	ENDDO
	1247
	1248	#if defined( __parallel )
	1249	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1250	CALL MPI_ALLREDUCE( c_u_m_l(nzb+1), c_u_m(nzb+1), nzt-nzb, MPI_REAL, &
	1251	MPI_SUM, comm1dx, ierr )
	1252	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1253	CALL MPI_ALLREDUCE( c_v_m_l(nzb+1), c_v_m(nzb+1), nzt-nzb, MPI_REAL, &
	1254	MPI_SUM, comm1dx, ierr )
	1255	IF ( collective_wait ) CALL MPI_BARRIER( comm1dx, ierr )
	1256	CALL MPI_ALLREDUCE( c_w_m_l(nzb+1), c_w_m(nzb+1), nzt-nzb, MPI_REAL, &
	1257	MPI_SUM, comm1dx, ierr )
	1258	#else
	1259	c_u_m = c_u_m_l
	1260	c_v_m = c_v_m_l
	1261	c_w_m = c_w_m_l
	1262	#endif
	1263
	1264	c_u_m = c_u_m / (nx+1)
	1265	c_v_m = c_v_m / (nx+1)
	1266	c_w_m = c_w_m / (nx+1)
	1267
	1268	!
	1269	!-- Save old timelevels for the next timestep
	1270	IF ( intermediate_timestep_count == 1 ) THEN
	1271	u_m_n(:,:,:) = u(:,ny-1:ny,:)
	1272	v_m_n(:,:,:) = v(:,ny-1:ny,:)
	1273	w_m_n(:,:,:) = w(:,ny-1:ny,:)
	1274	ENDIF
	1275
	1276	!
	1277	!-- Calculate the new velocities
	1278	DO k = nzb+1, nzt+1
	1279	DO i = nxlg, nxrg
	1280	u_p(k,ny+1,i) = u(k,ny+1,i) - dt_3d * tsc(2) * c_u_m(k) * &
	1281	( u(k,ny+1,i) - u(k,ny,i) ) * ddy
	1282
	1283	v_p(k,ny+1,i) = v(k,ny+1,i) - dt_3d * tsc(2) * c_v_m(k) * &
	1284	( v(k,ny+1,i) - v(k,ny,i) ) * ddy
	1285
	1286	w_p(k,ny+1,i) = w(k,ny+1,i) - dt_3d * tsc(2) * c_w_m(k) * &
	1287	( w(k,ny+1,i) - w(k,ny,i) ) * ddy
	1288	ENDDO
	1289	ENDDO
	1290
	1291	!
	1292	!-- Bottom boundary at the outflow
	1293	IF ( ibc_uv_b == 0 ) THEN
	1294	u_p(nzb,ny+1,:) = 0.0_wp
	1295	v_p(nzb,ny+1,:) = 0.0_wp
	1296	ELSE
	1297	u_p(nzb,ny+1,:) = u_p(nzb+1,ny+1,:)
	1298	v_p(nzb,ny+1,:) = v_p(nzb+1,ny+1,:)
	1299	ENDIF
	1300	w_p(nzb,ny+1,:) = 0.0_wp
	1301
	1302	!
	1303	!-- Top boundary at the outflow
	1304	IF ( ibc_uv_t == 0 ) THEN
	1305	u_p(nzt+1,ny+1,:) = u_init(nzt+1)
	1306	v_p(nzt+1,ny+1,:) = v_init(nzt+1)
	1307	ELSE
	1308	u_p(nzt+1,ny+1,:) = u_p(nzt,nyn+1,:)
	1309	v_p(nzt+1,ny+1,:) = v_p(nzt,nyn+1,:)
	1310	ENDIF
	1311	w_p(nzt:nzt+1,ny+1,:) = 0.0_wp
	1312
	1313	ENDIF
	1314
	1315	ENDIF
	1316
	1317	IF ( bc_radiation_l ) THEN
	1318
	1319	IF ( use_cmax ) THEN
	1320	u_p(:,:,0) = u(:,:,1)
	1321	v_p(:,:,-1) = v(:,:,0)
	1322	w_p(:,:,-1) = w(:,:,0)
	1323	ELSEIF ( .NOT. use_cmax ) THEN
	1324
	1325	c_max = dx / dt_3d
	1326
	1327	c_u_m_l = 0.0_wp
	1328	c_v_m_l = 0.0_wp
	1329	c_w_m_l = 0.0_wp
	1330
	1331	c_u_m = 0.0_wp
	1332	c_v_m = 0.0_wp
	1333	c_w_m = 0.0_wp
	1334
	1335	!
	1336	!-- Calculate the phase speeds for u, v, and w, first local and then
	1337	!-- average along the outflow boundary.
	1338	DO k = nzb+1, nzt+1
	1339	DO j = nys, nyn
	1340
	1341	denom = u_m_l(k,j,1) - u_m_l(k,j,2)
	1342
	1343	IF ( denom /= 0.0_wp ) THEN
	1344	c_u(k,j) = -c_max * ( u(k,j,1) - u_m_l(k,j,1) ) / ( denom * tsc(2) )
	1345	IF ( c_u(k,j) < 0.0_wp ) THEN
	1346	c_u(k,j) = 0.0_wp
	1347	ELSEIF ( c_u(k,j) > c_max ) THEN
	1348	c_u(k,j) = c_max
	1349	ENDIF
	1350	ELSE
	1351	c_u(k,j) = c_max
	1352	ENDIF
	1353
	1354	denom = v_m_l(k,j,0) - v_m_l(k,j,1)
	1355
	1356	IF ( denom /= 0.0_wp ) THEN
	1357	c_v(k,j) = -c_max * ( v(k,j,0) - v_m_l(k,j,0) ) / ( denom * tsc(2) )
	1358	IF ( c_v(k,j) < 0.0_wp ) THEN
	1359	c_v(k,j) = 0.0_wp
	1360	ELSEIF ( c_v(k,j) > c_max ) THEN
	1361	c_v(k,j) = c_max
	1362	ENDIF
	1363	ELSE
	1364	c_v(k,j) = c_max
	1365	ENDIF
	1366
	1367	denom = w_m_l(k,j,0) - w_m_l(k,j,1)
	1368
	1369	IF ( denom /= 0.0_wp ) THEN
	1370	c_w(k,j) = -c_max * ( w(k,j,0) - w_m_l(k,j,0) ) / ( denom * tsc(2) )
	1371	IF ( c_w(k,j) < 0.0_wp ) THEN
	1372	c_w(k,j) = 0.0_wp
	1373	ELSEIF ( c_w(k,j) > c_max ) THEN
	1374	c_w(k,j) = c_max
	1375	ENDIF
	1376	ELSE
	1377	c_w(k,j) = c_max
	1378	ENDIF
	1379
	1380	c_u_m_l(k) = c_u_m_l(k) + c_u(k,j)
	1381	c_v_m_l(k) = c_v_m_l(k) + c_v(k,j)
	1382	c_w_m_l(k) = c_w_m_l(k) + c_w(k,j)
	1383
	1384	ENDDO
	1385	ENDDO
	1386
	1387	#if defined( __parallel )
	1388	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1389	CALL MPI_ALLREDUCE( c_u_m_l(nzb+1), c_u_m(nzb+1), nzt-nzb, MPI_REAL, &
	1390	MPI_SUM, comm1dy, ierr )
	1391	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1392	CALL MPI_ALLREDUCE( c_v_m_l(nzb+1), c_v_m(nzb+1), nzt-nzb, MPI_REAL, &
	1393	MPI_SUM, comm1dy, ierr )
	1394	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1395	CALL MPI_ALLREDUCE( c_w_m_l(nzb+1), c_w_m(nzb+1), nzt-nzb, MPI_REAL, &
	1396	MPI_SUM, comm1dy, ierr )
	1397	#else
	1398	c_u_m = c_u_m_l
	1399	c_v_m = c_v_m_l
	1400	c_w_m = c_w_m_l
	1401	#endif
	1402
	1403	c_u_m = c_u_m / (ny+1)
	1404	c_v_m = c_v_m / (ny+1)
	1405	c_w_m = c_w_m / (ny+1)
	1406
	1407	!
	1408	!-- Save old timelevels for the next timestep
	1409	IF ( intermediate_timestep_count == 1 ) THEN
	1410	u_m_l(:,:,:) = u(:,:,1:2)
	1411	v_m_l(:,:,:) = v(:,:,0:1)
	1412	w_m_l(:,:,:) = w(:,:,0:1)
	1413	ENDIF
	1414
	1415	!
	1416	!-- Calculate the new velocities
	1417	DO k = nzb+1, nzt+1
	1418	DO j = nysg, nyng
	1419	u_p(k,j,0) = u(k,j,0) - dt_3d * tsc(2) * c_u_m(k) * &
	1420	( u(k,j,0) - u(k,j,1) ) * ddx
	1421
	1422	v_p(k,j,-1) = v(k,j,-1) - dt_3d * tsc(2) * c_v_m(k) * &
	1423	( v(k,j,-1) - v(k,j,0) ) * ddx
	1424
	1425	w_p(k,j,-1) = w(k,j,-1) - dt_3d * tsc(2) * c_w_m(k) * &
	1426	( w(k,j,-1) - w(k,j,0) ) * ddx
	1427	ENDDO
	1428	ENDDO
	1429
	1430	!
	1431	!-- Bottom boundary at the outflow
	1432	IF ( ibc_uv_b == 0 ) THEN
	1433	u_p(nzb,:,0) = 0.0_wp
	1434	v_p(nzb,:,-1) = 0.0_wp
	1435	ELSE
	1436	u_p(nzb,:,0) = u_p(nzb+1,:,0)
	1437	v_p(nzb,:,-1) = v_p(nzb+1,:,-1)
	1438	ENDIF
	1439	w_p(nzb,:,-1) = 0.0_wp
	1440
	1441	!
	1442	!-- Top boundary at the outflow
	1443	IF ( ibc_uv_t == 0 ) THEN
	1444	u_p(nzt+1,:,0) = u_init(nzt+1)
	1445	v_p(nzt+1,:,-1) = v_init(nzt+1)
	1446	ELSE
	1447	u_p(nzt+1,:,0) = u_p(nzt,:,0)
	1448	v_p(nzt+1,:,-1) = v_p(nzt,:,-1)
	1449	ENDIF
	1450	w_p(nzt:nzt+1,:,-1) = 0.0_wp
	1451
	1452	ENDIF
	1453
	1454	ENDIF
	1455
	1456	IF ( bc_radiation_r ) THEN
	1457
	1458	IF ( use_cmax ) THEN
	1459	u_p(:,:,nx+1) = u(:,:,nx)
	1460	v_p(:,:,nx+1) = v(:,:,nx)
	1461	w_p(:,:,nx+1) = w(:,:,nx)
	1462	ELSEIF ( .NOT. use_cmax ) THEN
	1463
	1464	c_max = dx / dt_3d
	1465
	1466	c_u_m_l = 0.0_wp
	1467	c_v_m_l = 0.0_wp
	1468	c_w_m_l = 0.0_wp
	1469
	1470	c_u_m = 0.0_wp
	1471	c_v_m = 0.0_wp
	1472	c_w_m = 0.0_wp
	1473
	1474	!
	1475	!-- Calculate the phase speeds for u, v, and w, first local and then
	1476	!-- average along the outflow boundary.
	1477	DO k = nzb+1, nzt+1
	1478	DO j = nys, nyn
	1479
	1480	denom = u_m_r(k,j,nx) - u_m_r(k,j,nx-1)
	1481
	1482	IF ( denom /= 0.0_wp ) THEN
	1483	c_u(k,j) = -c_max * ( u(k,j,nx) - u_m_r(k,j,nx) ) / ( denom * tsc(2) )
	1484	IF ( c_u(k,j) < 0.0_wp ) THEN
	1485	c_u(k,j) = 0.0_wp
	1486	ELSEIF ( c_u(k,j) > c_max ) THEN
	1487	c_u(k,j) = c_max
	1488	ENDIF
	1489	ELSE
	1490	c_u(k,j) = c_max
	1491	ENDIF
	1492
	1493	denom = v_m_r(k,j,nx) - v_m_r(k,j,nx-1)
	1494
	1495	IF ( denom /= 0.0_wp ) THEN
	1496	c_v(k,j) = -c_max * ( v(k,j,nx) - v_m_r(k,j,nx) ) / ( denom * tsc(2) )
	1497	IF ( c_v(k,j) < 0.0_wp ) THEN
	1498	c_v(k,j) = 0.0_wp
	1499	ELSEIF ( c_v(k,j) > c_max ) THEN
	1500	c_v(k,j) = c_max
	1501	ENDIF
	1502	ELSE
	1503	c_v(k,j) = c_max
	1504	ENDIF
	1505
	1506	denom = w_m_r(k,j,nx) - w_m_r(k,j,nx-1)
	1507
	1508	IF ( denom /= 0.0_wp ) THEN
	1509	c_w(k,j) = -c_max * ( w(k,j,nx) - w_m_r(k,j,nx) ) / ( denom * tsc(2) )
	1510	IF ( c_w(k,j) < 0.0_wp ) THEN
	1511	c_w(k,j) = 0.0_wp
	1512	ELSEIF ( c_w(k,j) > c_max ) THEN
	1513	c_w(k,j) = c_max
	1514	ENDIF
	1515	ELSE
	1516	c_w(k,j) = c_max
	1517	ENDIF
	1518
	1519	c_u_m_l(k) = c_u_m_l(k) + c_u(k,j)
	1520	c_v_m_l(k) = c_v_m_l(k) + c_v(k,j)
	1521	c_w_m_l(k) = c_w_m_l(k) + c_w(k,j)
	1522
	1523	ENDDO
	1524	ENDDO
	1525
	1526	#if defined( __parallel )
	1527	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1528	CALL MPI_ALLREDUCE( c_u_m_l(nzb+1), c_u_m(nzb+1), nzt-nzb, MPI_REAL, &
	1529	MPI_SUM, comm1dy, ierr )
	1530	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1531	CALL MPI_ALLREDUCE( c_v_m_l(nzb+1), c_v_m(nzb+1), nzt-nzb, MPI_REAL, &
	1532	MPI_SUM, comm1dy, ierr )
	1533	IF ( collective_wait ) CALL MPI_BARRIER( comm1dy, ierr )
	1534	CALL MPI_ALLREDUCE( c_w_m_l(nzb+1), c_w_m(nzb+1), nzt-nzb, MPI_REAL, &
	1535	MPI_SUM, comm1dy, ierr )
	1536	#else
	1537	c_u_m = c_u_m_l
	1538	c_v_m = c_v_m_l
	1539	c_w_m = c_w_m_l
	1540	#endif
	1541
	1542	c_u_m = c_u_m / (ny+1)
	1543	c_v_m = c_v_m / (ny+1)
	1544	c_w_m = c_w_m / (ny+1)
	1545
	1546	!
	1547	!-- Save old timelevels for the next timestep
	1548	IF ( intermediate_timestep_count == 1 ) THEN
	1549	u_m_r(:,:,:) = u(:,:,nx-1:nx)
	1550	v_m_r(:,:,:) = v(:,:,nx-1:nx)
	1551	w_m_r(:,:,:) = w(:,:,nx-1:nx)
	1552	ENDIF
	1553
	1554	!
	1555	!-- Calculate the new velocities
	1556	DO k = nzb+1, nzt+1
	1557	DO j = nysg, nyng
	1558	u_p(k,j,nx+1) = u(k,j,nx+1) - dt_3d * tsc(2) * c_u_m(k) * &
	1559	( u(k,j,nx+1) - u(k,j,nx) ) * ddx
	1560
	1561	v_p(k,j,nx+1) = v(k,j,nx+1) - dt_3d * tsc(2) * c_v_m(k) * &
	1562	( v(k,j,nx+1) - v(k,j,nx) ) * ddx
	1563
	1564	w_p(k,j,nx+1) = w(k,j,nx+1) - dt_3d * tsc(2) * c_w_m(k) * &
	1565	( w(k,j,nx+1) - w(k,j,nx) ) * ddx
	1566	ENDDO
	1567	ENDDO
	1568
	1569	!
	1570	!-- Bottom boundary at the outflow
	1571	IF ( ibc_uv_b == 0 ) THEN
	1572	u_p(nzb,:,nx+1) = 0.0_wp
	1573	v_p(nzb,:,nx+1) = 0.0_wp
	1574	ELSE
	1575	u_p(nzb,:,nx+1) = u_p(nzb+1,:,nx+1)
	1576	v_p(nzb,:,nx+1) = v_p(nzb+1,:,nx+1)
	1577	ENDIF
	1578	w_p(nzb,:,nx+1) = 0.0_wp
	1579
	1580	!
	1581	!-- Top boundary at the outflow
	1582	IF ( ibc_uv_t == 0 ) THEN
	1583	u_p(nzt+1,:,nx+1) = u_init(nzt+1)
	1584	v_p(nzt+1,:,nx+1) = v_init(nzt+1)
	1585	ELSE
	1586	u_p(nzt+1,:,nx+1) = u_p(nzt,:,nx+1)
	1587	v_p(nzt+1,:,nx+1) = v_p(nzt,:,nx+1)
	1588	ENDIF
	1589	w_p(nzt:nzt+1,:,nx+1) = 0.0_wp
	1590
	1591	ENDIF
	1592
	1593	ENDIF
	1594
	1595	END SUBROUTINE dynamics_boundary_conditions
[4047]	1596	!------------------------------------------------------------------------------!
	1597	! Description:
	1598	! ------------
	1599	!> Swap timelevels of module-specific array pointers
	1600	!------------------------------------------------------------------------------!
	1601	SUBROUTINE dynamics_swap_timelevel ( mod_count )
	1602
	1603
	1604	INTEGER, INTENT(IN) :: mod_count
	1605
	1606
	1607	SELECT CASE ( mod_count )
	1608
	1609	CASE ( 0 )
	1610
	1611	u => u_1; u_p => u_2
	1612	v => v_1; v_p => v_2
	1613	w => w_1; w_p => w_2
	1614	IF ( .NOT. neutral ) THEN
	1615	pt => pt_1; pt_p => pt_2
	1616	ENDIF
	1617	IF ( humidity ) THEN
	1618	q => q_1; q_p => q_2
	1619	ENDIF
	1620	IF ( passive_scalar ) THEN
	1621	s => s_1; s_p => s_2
	1622	ENDIF
	1623
	1624	CASE ( 1 )
	1625
	1626	u => u_2; u_p => u_1
	1627	v => v_2; v_p => v_1
	1628	w => w_2; w_p => w_1
	1629	IF ( .NOT. neutral ) THEN
	1630	pt => pt_2; pt_p => pt_1
	1631	ENDIF
	1632	IF ( humidity ) THEN
	1633	q => q_2; q_p => q_1
	1634	ENDIF
	1635	IF ( passive_scalar ) THEN
	1636	s => s_2; s_p => s_1
	1637	ENDIF
	1638
	1639	END SELECT
	1640
	1641	END SUBROUTINE dynamics_swap_timelevel
	1642
	1643
	1644	!--------------------------------------------------------------------------------------------------!
	1645	! Description:
	1646	! ------------
	1647	!> Sum up and time-average module-specific output quantities
	1648	!> as well as allocate the array necessary for storing the average.
	1649	!--------------------------------------------------------------------------------------------------!
	1650	SUBROUTINE dynamics_3d_data_averaging( mode, variable )
	1651
	1652
	1653	CHARACTER (LEN=*) :: mode !<
	1654	CHARACTER (LEN=*) :: variable !<
	1655
	1656
	1657	IF ( mode == 'allocate' ) THEN
	1658
	1659	SELECT CASE ( TRIM( variable ) )
	1660
	1661	! CASE ( 'u2' )
	1662
	1663	CASE DEFAULT
	1664	CONTINUE
	1665
	1666	END SELECT
	1667
	1668	ELSEIF ( mode == 'sum' ) THEN
	1669
	1670	SELECT CASE ( TRIM( variable ) )
	1671
	1672	! CASE ( 'u2' )
	1673
	1674	CASE DEFAULT
	1675	CONTINUE
	1676
	1677	END SELECT
	1678
	1679	ELSEIF ( mode == 'average' ) THEN
	1680
	1681	SELECT CASE ( TRIM( variable ) )
	1682
	1683	! CASE ( 'u2' )
	1684
	1685	END SELECT
	1686
	1687	ENDIF
	1688
	1689
	1690	END SUBROUTINE dynamics_3d_data_averaging
	1691
	1692
	1693	!--------------------------------------------------------------------------------------------------!
	1694	! Description:
	1695	! ------------
	1696	!> Resorts the module-specific output quantity with indices (k,j,i) to a
	1697	!> temporary array with indices (i,j,k) and sets the grid on which it is defined.
	1698	!> Allowed values for grid are "zu" and "zw".
	1699	!--------------------------------------------------------------------------------------------------!
	1700	SUBROUTINE dynamics_data_output_2d( av, variable, found, grid, mode, local_pf, &
	1701	two_d, nzb_do, nzt_do, fill_value )
	1702
	1703
	1704	CHARACTER (LEN=*) :: grid !<
	1705	CHARACTER (LEN=*), INTENT(IN) :: mode !< either 'xy', 'xz' or 'yz'
	1706	CHARACTER (LEN=*) :: variable !<
	1707
	1708	INTEGER(iwp) :: av !< flag to control data output of instantaneous or time-averaged data
	1709	! INTEGER(iwp) :: i !< grid index along x-direction
	1710	! INTEGER(iwp) :: j !< grid index along y-direction
	1711	! INTEGER(iwp) :: k !< grid index along z-direction
	1712	! INTEGER(iwp) :: m !< running index surface elements
	1713	INTEGER(iwp) :: nzb_do !< lower limit of the domain (usually nzb)
	1714	INTEGER(iwp) :: nzt_do !< upper limit of the domain (usually nzt+1)
	1715
	1716	LOGICAL :: found !<
	1717	LOGICAL :: two_d !< flag parameter that indicates 2D variables (horizontal cross sections)
	1718
	1719	REAL(wp), INTENT(IN) :: fill_value
	1720
	1721	REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
	1722
	1723	!
	1724	!-- Next line is just to avoid compiler warnings about unused variables. You may remove it.
	1725	IF ( two_d .AND. av + LEN( mode ) + local_pf(nxl,nys,nzb_do) + fill_value < 0.0 ) CONTINUE
	1726
	1727	found = .TRUE.
	1728
	1729	SELECT CASE ( TRIM( variable ) )
	1730
	1731	! CASE ( 'u2_xy', 'u2_xz', 'u2_yz' )
	1732
	1733	CASE DEFAULT
	1734	found = .FALSE.
	1735	grid = 'none'
	1736
	1737	END SELECT
	1738
	1739
	1740	END SUBROUTINE dynamics_data_output_2d
	1741
	1742
	1743	!--------------------------------------------------------------------------------------------------!
	1744	! Description:
	1745	! ------------
	1746	!> Resorts the module-specific output quantity with indices (k,j,i)
	1747	!> to a temporary array with indices (i,j,k).
	1748	!--------------------------------------------------------------------------------------------------!
	1749	SUBROUTINE dynamics_data_output_3d( av, variable, found, local_pf, fill_value, nzb_do, nzt_do )
	1750
	1751
	1752	CHARACTER (LEN=*) :: variable !<
	1753
	1754	INTEGER(iwp) :: av !<
	1755	! INTEGER(iwp) :: i !<
	1756	! INTEGER(iwp) :: j !<
	1757	! INTEGER(iwp) :: k !<
	1758	INTEGER(iwp) :: nzb_do !< lower limit of the data output (usually 0)
	1759	INTEGER(iwp) :: nzt_do !< vertical upper limit of the data output (usually nz_do3d)
	1760
	1761	LOGICAL :: found !<
	1762
	1763	REAL(wp), INTENT(IN) :: fill_value !< value for the _FillValue attribute
	1764
	1765	REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !<
	1766
	1767	!
	1768	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	1769	IF ( av + local_pf(nxl,nys,nzb_do) + fill_value < 0.0 ) CONTINUE
	1770
	1771
	1772	found = .TRUE.
	1773
	1774	SELECT CASE ( TRIM( variable ) )
	1775
	1776	! CASE ( 'u2' )
	1777
	1778	CASE DEFAULT
	1779	found = .FALSE.
	1780
	1781	END SELECT
	1782
	1783
	1784	END SUBROUTINE dynamics_data_output_3d
	1785
	1786
	1787	!--------------------------------------------------------------------------------------------------!
	1788	! Description:
	1789	! ------------
	1790	!> Calculation of module-specific statistics, i.e. horizontally averaged profiles and time series.
	1791	!> This is called for every statistic region sr, but at least for the region "total domain" (sr=0).
	1792	!--------------------------------------------------------------------------------------------------!
	1793	SUBROUTINE dynamics_statistics( mode, sr, tn )
	1794
	1795
	1796	CHARACTER (LEN=*) :: mode !<
	1797	! INTEGER(iwp) :: i !<
	1798	! INTEGER(iwp) :: j !<
	1799	! INTEGER(iwp) :: k !<
	1800	INTEGER(iwp) :: sr !<
	1801	INTEGER(iwp) :: tn !<
	1802
	1803	!
	1804	!-- Next line is to avoid compiler warning about unused variables. Please remove.
	1805	IF ( sr == 0 .OR. tn == 0 ) CONTINUE
	1806
	1807	IF ( mode == 'profiles' ) THEN
	1808
	1809	ELSEIF ( mode == 'time_series' ) THEN
	1810
	1811	ENDIF
	1812
	1813	END SUBROUTINE dynamics_statistics
	1814
	1815
	1816	!--------------------------------------------------------------------------------------------------!
	1817	! Description:
	1818	! ------------
	1819	!> Read module-specific global restart data.
	1820	!--------------------------------------------------------------------------------------------------!
	1821	SUBROUTINE dynamics_rrd_global( found )
	1822
	1823
	1824	LOGICAL, INTENT(OUT) :: found
	1825
	1826
	1827	found = .TRUE.
	1828
	1829
	1830	SELECT CASE ( restart_string(1:length) )
	1831
	1832	CASE ( 'global_paramter' )
	1833	! READ ( 13 ) global_parameter
	1834
	1835	CASE DEFAULT
	1836
	1837	found = .FALSE.
	1838
	1839	END SELECT
	1840
	1841
	1842	END SUBROUTINE dynamics_rrd_global
	1843
	1844
	1845	!--------------------------------------------------------------------------------------------------!
	1846	! Description:
	1847	! ------------
	1848	!> Read module-specific processor specific restart data from file(s).
	1849	!> Subdomain index limits on file are given by nxl_on_file, etc.
	1850	!> Indices nxlc, etc. indicate the range of gridpoints to be mapped from the subdomain on file (f)
	1851	!> to the subdomain of the current PE (c). They have been calculated in routine rrd_local.
	1852	!--------------------------------------------------------------------------------------------------!
	1853	SUBROUTINE dynamics_rrd_local( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, nxr_on_file, nynf, nync, &
	1854	nyn_on_file, nysf, nysc, nys_on_file, tmp_2d, tmp_3d, found )
	1855
	1856
	1857	INTEGER(iwp) :: k !<
	1858	INTEGER(iwp) :: nxlc !<
	1859	INTEGER(iwp) :: nxlf !<
	1860	INTEGER(iwp) :: nxl_on_file !<
	1861	INTEGER(iwp) :: nxrc !<
	1862	INTEGER(iwp) :: nxrf !<
	1863	INTEGER(iwp) :: nxr_on_file !<
	1864	INTEGER(iwp) :: nync !<
	1865	INTEGER(iwp) :: nynf !<
	1866	INTEGER(iwp) :: nyn_on_file !<
	1867	INTEGER(iwp) :: nysc !<
	1868	INTEGER(iwp) :: nysf !<
	1869	INTEGER(iwp) :: nys_on_file !<
	1870
	1871	LOGICAL, INTENT(OUT) :: found
	1872
	1873	REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d !<
	1874	REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d !<
	1875
	1876	!
	1877	!-- Next line is to avoid compiler warning about unused variables. Please remove.
[4097]	1878	IF ( k + nxlc + nxlf + nxrc + nxrf + nync + nynf + nysc + nysf + &
	1879	tmp_2d(nys_on_file,nxl_on_file) + &
	1880	tmp_3d(nzb,nys_on_file,nxl_on_file) < 0.0 ) CONTINUE
[4047]	1881	!
	1882	!-- Here the reading of user-defined restart data follows:
	1883	!-- Sample for user-defined output
	1884
	1885	found = .TRUE.
	1886
	1887	SELECT CASE ( restart_string(1:length) )
	1888
	1889	! CASE ( 'u2_av' )
	1890
	1891	CASE DEFAULT
	1892
	1893	found = .FALSE.
	1894
	1895	END SELECT
	1896
	1897	END SUBROUTINE dynamics_rrd_local
	1898
	1899
	1900	!--------------------------------------------------------------------------------------------------!
	1901	! Description:
	1902	! ------------
	1903	!> Writes global module-specific restart data into binary file(s) for restart runs.
	1904	!--------------------------------------------------------------------------------------------------!
	1905	SUBROUTINE dynamics_wrd_global
	1906
	1907
	1908	END SUBROUTINE dynamics_wrd_global
	1909
	1910
	1911	!--------------------------------------------------------------------------------------------------!
	1912	! Description:
	1913	! ------------
	1914	!> Writes processor specific and module-specific restart data into binary file(s) for restart runs.
	1915	!--------------------------------------------------------------------------------------------------!
	1916	SUBROUTINE dynamics_wrd_local
	1917
	1918
	1919	END SUBROUTINE dynamics_wrd_local
	1920
	1921
	1922	!--------------------------------------------------------------------------------------------------!
	1923	! Description:
	1924	! ------------
	1925	!> Execute module-specific actions at the very end of the program.
	1926	!--------------------------------------------------------------------------------------------------!
	1927	SUBROUTINE dynamics_last_actions
	1928
	1929
	1930	END SUBROUTINE dynamics_last_actions
	1931
	1932	END MODULE dynamics_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |