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

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

your message

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