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

Last change on this file since 2000 was 2000, checked in by knoop, 8 years ago

Forced header and separation lines into 80 columns

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