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

Last change on this file since 1928 was 1928, checked in by hellstea, 8 years ago

last commit documented

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