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

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

Error check for overlapping parallel nests added

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