Home

Context Navigation

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

Last change on this file since 4359 was 4359, checked in by suehring, 4 years ago
Fix for last commit
Property svn:keywords set to `Id`
File size: 59.4 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |