source: palm/trunk/SOURCE/time_integration.f90 @ 1976

Last change on this file since 1976 was 1976, checked in by maronga, 8 years ago

further modularization of land surface model (2D/3D output and restart data). Bugfix for restart runs without land surface model

  • Property svn:keywords set to Id
File size: 46.6 KB
Line 
1!> @file time_integration.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms
6! of the GNU General Public License as published by the Free Software Foundation,
7! either version 3 of the License, or (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with
14! PALM. If not, see <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2016 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------!
18!
19! Current revisions:
20! ------------------
21! Simplified calls to radiation model
22!
23! Former revisions:
24! -----------------
25! $Id: time_integration.f90 1976 2016-07-27 13:28:04Z maronga $
26!
27! 1960 2016-07-12 16:34:24Z suehring
28! Separate humidity and passive scalar
29!
30! 1957 2016-07-07 10:43:48Z suehring
31! flight module added
32!
33! 1919 2016-05-27 14:51:23Z raasch
34! Initial version of purely vertical nesting introduced.
35!
36! 1918 2016-05-27 14:35:57Z raasch
37! determination of time step moved to the end of the time step loop,
38! the first time step is now always calculated before the time step loop (i.e.
39! also in case of restart runs)
40!
41! 1914 2016-05-26 14:44:07Z witha
42! Added call for wind turbine model
43!
44! 1878 2016-04-19 12:30:36Z hellstea
45! Synchronization for nested runs rewritten
46!
47! 1853 2016-04-11 09:00:35Z maronga
48! Adjusted for use with radiation_scheme = constant
49!
50! 1849 2016-04-08 11:33:18Z hoffmann
51! Adapted for modularization of microphysics
52!
53! 1833 2016-04-07 14:23:03Z raasch
54! spectrum renamed spectra_mod, spectra related variables moved to spectra_mod
55!
56! 1831 2016-04-07 13:15:51Z hoffmann
57! turbulence renamed collision_turbulence
58!
59! 1822 2016-04-07 07:49:42Z hoffmann
60! icloud_scheme replaced by microphysics_*
61!
62! 1808 2016-04-05 19:44:00Z raasch
63! output message in case unscheduled radiation calls removed
64!
65! 1797 2016-03-21 16:50:28Z raasch
66! introduction of different datatransfer modes
67!
68! 1791 2016-03-11 10:41:25Z raasch
69! call of pmci_update_new removed
70!
71! 1786 2016-03-08 05:49:27Z raasch
72! +module spectrum
73!
74! 1783 2016-03-06 18:36:17Z raasch
75! switch back of netcdf data format for mask output moved to the mask output
76! routine
77!
78! 1781 2016-03-03 15:12:23Z raasch
79! some pmc calls removed at the beginning (before timeloop),
80! pmc initialization moved to the main program
81!
82! 1764 2016-02-28 12:45:19Z raasch
83! PMC_ACTIVE flags removed,
84! bugfix: nest synchronization after first call of timestep
85!
86! 1762 2016-02-25 12:31:13Z hellstea
87! Introduction of nested domain feature
88!
89! 1736 2015-12-04 08:56:33Z raasch
90! no perturbations added to total domain if energy limit has been set zero
91!
92! 1691 2015-10-26 16:17:44Z maronga
93! Added option for spin-ups without land surface and radiation models. Moved calls
94! for radiation and lan surface schemes.
95!
96! 1682 2015-10-07 23:56:08Z knoop
97! Code annotations made doxygen readable
98!
99! 1671 2015-09-25 03:29:37Z raasch
100! bugfix: ghostpoint exchange for array diss in case that sgs velocities are used
101! for particles
102!
103! 1585 2015-04-30 07:05:52Z maronga
104! Moved call of radiation scheme. Added support for RRTM
105!
106! 1551 2015-03-03 14:18:16Z maronga
107! Added interface for different radiation schemes.
108!
109! 1496 2014-12-02 17:25:50Z maronga
110! Added calls for the land surface model and radiation scheme
111!
112! 1402 2014-05-09 14:25:13Z raasch
113! location messages modified
114!
115! 1384 2014-05-02 14:31:06Z raasch
116! location messages added
117!
118! 1380 2014-04-28 12:40:45Z heinze
119! CALL of nudge_ref added
120! bc_pt_t_val and bc_q_t_val are updated in case nudging is used
121!
122! 1365 2014-04-22 15:03:56Z boeske
123! Reset sums_ls_l to zero at each timestep
124! +sums_ls_l
125! Calculation of reference state (previously in subroutine calc_mean_profile)
126
127! 1342 2014-03-26 17:04:47Z kanani
128! REAL constants defined as wp-kind
129!
130! 1320 2014-03-20 08:40:49Z raasch
131! ONLY-attribute added to USE-statements,
132! kind-parameters added to all INTEGER and REAL declaration statements,
133! kinds are defined in new module kinds,
134! old module precision_kind is removed,
135! revision history before 2012 removed,
136! comment fields (!:) to be used for variable explanations added to
137! all variable declaration statements
138! 1318 2014-03-17 13:35:16Z raasch
139! module interfaces removed
140!
141! 1308 2014-03-13 14:58:42Z fricke
142! +netcdf_data_format_save
143! For masked data, parallel netcdf output is not tested so far, hence
144! netcdf_data_format is switched back to non-paralell output.
145!
146! 1276 2014-01-15 13:40:41Z heinze
147! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
148!
149! 1257 2013-11-08 15:18:40Z raasch
150! acc-update-host directive for timestep removed
151!
152! 1241 2013-10-30 11:36:58Z heinze
153! Generalize calc_mean_profile for wider use
154! Determine shf and qsws in dependence on data from LSF_DATA
155! Determine ug and vg in dependence on data from LSF_DATA
156! 1221 2013-09-10 08:59:13Z raasch
157! host update of arrays before timestep is called
158!
159! 1179 2013-06-14 05:57:58Z raasch
160! mean profiles for reference state are only calculated if required,
161! small bugfix for background communication
162!
163! 1171 2013-05-30 11:27:45Z raasch
164! split of prognostic_equations deactivated (comment lines), for the time being
165!
166! 1128 2013-04-12 06:19:32Z raasch
167! asynchronous transfer of ghost point data realized for acc-optimized version:
168! prognostic_equations are first called two times for those points required for
169! the left-right and north-south exchange, respectively, and then for the
170! remaining points,
171! those parts requiring global communication moved from prognostic_equations to
172! here
173!
174! 1115 2013-03-26 18:16:16Z hoffmann
175! calculation of qr and nr is restricted to precipitation
176!
177! 1113 2013-03-10 02:48:14Z raasch
178! GPU-porting of boundary conditions,
179! openACC directives updated
180! formal parameter removed from routine boundary_conds
181!
182! 1111 2013-03-08 23:54:10Z raasch
183! +internal timestep counter for cpu statistics added,
184! openACC directives updated
185!
186! 1092 2013-02-02 11:24:22Z raasch
187! unused variables removed
188!
189! 1065 2012-11-22 17:42:36Z hoffmann
190! exchange of diss (dissipation rate) in case of turbulence = .TRUE. added
191!
192! 1053 2012-11-13 17:11:03Z hoffmann
193! exchange of ghost points for nr, qr added
194!
195! 1036 2012-10-22 13:43:42Z raasch
196! code put under GPL (PALM 3.9)
197!
198! 1019 2012-09-28 06:46:45Z raasch
199! non-optimized version of prognostic_equations removed
200!
201! 1015 2012-09-27 09:23:24Z raasch
202! +call of prognostic_equations_acc
203!
204! 1001 2012-09-13 14:08:46Z raasch
205! all actions concerning leapfrog- and upstream-spline-scheme removed
206!
207! 849 2012-03-15 10:35:09Z raasch
208! advec_particles renamed lpm, first_call_advec_particles renamed first_call_lpm
209!
210! 825 2012-02-19 03:03:44Z raasch
211! wang_collision_kernel renamed wang_kernel
212!
213! Revision 1.1  1997/08/11 06:19:04  raasch
214! Initial revision
215!
216!
217! Description:
218! ------------
219!> Integration in time of the model equations, statistical analysis and graphic
220!> output
221!------------------------------------------------------------------------------!
222 SUBROUTINE time_integration
223 
224
225    USE advec_ws,                                                              &
226        ONLY:  ws_statistics
227
228    USE arrays_3d,                                                             &
229        ONLY:  diss, dzu, e, e_p, nr_p, prho, pt, pt_p, pt_init, q_init, q,    &
230               ql, ql_c, ql_v, ql_vp, qr_p, q_p, ref_state, rho, s, s_p, sa_p, &
231               tend, u, u_p, v, vpt, v_p, w, w_p
232
233    USE calc_mean_profile_mod,                                                 &
234        ONLY:  calc_mean_profile
235
236    USE control_parameters,                                                    &
237        ONLY:  advected_distance_x, advected_distance_y, average_count_3d,     &
238               averaging_interval, averaging_interval_pr,                      &
239               bc_lr_cyc, bc_ns_cyc, bc_pt_t_val,                              &
240               bc_q_t_val, call_psolver_at_all_substeps, cloud_droplets,       &
241               cloud_physics, constant_flux_layer, constant_heatflux,          &
242               create_disturbances, dopr_n, constant_diffusion, coupling_mode, &
243               coupling_start_time, current_timestep_number,                   &
244               disturbance_created, disturbance_energy_limit, dist_range,      &
245               do_sum, dt_3d, dt_averaging_input, dt_averaging_input_pr,       &
246               dt_coupling, dt_data_output_av, dt_disturb, dt_do2d_xy,         &
247               dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_domask,dt_dopts, dt_dopr,   &
248               dt_dopr_listing, dt_dots, dt_dvrp, dt_run_control,              &
249               end_time, first_call_lpm, galilei_transformation, humidity,     &
250               intermediate_timestep_count,                                    &
251               intermediate_timestep_count_max, large_scale_forcing,           &
252               loop_optimization, lsf_surf, lsf_vert, masks,                   &
253               microphysics_seifert, mid, nest_domain,                         &
254               neutral, nr_timesteps_this_run, nudging,                        &
255               ocean, on_device, passive_scalar,                               &
256               prho_reference, pt_reference, pt_slope_offset, random_heatflux, &
257               run_coupled, simulated_time, simulated_time_chr,                &
258               skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz,        &
259               skip_time_do3d, skip_time_domask, skip_time_dopr,               &
260               skip_time_data_output_av, sloping_surface,                      &
261               stop_dt, terminate_coupled, terminate_run, timestep_scheme,     &
262               time_coupling, time_do2d_xy, time_do2d_xz, time_do2d_yz,        &
263               time_do3d, time_domask, time_dopr, time_dopr_av,                &
264               time_dopr_listing, time_dopts, time_dosp, time_dosp_av,         &
265               time_dots, time_do_av, time_do_sla, time_disturb, time_dvrp,    &
266               time_run_control, time_since_reference_point,                   &
267               turbulent_inflow, use_initial_profile_as_reference,             &
268               use_single_reference_value, u_gtrans, v_gtrans, virtual_flight, &
269               ws_scheme_mom, ws_scheme_sca
270
271    USE cpulog,                                                                &
272        ONLY:  cpu_log, log_point, log_point_s
273
274    USE flight_mod,                                                            &
275        ONLY:  flight_measurement
276
277
278    USE indices,                                                               &
279        ONLY:  i_left, i_right, j_north, j_south, nbgp, nx, nxl, nxlg, nxr,    &
280               nxrg, nyn, nyng, nys, nysg, nzb, nzt, nzb_u_inner, nzb_v_inner
281
282    USE interaction_droplets_ptq_mod,                                          &
283        ONLY:  interaction_droplets_ptq
284
285    USE interfaces
286
287    USE kinds
288
289    USE land_surface_model_mod,                                                &
290        ONLY:  land_surface, lsm_energy_balance, lsm_soil_model,               &
291               skip_time_do_lsm
292
293    USE ls_forcing_mod,                                                        &
294        ONLY:  ls_forcing_surf, ls_forcing_vert
295
296    USE microphysics_mod,                                                      &
297        ONLY: collision_turbulence
298
299    USE nudge_mod,                                                             &
300        ONLY:  calc_tnudge, nudge_ref
301
302    USE particle_attributes,                                                   &
303        ONLY:  particle_advection, particle_advection_start,                   &
304               use_sgs_for_particles, wang_kernel
305
306    USE pegrid
307
308    USE pmc_interface,                                                         &
309        ONLY:  nested_run, nesting_mode, pmci_datatrans,                       &
310               pmci_ensure_nest_mass_conservation, pmci_synchronize
311
312    USE production_e_mod,                                                      &
313        ONLY:  production_e_init
314
315    USE progress_bar,                                                          &
316        ONLY:  finish_progress_bar, output_progress_bar
317
318    USE prognostic_equations_mod,                                              &
319        ONLY:  prognostic_equations_acc, prognostic_equations_cache,           &
320               prognostic_equations_vector
321
322    USE radiation_model_mod,                                                   &
323        ONLY: dt_radiation, force_radiation_call, radiation, radiation_control,&
324              skip_time_do_radiation, time_radiation
325
326    USE spectra_mod,                                                           &
327        ONLY: average_count_sp, averaging_interval_sp, calc_spectra, dt_dosp,  &
328              skip_time_dosp
329
330    USE statistics,                                                            &
331        ONLY:  flow_statistics_called, hom, pr_palm, sums_ls_l, u_max,         &
332               u_max_ijk, v_max, v_max_ijk, w_max, w_max_ijk
333
334    USE surface_layer_fluxes_mod,                                              &
335        ONLY:  surface_layer_fluxes
336
337    USE user_actions_mod,                                                      &
338        ONLY:  user_actions
339
340    USE wind_turbine_model_mod,                                                &
341        ONLY:  wind_turbine, wtm_forces
342
343    IMPLICIT NONE
344
345    CHARACTER (LEN=9) ::  time_to_string          !<
346
347    REAL(wp) ::  dt_3d_old  !< temporary storage of timestep to be used for
348                            !< steering of run control output interval
349
350!
351!-- At beginning determine the first time step
352    CALL timestep
353
354!
355!-- Synchronize the timestep in case of nested run.
356    IF ( nested_run )  THEN
357!
358!--    Synchronization by unifying the time step.
359!--    Global minimum of all time-steps is used for all.
360       CALL pmci_synchronize
361    ENDIF
362
363!
364!-- Determine and print out the run control quantities before the first time
365!-- step of this run. For the initial run, some statistics (e.g. divergence)
366!-- need to be determined first.
367    IF ( simulated_time == 0.0_wp )  CALL flow_statistics
368    CALL run_control
369
370!
371!-- Data exchange between coupled models in case that a call has been omitted
372!-- at the end of the previous run of a job chain.
373    IF ( coupling_mode /= 'uncoupled'  .AND.  run_coupled )  THEN
374!
375!--    In case of model termination initiated by the local model the coupler
376!--    must not be called because this would again cause an MPI hang.
377       DO WHILE ( time_coupling >= dt_coupling  .AND.  terminate_coupled == 0 )
378          CALL surface_coupler
379          time_coupling = time_coupling - dt_coupling
380       ENDDO
381       IF (time_coupling == 0.0_wp  .AND.                                      &
382           time_since_reference_point < dt_coupling )                          &
383       THEN
384          time_coupling = time_since_reference_point
385       ENDIF
386    ENDIF
387
388#if defined( __dvrp_graphics )
389!
390!-- Time measurement with dvrp software 
391    CALL DVRP_LOG_EVENT( 2, current_timestep_number )
392#endif
393
394    CALL location_message( 'start with time-stepping', .TRUE. )
395!
396!-- Start of the time loop
397    DO  WHILE ( simulated_time < end_time  .AND.  .NOT. stop_dt  .AND. &
398                .NOT. terminate_run )
399
400       CALL cpu_log( log_point_s(10), 'timesteps', 'start' )
401
402!
403!--    Determine ug, vg and w_subs in dependence on data from external file
404!--    LSF_DATA
405       IF ( large_scale_forcing .AND. lsf_vert )  THEN
406           CALL ls_forcing_vert ( simulated_time )
407           sums_ls_l = 0.0_wp
408       ENDIF
409
410!
411!--    Set pt_init and q_init to the current profiles taken from
412!--    NUDGING_DATA
413       IF ( nudging )  THEN
414           CALL nudge_ref ( simulated_time )
415!
416!--        Store temperature gradient at the top boundary for possible Neumann
417!--        boundary condition
418           bc_pt_t_val = ( pt_init(nzt+1) - pt_init(nzt) ) / dzu(nzt+1)
419           bc_q_t_val  = ( q_init(nzt+1) - q_init(nzt) ) / dzu(nzt+1)
420       ENDIF
421
422!
423!--    Execute the user-defined actions
424       CALL user_actions( 'before_timestep' )
425
426!
427!--    Calculate forces by wind turbines
428       IF ( wind_turbine )  THEN
429
430          CALL cpu_log( log_point(55), 'wind_turbine', 'start' )
431
432          CALL wtm_forces
433
434          CALL cpu_log( log_point(55), 'wind_turbine', 'stop' )
435
436       ENDIF       
437       
438!
439!--    Start of intermediate step loop
440       intermediate_timestep_count = 0
441       DO  WHILE ( intermediate_timestep_count < &
442                   intermediate_timestep_count_max )
443
444          intermediate_timestep_count = intermediate_timestep_count + 1
445
446!
447!--       Set the steering factors for the prognostic equations which depend
448!--       on the timestep scheme
449          CALL timestep_scheme_steering
450
451!
452!--       Calculate those variables needed in the tendency terms which need
453!--       global communication
454          IF ( .NOT. use_single_reference_value  .AND. &
455               .NOT. use_initial_profile_as_reference )  THEN
456!
457!--          Horizontally averaged profiles to be used as reference state in
458!--          buoyancy terms (WARNING: only the respective last call of
459!--          calc_mean_profile defines the reference state!)
460             IF ( .NOT. neutral )  THEN
461                CALL calc_mean_profile( pt, 4 )
462                ref_state(:)  = hom(:,1,4,0) ! this is used in the buoyancy term
463             ENDIF
464             IF ( ocean )  THEN
465                CALL calc_mean_profile( rho, 64 )
466                ref_state(:)  = hom(:,1,64,0)
467             ENDIF
468             IF ( humidity )  THEN
469                CALL calc_mean_profile( vpt, 44 )
470                ref_state(:)  = hom(:,1,44,0)
471             ENDIF
472
473          ENDIF
474
475          IF ( .NOT. constant_diffusion )  CALL production_e_init
476          IF ( ( ws_scheme_mom .OR. ws_scheme_sca )  .AND.  &
477               intermediate_timestep_count == 1 )  CALL ws_statistics
478!
479!--       In case of nudging calculate current nudging time scale and horizontal
480!--       means of u, v, pt and q
481          IF ( nudging )  THEN
482             CALL calc_tnudge( simulated_time )
483             CALL calc_mean_profile( u, 1 )
484             CALL calc_mean_profile( v, 2 )
485             CALL calc_mean_profile( pt, 4 )
486             CALL calc_mean_profile( q, 41 )
487          ENDIF
488
489!
490!--       Solve the prognostic equations. A fast cache optimized version with
491!--       only one single loop is used in case of Piascek-Williams advection
492!--       scheme. NEC vector machines use a different version, because
493!--       in the other versions a good vectorization is prohibited due to
494!--       inlining problems.
495          IF ( loop_optimization == 'cache' )  THEN
496             CALL prognostic_equations_cache
497          ELSEIF ( loop_optimization == 'vector' )  THEN
498             CALL prognostic_equations_vector
499          ELSEIF ( loop_optimization == 'acc' )  THEN
500             i_left  = nxl;         i_right = nxr
501             j_south = nys;         j_north = nyn
502             CALL prognostic_equations_acc
503
504!             i_left  = nxl;         i_right = nxl+nbgp-1
505!             j_south = nys;         j_north = nyn
506!             CALL prognostic_equations_acc
507!             i_left  = nxr-nbgp+1;  i_right = nxr
508!             j_south = nys;         j_north = nyn
509!             CALL prognostic_equations_acc
510
511!
512!--          Exchange of ghost points (lateral boundary conditions)
513             IF ( background_communication )  THEN
514
515                CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'start' )
516               
517                send_receive = 'lr'
518                sendrecv_in_background = .TRUE.
519                req          = 0
520                req_count    = 0
521
522                IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
523                   on_device = .TRUE.         ! to be removed after complete porting
524                ELSE                          ! of ghost point exchange
525                   !$acc update host( e_p, pt_p, u_p, v_p, w_p )
526                ENDIF
527
528                CALL exchange_horiz( u_p, nbgp )
529                CALL exchange_horiz( v_p, nbgp )
530                CALL exchange_horiz( w_p, nbgp )
531                CALL exchange_horiz( pt_p, nbgp )
532                IF ( .NOT. constant_diffusion )  CALL exchange_horiz( e_p, nbgp )
533                IF ( ocean )  THEN
534                   CALL exchange_horiz( sa_p, nbgp )
535                   CALL exchange_horiz( rho, nbgp )
536                   CALL exchange_horiz( prho, nbgp )
537                ENDIF
538                IF ( humidity )  THEN
539                   CALL exchange_horiz( q_p, nbgp )
540                   IF ( cloud_physics .AND. microphysics_seifert )  THEN
541                      CALL exchange_horiz( qr_p, nbgp )
542                      CALL exchange_horiz( nr_p, nbgp )
543                   ENDIF
544                ENDIF
545                IF ( cloud_droplets )  THEN
546                   CALL exchange_horiz( ql, nbgp )
547                   CALL exchange_horiz( ql_c, nbgp )
548                   CALL exchange_horiz( ql_v, nbgp )
549                   CALL exchange_horiz( ql_vp, nbgp )
550                ENDIF
551                IF ( wang_kernel  .OR.  collision_turbulence  .OR.             &
552                     use_sgs_for_particles )  THEN
553                   CALL exchange_horiz( diss, nbgp )
554                ENDIF
555                IF ( passive_scalar )  CALL exchange_horiz( s_p, nbgp )
556
557                IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
558                   on_device = .FALSE.        ! to be removed after complete porting
559                ELSE                          ! of ghost point exchange
560                   !$acc update device( e_p, pt_p, u_p, v_p, w_p )
561                ENDIF
562
563                sendrecv_in_background = .FALSE.
564
565                CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'pause' )
566
567             ENDIF
568
569!             i_left  = nxl+nbgp;    i_right = nxr-nbgp
570!             j_south = nys;         j_north = nys+nbgp-1
571!             CALL prognostic_equations_acc
572!             i_left  = nxl+nbgp;    i_right = nxr-nbgp
573!             j_south = nyn-nbgp+1;  j_north = nyn
574!             CALL prognostic_equations_acc
575
576             IF ( background_communication )  THEN
577                CALL cpu_log( log_point(41), 'exchange-horiz-wait', 'start' )
578#if defined( __parallel )
579                CALL MPI_WAITALL( req_count, req, wait_stat, ierr )
580#endif
581                CALL cpu_log( log_point(41), 'exchange-horiz-wait', 'pause' )
582
583                CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'continue' )
584
585                send_receive = 'ns'
586                sendrecv_in_background = .TRUE.
587                req          = 0
588                req_count    = 0
589
590                IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
591                   on_device = .TRUE.         ! to be removed after complete porting
592                ELSE                          ! of ghost point exchange
593                   !$acc update host( e_p, pt_p, u_p, v_p, w_p )
594                ENDIF
595
596                CALL exchange_horiz( u_p, nbgp )
597                CALL exchange_horiz( v_p, nbgp )
598                CALL exchange_horiz( w_p, nbgp )
599                CALL exchange_horiz( pt_p, nbgp )
600                IF ( .NOT. constant_diffusion )  CALL exchange_horiz( e_p, nbgp )
601                IF ( ocean )  THEN
602                   CALL exchange_horiz( sa_p, nbgp )
603                   CALL exchange_horiz( rho, nbgp )
604                  CALL exchange_horiz( prho, nbgp )
605                ENDIF
606                IF ( humidity )  THEN
607                   CALL exchange_horiz( q_p, nbgp )
608                   IF ( cloud_physics .AND. microphysics_seifert )  THEN
609                      CALL exchange_horiz( qr_p, nbgp )
610                      CALL exchange_horiz( nr_p, nbgp )
611                   ENDIF
612                ENDIF
613                IF ( cloud_droplets )  THEN
614                   CALL exchange_horiz( ql, nbgp )
615                   CALL exchange_horiz( ql_c, nbgp )
616                   CALL exchange_horiz( ql_v, nbgp )
617                   CALL exchange_horiz( ql_vp, nbgp )
618                ENDIF
619                IF ( wang_kernel  .OR.  collision_turbulence  .OR.             &
620                     use_sgs_for_particles )  THEN
621                   CALL exchange_horiz( diss, nbgp )
622                ENDIF
623                IF ( passive_scalar )  CALL exchange_horiz( s_p, nbgp )
624
625                IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
626                   on_device = .FALSE.        ! to be removed after complete porting
627                ELSE                          ! of ghost point exchange
628                   !$acc update device( e_p, pt_p, u_p, v_p, w_p )
629                ENDIF
630
631                sendrecv_in_background = .FALSE.
632
633                CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'stop' )
634
635             ENDIF
636
637!             i_left  = nxl+nbgp;    i_right = nxr-nbgp
638!             j_south = nys+nbgp;    j_north = nyn-nbgp
639!             CALL prognostic_equations_acc
640
641             IF ( background_communication )  THEN
642                CALL cpu_log( log_point(41), 'exchange-horiz-wait', 'continue' )
643#if defined( __parallel )
644                CALL MPI_WAITALL( req_count, req, wait_stat, ierr )
645#endif
646                send_receive = 'al'
647                CALL cpu_log( log_point(41), 'exchange-horiz-wait', 'stop' )
648             ENDIF
649
650          ENDIF
651
652!
653!--       Particle transport/physics with the Lagrangian particle model
654!--       (only once during intermediate steps, because it uses an Euler-step)
655!--       ### particle model should be moved before prognostic_equations, in order
656!--       to regard droplet interactions directly
657          IF ( particle_advection  .AND.                         &
658               simulated_time >= particle_advection_start  .AND. &
659               intermediate_timestep_count == 1 )  THEN
660             CALL lpm
661             first_call_lpm = .FALSE.
662          ENDIF
663
664!
665!--       Interaction of droplets with temperature and specific humidity.
666!--       Droplet condensation and evaporation is calculated within
667!--       advec_particles.
668          IF ( cloud_droplets  .AND.  &
669               intermediate_timestep_count == intermediate_timestep_count_max )&
670          THEN
671             CALL interaction_droplets_ptq
672          ENDIF
673
674!
675!--       Exchange of ghost points (lateral boundary conditions)
676          IF ( .NOT. background_communication )  THEN
677
678             CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'start' )
679
680             IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
681                on_device = .TRUE.         ! to be removed after complete porting
682             ELSE                          ! of ghost point exchange
683                !$acc update host( e_p, pt_p, u_p, v_p, w_p )
684             ENDIF
685
686             CALL exchange_horiz( u_p, nbgp )
687             CALL exchange_horiz( v_p, nbgp )
688             CALL exchange_horiz( w_p, nbgp )
689             CALL exchange_horiz( pt_p, nbgp )
690             IF ( .NOT. constant_diffusion )  CALL exchange_horiz( e_p, nbgp )
691             IF ( ocean )  THEN
692                CALL exchange_horiz( sa_p, nbgp )
693                CALL exchange_horiz( rho, nbgp )
694                CALL exchange_horiz( prho, nbgp )
695             ENDIF
696             IF ( humidity )  THEN
697                CALL exchange_horiz( q_p, nbgp )
698                IF ( cloud_physics .AND. microphysics_seifert )  THEN
699                   CALL exchange_horiz( qr_p, nbgp )
700                   CALL exchange_horiz( nr_p, nbgp )
701                ENDIF
702             ENDIF
703             IF ( cloud_droplets )  THEN
704                CALL exchange_horiz( ql, nbgp )
705                CALL exchange_horiz( ql_c, nbgp )
706                CALL exchange_horiz( ql_v, nbgp )
707                CALL exchange_horiz( ql_vp, nbgp )
708             ENDIF
709             IF ( wang_kernel  .OR.  collision_turbulence  .OR.                &
710                  use_sgs_for_particles )  THEN
711                CALL exchange_horiz( diss, nbgp )
712             ENDIF
713             IF ( passive_scalar )  CALL exchange_horiz( s_p, nbgp )     
714
715             IF ( numprocs == 1 )  THEN    ! workaround for single-core GPU runs
716                on_device = .FALSE.        ! to be removed after complete porting
717             ELSE                          ! of ghost point exchange
718                !$acc update device( e_p, pt_p, u_p, v_p, w_p )
719             ENDIF
720
721             CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'stop' )
722
723          ENDIF
724
725!
726!--       Boundary conditions for the prognostic quantities (except of the
727!--       velocities at the outflow in case of a non-cyclic lateral wall)
728          CALL boundary_conds
729
730!
731!--       Swap the time levels in preparation for the next time step.
732          CALL swap_timelevel
733
734          IF ( nested_run )  THEN
735
736             CALL cpu_log( log_point(60), 'nesting', 'start' )
737!
738!--          Domain nesting. The data transfer subroutines pmci_parent_datatrans
739!--          and pmci_child_datatrans are called inside the wrapper
740!--          subroutine pmci_datatrans according to the control parameters
741!--          nesting_mode and nesting_datatransfer_mode.
742!--          TO_DO: why is nesting_mode given as a parameter here?
743             CALL pmci_datatrans( nesting_mode )
744
745             IF ( TRIM( nesting_mode ) == 'two-way' .OR.                               &
746                  nesting_mode == 'vertical' )  THEN
747!
748!--             Exchange_horiz is needed for all parent-domains after the
749!--             anterpolation
750                CALL exchange_horiz( u, nbgp )
751                CALL exchange_horiz( v, nbgp )
752                CALL exchange_horiz( w, nbgp )
753                IF ( .NOT. neutral            )  CALL exchange_horiz( pt, nbgp )
754                IF ( humidity                 )  CALL exchange_horiz( q, nbgp  )
755                IF ( passive_scalar           )  CALL exchange_horiz( s, nbgp  )
756                IF ( .NOT. constant_diffusion )  CALL exchange_horiz( e, nbgp  )
757             ENDIF
758!
759!--          Correct the w top-BC in nest domains to ensure mass conservation.
760!--          This action must never be done for the root domain. Vertical
761!--          nesting implies mass conservation.
762             IF ( nest_domain )  THEN
763                CALL pmci_ensure_nest_mass_conservation
764             ENDIF
765
766             CALL cpu_log( log_point(60), 'nesting', 'stop' )
767
768          ENDIF
769
770!
771!--       Temperature offset must be imposed at cyclic boundaries in x-direction
772!--       when a sloping surface is used
773          IF ( sloping_surface )  THEN
774             IF ( nxl ==  0 )  pt(:,:,nxlg:nxl-1) = pt(:,:,nxlg:nxl-1) - &
775                                                    pt_slope_offset
776             IF ( nxr == nx )  pt(:,:,nxr+1:nxrg) = pt(:,:,nxr+1:nxrg) + &
777                                                    pt_slope_offset
778          ENDIF
779
780!
781!--       Impose a turbulent inflow using the recycling method
782          IF ( turbulent_inflow )  CALL  inflow_turbulence
783
784!
785!--       Impose a random perturbation on the horizontal velocity field
786          IF ( create_disturbances  .AND.                                      &
787               ( call_psolver_at_all_substeps  .AND.                           &
788               intermediate_timestep_count == intermediate_timestep_count_max )&
789          .OR. ( .NOT. call_psolver_at_all_substeps  .AND.                     &
790               intermediate_timestep_count == 1 ) )                            &
791          THEN
792             time_disturb = time_disturb + dt_3d
793             IF ( time_disturb >= dt_disturb )  THEN
794                !$acc update host( u, v )
795                IF ( numprocs == 1 )  on_device = .FALSE.  ! workaround, remove later
796                IF ( disturbance_energy_limit /= 0.0_wp  .AND.                 &
797                     hom(nzb+5,1,pr_palm,0) < disturbance_energy_limit )  THEN
798                   CALL disturb_field( nzb_u_inner, tend, u )
799                   CALL disturb_field( nzb_v_inner, tend, v )
800                ELSEIF ( .NOT. bc_lr_cyc  .OR.  .NOT. bc_ns_cyc )  THEN
801!
802!--                Runs with a non-cyclic lateral wall need perturbations
803!--                near the inflow throughout the whole simulation
804                   dist_range = 1
805                   CALL disturb_field( nzb_u_inner, tend, u )
806                   CALL disturb_field( nzb_v_inner, tend, v )
807                   dist_range = 0
808                ENDIF
809                IF ( numprocs == 1 )  on_device = .TRUE.  ! workaround, remove later
810                !$acc update device( u, v )
811                time_disturb = time_disturb - dt_disturb
812             ENDIF
813          ENDIF
814
815!
816!--       Reduce the velocity divergence via the equation for perturbation
817!--       pressure.
818          IF ( intermediate_timestep_count == 1  .OR. &
819                call_psolver_at_all_substeps )  THEN
820             CALL pres
821          ENDIF
822
823!
824!--       If required, compute liquid water content
825          IF ( cloud_physics )  THEN
826             CALL calc_liquid_water_content
827             !$acc update device( ql )
828          ENDIF
829!
830!--       If required, compute virtual potential temperature
831          IF ( humidity )  THEN
832             CALL compute_vpt 
833             !$acc update device( vpt )
834          ENDIF 
835
836!
837!--       Compute the diffusion quantities
838          IF ( .NOT. constant_diffusion )  THEN
839
840!
841!--          Determine surface fluxes shf and qsws and surface values
842!--          pt_surface and q_surface in dependence on data from external
843!--          file LSF_DATA respectively
844             IF ( ( large_scale_forcing .AND. lsf_surf ) .AND. &
845                 intermediate_timestep_count == intermediate_timestep_count_max )&
846             THEN
847                CALL ls_forcing_surf ( simulated_time )
848             ENDIF
849
850!
851!--          First the vertical fluxes in the surface (constant flux) layer are computed
852             IF ( constant_flux_layer )  THEN
853                CALL cpu_log( log_point(19), 'surface_layer_fluxes', 'start' )
854                CALL surface_layer_fluxes
855                CALL cpu_log( log_point(19), 'surface_layer_fluxes', 'stop' )
856             ENDIF
857
858!
859!--          If required, solve the energy balance for the surface and run soil
860!--          model
861             IF ( land_surface .AND. simulated_time > skip_time_do_lsm)  THEN
862
863                CALL cpu_log( log_point(54), 'land_surface', 'start' )
864                CALL lsm_energy_balance
865                CALL lsm_soil_model
866                CALL cpu_log( log_point(54), 'land_surface', 'stop' )
867             ENDIF
868!
869!--          Compute the diffusion coefficients
870             CALL cpu_log( log_point(17), 'diffusivities', 'start' )
871             IF ( .NOT. humidity ) THEN
872                IF ( ocean )  THEN
873                   CALL diffusivities( prho, prho_reference )
874                ELSE
875                   CALL diffusivities( pt, pt_reference )
876                ENDIF
877             ELSE
878                CALL diffusivities( vpt, pt_reference )
879             ENDIF
880             CALL cpu_log( log_point(17), 'diffusivities', 'stop' )
881
882          ENDIF
883
884!
885!--       If required, calculate radiative fluxes and heating rates
886          IF ( radiation .AND. intermediate_timestep_count                     &
887               == intermediate_timestep_count_max .AND. simulated_time >    &
888               skip_time_do_radiation )  THEN
889
890               time_radiation = time_radiation + dt_3d
891
892             IF ( time_radiation >= dt_radiation .OR. force_radiation_call )   &
893             THEN
894
895                CALL cpu_log( log_point(50), 'radiation', 'start' )
896
897                IF ( .NOT. force_radiation_call )  THEN
898                   time_radiation = time_radiation - dt_radiation
899                ENDIF
900
901                CALL radiation_control
902
903                CALL cpu_log( log_point(50), 'radiation', 'stop' )
904             ENDIF
905          ENDIF
906
907       ENDDO   ! Intermediate step loop
908
909!
910!--    Increase simulation time and output times
911       nr_timesteps_this_run      = nr_timesteps_this_run + 1
912       current_timestep_number    = current_timestep_number + 1
913       simulated_time             = simulated_time   + dt_3d
914       simulated_time_chr         = time_to_string( simulated_time )
915       time_since_reference_point = simulated_time - coupling_start_time
916
917
918
919       IF ( simulated_time >= skip_time_data_output_av )  THEN
920          time_do_av         = time_do_av       + dt_3d
921       ENDIF
922       IF ( simulated_time >= skip_time_do2d_xy )  THEN
923          time_do2d_xy       = time_do2d_xy     + dt_3d
924       ENDIF
925       IF ( simulated_time >= skip_time_do2d_xz )  THEN
926          time_do2d_xz       = time_do2d_xz     + dt_3d
927       ENDIF
928       IF ( simulated_time >= skip_time_do2d_yz )  THEN
929          time_do2d_yz       = time_do2d_yz     + dt_3d
930       ENDIF
931       IF ( simulated_time >= skip_time_do3d    )  THEN
932          time_do3d          = time_do3d        + dt_3d
933       ENDIF
934       DO  mid = 1, masks
935          IF ( simulated_time >= skip_time_domask(mid) )  THEN
936             time_domask(mid)= time_domask(mid) + dt_3d
937          ENDIF
938       ENDDO
939       time_dvrp          = time_dvrp        + dt_3d
940       IF ( simulated_time >= skip_time_dosp )  THEN
941          time_dosp       = time_dosp        + dt_3d
942       ENDIF
943       time_dots          = time_dots        + dt_3d
944       IF ( .NOT. first_call_lpm )  THEN
945          time_dopts      = time_dopts       + dt_3d
946       ENDIF
947       IF ( simulated_time >= skip_time_dopr )  THEN
948          time_dopr       = time_dopr        + dt_3d
949       ENDIF
950       time_dopr_listing          = time_dopr_listing        + dt_3d
951       time_run_control   = time_run_control + dt_3d
952
953!
954!--    Data exchange between coupled models
955       IF ( coupling_mode /= 'uncoupled'  .AND.  run_coupled )  THEN
956          time_coupling = time_coupling + dt_3d
957
958!
959!--       In case of model termination initiated by the local model
960!--       (terminate_coupled > 0), the coupler must be skipped because it would
961!--       cause an MPI intercomminucation hang.
962!--       If necessary, the coupler will be called at the beginning of the
963!--       next restart run.
964          DO WHILE ( time_coupling >= dt_coupling .AND. terminate_coupled == 0 )
965             CALL surface_coupler
966             time_coupling = time_coupling - dt_coupling
967          ENDDO
968       ENDIF
969
970!
971!--    Execute user-defined actions
972       CALL user_actions( 'after_integration' )
973
974!
975!--    If Galilei transformation is used, determine the distance that the
976!--    model has moved so far
977       IF ( galilei_transformation )  THEN
978          advected_distance_x = advected_distance_x + u_gtrans * dt_3d
979          advected_distance_y = advected_distance_y + v_gtrans * dt_3d
980       ENDIF
981
982!
983!--    Check, if restart is necessary (because cpu-time is expiring or
984!--    because it is forced by user) and set stop flag
985!--    This call is skipped if the remote model has already initiated a restart.
986       IF ( .NOT. terminate_run )  CALL check_for_restart
987
988!
989!--    Carry out statistical analysis and output at the requested output times.
990!--    The MOD function is used for calculating the output time counters (like
991!--    time_dopr) in order to regard a possible decrease of the output time
992!--    interval in case of restart runs
993
994!
995!--    Set a flag indicating that so far no statistics have been created
996!--    for this time step
997       flow_statistics_called = .FALSE.
998
999!
1000!--    If required, call flow_statistics for averaging in time
1001       IF ( averaging_interval_pr /= 0.0_wp  .AND.  &
1002            ( dt_dopr - time_dopr ) <= averaging_interval_pr  .AND.  &
1003            simulated_time >= skip_time_dopr )  THEN
1004          time_dopr_av = time_dopr_av + dt_3d
1005          IF ( time_dopr_av >= dt_averaging_input_pr )  THEN
1006             do_sum = .TRUE.
1007             time_dopr_av = MOD( time_dopr_av, &
1008                                    MAX( dt_averaging_input_pr, dt_3d ) )
1009          ENDIF
1010       ENDIF
1011       IF ( do_sum )  CALL flow_statistics
1012
1013!
1014!--    Sum-up 3d-arrays for later output of time-averaged 2d/3d/masked data
1015       IF ( averaging_interval /= 0.0_wp  .AND.                                &
1016            ( dt_data_output_av - time_do_av ) <= averaging_interval  .AND. &
1017            simulated_time >= skip_time_data_output_av )                    &
1018       THEN
1019          time_do_sla = time_do_sla + dt_3d
1020          IF ( time_do_sla >= dt_averaging_input )  THEN
1021             CALL sum_up_3d_data
1022             average_count_3d = average_count_3d + 1
1023             time_do_sla = MOD( time_do_sla, MAX( dt_averaging_input, dt_3d ) )
1024          ENDIF
1025       ENDIF
1026
1027!
1028!--    Calculate spectra for time averaging
1029       IF ( averaging_interval_sp /= 0.0_wp  .AND.  &
1030            ( dt_dosp - time_dosp ) <= averaging_interval_sp  .AND.  &
1031            simulated_time >= skip_time_dosp )  THEN
1032          time_dosp_av = time_dosp_av + dt_3d
1033          IF ( time_dosp_av >= dt_averaging_input_pr )  THEN
1034             CALL calc_spectra
1035             time_dosp_av = MOD( time_dosp_av, &
1036                                 MAX( dt_averaging_input_pr, dt_3d ) )
1037          ENDIF
1038       ENDIF
1039
1040!
1041!--    Call flight module and output data
1042       IF ( virtual_flight )  THEN
1043          CALL flight_measurement
1044          CALL data_output_flight
1045       ENDIF
1046
1047!
1048!--    Profile output (ASCII) on file
1049       IF ( time_dopr_listing >= dt_dopr_listing )  THEN
1050          CALL print_1d
1051          time_dopr_listing = MOD( time_dopr_listing, MAX( dt_dopr_listing, &
1052                                                           dt_3d ) )
1053       ENDIF
1054
1055!
1056!--    Graphic output for PROFIL
1057       IF ( time_dopr >= dt_dopr )  THEN
1058          IF ( dopr_n /= 0 )  CALL data_output_profiles
1059          time_dopr = MOD( time_dopr, MAX( dt_dopr, dt_3d ) )
1060          time_dopr_av = 0.0_wp    ! due to averaging (see above)
1061       ENDIF
1062
1063!
1064!--    Graphic output for time series
1065       IF ( time_dots >= dt_dots )  THEN
1066          CALL data_output_tseries
1067          time_dots = MOD( time_dots, MAX( dt_dots, dt_3d ) )
1068       ENDIF
1069
1070!
1071!--    Output of spectra (formatted for use with PROFIL), in case of no
1072!--    time averaging, spectra has to be calculated before
1073       IF ( time_dosp >= dt_dosp )  THEN
1074          IF ( average_count_sp == 0 )  CALL calc_spectra
1075          CALL data_output_spectra
1076          time_dosp = MOD( time_dosp, MAX( dt_dosp, dt_3d ) )
1077       ENDIF
1078
1079!
1080!--    2d-data output (cross-sections)
1081       IF ( time_do2d_xy >= dt_do2d_xy )  THEN
1082          CALL data_output_2d( 'xy', 0 )
1083          time_do2d_xy = MOD( time_do2d_xy, MAX( dt_do2d_xy, dt_3d ) )
1084       ENDIF
1085       IF ( time_do2d_xz >= dt_do2d_xz )  THEN
1086          CALL data_output_2d( 'xz', 0 )
1087          time_do2d_xz = MOD( time_do2d_xz, MAX( dt_do2d_xz, dt_3d ) )
1088       ENDIF
1089       IF ( time_do2d_yz >= dt_do2d_yz )  THEN
1090          CALL data_output_2d( 'yz', 0 )
1091          time_do2d_yz = MOD( time_do2d_yz, MAX( dt_do2d_yz, dt_3d ) )
1092       ENDIF
1093
1094!
1095!--    3d-data output (volume data)
1096       IF ( time_do3d >= dt_do3d )  THEN
1097          CALL data_output_3d( 0 )
1098          time_do3d = MOD( time_do3d, MAX( dt_do3d, dt_3d ) )
1099       ENDIF
1100
1101!
1102!--    Masked data output
1103       DO  mid = 1, masks
1104          IF ( time_domask(mid) >= dt_domask(mid) )  THEN
1105             CALL data_output_mask( 0 )
1106             time_domask(mid) = MOD( time_domask(mid),  &
1107                                     MAX( dt_domask(mid), dt_3d ) )
1108          ENDIF
1109       ENDDO
1110
1111!
1112!--    Output of time-averaged 2d/3d/masked data
1113       IF ( time_do_av >= dt_data_output_av )  THEN
1114          CALL average_3d_data
1115          CALL data_output_2d( 'xy', 1 )
1116          CALL data_output_2d( 'xz', 1 )
1117          CALL data_output_2d( 'yz', 1 )
1118          CALL data_output_3d( 1 )
1119          DO  mid = 1, masks
1120             CALL data_output_mask( 1 )
1121          ENDDO
1122          time_do_av = MOD( time_do_av, MAX( dt_data_output_av, dt_3d ) )
1123       ENDIF
1124
1125!
1126!--    Output of particle time series
1127       IF ( particle_advection )  THEN
1128          IF ( time_dopts >= dt_dopts  .OR. &
1129               ( simulated_time >= particle_advection_start  .AND. &
1130                 first_call_lpm ) )  THEN
1131             CALL data_output_ptseries
1132             time_dopts = MOD( time_dopts, MAX( dt_dopts, dt_3d ) )
1133          ENDIF
1134       ENDIF
1135
1136!
1137!--    Output of dvrp-graphics (isosurface, particles, slicer)
1138#if defined( __dvrp_graphics )
1139       CALL DVRP_LOG_EVENT( -2, current_timestep_number-1 )
1140#endif
1141       IF ( time_dvrp >= dt_dvrp )  THEN
1142          CALL data_output_dvrp
1143          time_dvrp = MOD( time_dvrp, MAX( dt_dvrp, dt_3d ) )
1144       ENDIF
1145#if defined( __dvrp_graphics )
1146       CALL DVRP_LOG_EVENT( 2, current_timestep_number )
1147#endif
1148
1149!
1150!--    If required, set the heat flux for the next time step at a random value
1151       IF ( constant_heatflux  .AND.  random_heatflux )  CALL disturb_heatflux
1152
1153!
1154!--    Execute user-defined actions
1155       CALL user_actions( 'after_timestep' )
1156
1157!
1158!--    Determine size of next time step. Save timestep dt_3d because it is
1159!--    newly calculated in routine timestep, but required further below for
1160!--    steering the run control output interval
1161       dt_3d_old = dt_3d
1162       CALL timestep
1163
1164!
1165!--    Synchronize the timestep in case of nested run.
1166       IF ( nested_run )  THEN
1167!
1168!--       Synchronize by unifying the time step.
1169!--       Global minimum of all time-steps is used for all.
1170          CALL pmci_synchronize
1171       ENDIF
1172
1173!
1174!--    Computation and output of run control parameters.
1175!--    This is also done whenever perturbations have been imposed
1176       IF ( time_run_control >= dt_run_control  .OR.                     &
1177            timestep_scheme(1:5) /= 'runge'  .OR.  disturbance_created ) &
1178       THEN
1179!          IF ( current_timestep_number == 1 )  THEN
1180!             IF ( nxl < 7 .AND.  nxr > 7 .AND. nys < 7 .AND. nyn > 7 )  THEN
1181!                u(10,7,7) = 0.55
1182!             ENDIF
1183!             PRINT*, 'calculating minmax'
1184!             CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, u,       &
1185!                                  'abs', 0.0_wp, u_max, u_max_ijk )
1186!             CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, v,       &
1187!                                  'abs', 0.0_wp, v_max, v_max_ijk )
1188!             CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, w,       &
1189!                                  'abs', 0.0_wp, w_max, w_max_ijk )
1190!             PRINT*, 'calculated u_max = ', u_max, '  myid = ', myid
1191!          ENDIF
1192          CALL run_control
1193          IF ( time_run_control >= dt_run_control )  THEN
1194             time_run_control = MOD( time_run_control, &
1195                                     MAX( dt_run_control, dt_3d_old ) )
1196          ENDIF
1197       ENDIF
1198
1199!
1200!--    Output elapsed simulated time in form of a progress bar on stdout
1201       IF ( myid == 0 )  CALL output_progress_bar
1202
1203       CALL cpu_log( log_point_s(10), 'timesteps', 'stop' )
1204
1205
1206    ENDDO   ! time loop
1207
1208    IF ( myid == 0 )  CALL finish_progress_bar
1209
1210#if defined( __dvrp_graphics )
1211    CALL DVRP_LOG_EVENT( -2, current_timestep_number )
1212#endif
1213
1214    CALL location_message( 'finished time-stepping', .TRUE. )
1215
1216 END SUBROUTINE time_integration
Note: See TracBrowser for help on using the repository browser.