Home

Context Navigation

prognostic_equations.f90 @ 1484

Last change on this file since 1484 was 1484, checked in by kanani, 9 years ago

New:
---
Subroutine init_plant_canopy added to module plant_canopy_model_mod. (plant_canopy_model)
Alternative method for lad-profile construction added, also, new parameters added.
(header, package_parin, plant_canopy_model, read_var_list, write_var_list)
plant_canopy_model-dependency added to several subroutines. (Makefile)
New package/namelist canopy_par for canopy-related parameters added. (package_parin)

Changed:
---
Code structure of the plant canopy model changed, all canopy-model related code
combined to module plant_canopy_model_mod. (check_parameters, init_3d_model,
modules, timestep)
Module plant_canopy_model_mod added in USE-lists of some subroutines. (check_parameters,
header, init_3d_model, package_parin, read_var_list, user_init_plant_canopy, write_var_list)
Canopy initialization moved to new subroutine init_plant_canopy. (check_parameters,
init_3d_model, plant_canopy_model)
Calculation of canopy timestep-criterion removed, instead, the canopy
drag is now directly limited in the calculation of the canopy tendency terms.
(plant_canopy_model, timestep)
Some parameters renamed. (check_parameters, header, init_plant_canopy,
plant_canopy_model, read_var_list, write_var_list)
Unnecessary 3d-arrays removed. (init_plant_canopy, plant_canopy_model, user_init_plant_canopy)
Parameter checks regarding canopy initialization added. (check_parameters)
All canopy steering parameters moved from namelist inipar to canopy_par. (package_parin, parin)
Some redundant MPI communication removed. (init_plant_canopy)

Bugfix:
---
Missing KIND-attribute for REAL constant added. (check_parameters)
DO-WHILE-loop for lad-profile output restricted. (header)
Removed double-listing of use_upstream_for_tke in ONLY-list of module
control_parameters. (prognostic_equations)

Property svn:keywords set to Id

File size: 83.7 KB

Line
1	MODULE prognostic_equations_mod
2
3	!--------------------------------------------------------------------------------!
4	! This file is part of PALM.
5	!
6	! PALM is free software: you can redistribute it and/or modify it under the terms
7	! of the GNU General Public License as published by the Free Software Foundation,
8	! either version 3 of the License, or (at your option) any later 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-2014 Leibniz Universitaet Hannover
18	!--------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! ------------------
22	! Changes due to new module structure of the plant canopy model:
23	! parameters cthf and plant_canopy moved to module plant_canopy_model_mod.
24	! Removed double-listing of use_upstream_for_tke in ONLY-list of module
25	! control_parameters
26	!
27	! Former revisions:
28	! -----------------
29	! $Id: prognostic_equations.f90 1484 2014-10-21 10:53:05Z kanani $
30	!
31	! 1409 2014-05-23 12:11:32Z suehring
32	! Bugfix: i_omp_start changed for advec_u_ws at left inflow and outflow boundary.
33	! This ensures that left-hand side fluxes are also calculated for nxl in that
34	! case, even though the solution at nxl is overwritten in boundary_conds()
35	!
36	! 1398 2014-05-07 11:15:00Z heinze
37	! Rayleigh-damping for horizontal velocity components changed: instead of damping
38	! against ug and vg, damping against u_init and v_init is used to allow for a
39	! homogenized treatment in case of nudging
40	!
41	! 1380 2014-04-28 12:40:45Z heinze
42	! Change order of calls for scalar prognostic quantities:
43	! ls_advec -> nudging -> subsidence since initial profiles
44	!
45	! 1374 2014-04-25 12:55:07Z raasch
46	! missing variables added to ONLY lists
47	!
48	! 1365 2014-04-22 15:03:56Z boeske
49	! Calls of ls_advec for large scale advection added,
50	! subroutine subsidence is only called if use_subsidence_tendencies = .F.,
51	! new argument ls_index added to the calls of subsidence
52	! +ls_index
53	!
54	! 1361 2014-04-16 15:17:48Z hoffmann
55	! Two-moment microphysics moved to the start of prognostic equations. This makes
56	! the 3d arrays for tend_q, tend_qr, tend_pt and tend_pt redundant.
57	! Additionally, it is allowed to call the microphysics just once during the time
58	! step (not at each sub-time step).
59	!
60	! Two-moment cloud physics added for vector and accelerator optimization.
61	!
62	! 1353 2014-04-08 15:21:23Z heinze
63	! REAL constants provided with KIND-attribute
64	!
65	! 1337 2014-03-25 15:11:48Z heinze
66	! Bugfix: REAL constants provided with KIND-attribute
67	!
68	! 1332 2014-03-25 11:59:43Z suehring
69	! Bugfix: call advec_ws or advec_pw for TKE only if NOT use_upstream_for_tke
70	!
71	! 1330 2014-03-24 17:29:32Z suehring
72	! In case of SGS-particle velocity advection of TKE is also allowed with
73	! dissipative 5th-order scheme.
74	!
75	! 1320 2014-03-20 08:40:49Z raasch
76	! ONLY-attribute added to USE-statements,
77	! kind-parameters added to all INTEGER and REAL declaration statements,
78	! kinds are defined in new module kinds,
79	! old module precision_kind is removed,
80	! revision history before 2012 removed,
81	! comment fields (!:) to be used for variable explanations added to
82	! all variable declaration statements
83	!
84	! 1318 2014-03-17 13:35:16Z raasch
85	! module interfaces removed
86	!
87	! 1257 2013-11-08 15:18:40Z raasch
88	! openacc loop vector clauses removed, independent clauses added
89	!
90	! 1246 2013-11-01 08:59:45Z heinze
91	! enable nudging also for accelerator version
92	!
93	! 1241 2013-10-30 11:36:58Z heinze
94	! usage of nudging enabled (so far not implemented for accelerator version)
95	!
96	! 1179 2013-06-14 05:57:58Z raasch
97	! two arguments removed from routine buoyancy, ref_state updated on device
98	!
99	! 1128 2013-04-12 06:19:32Z raasch
100	! those parts requiring global communication moved to time_integration,
101	! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
102	! j_north
103	!
104	! 1115 2013-03-26 18:16:16Z hoffmann
105	! optimized cloud physics: calculation of microphysical tendencies transfered
106	! to microphysics.f90; qr and nr are only calculated if precipitation is required
107	!
108	! 1111 2013-03-08 23:54:10Z raasch
109	! update directives for prognostic quantities removed
110	!
111	! 1106 2013-03-04 05:31:38Z raasch
112	! small changes in code formatting
113	!
114	! 1092 2013-02-02 11:24:22Z raasch
115	! unused variables removed
116	!
117	! 1053 2012-11-13 17:11:03Z hoffmann
118	! implementation of two new prognostic equations for rain drop concentration (nr)
119	! and rain water content (qr)
120	!
121	! currently, only available for cache loop optimization
122	!
123	! 1036 2012-10-22 13:43:42Z raasch
124	! code put under GPL (PALM 3.9)
125	!
126	! 1019 2012-09-28 06:46:45Z raasch
127	! non-optimized version of prognostic_equations removed
128	!
129	! 1015 2012-09-27 09:23:24Z raasch
130	! new branch prognostic_equations_acc
131	! OpenACC statements added + code changes required for GPU optimization
132	!
133	! 1001 2012-09-13 14:08:46Z raasch
134	! all actions concerning leapfrog- and upstream-spline-scheme removed
135	!
136	! 978 2012-08-09 08:28:32Z fricke
137	! km_damp_x and km_damp_y removed in calls of diffusion_u and diffusion_v
138	! add ptdf_x, ptdf_y for damping the potential temperature at the inflow
139	! boundary in case of non-cyclic lateral boundaries
140	! Bugfix: first thread index changes for WS-scheme at the inflow
141	!
142	! 940 2012-07-09 14:31:00Z raasch
143	! temperature equation can be switched off
144	!
145	! Revision 1.1 2000/04/13 14:56:27 schroeter
146	! Initial revision
147	!
148	!
149	! Description:
150	! ------------
151	! Solving the prognostic equations.
152	!------------------------------------------------------------------------------!
153
154	USE arrays_3d, &
155	ONLY: diss_l_e, diss_l_nr, diss_l_pt, diss_l_q, diss_l_qr, &
156	diss_l_sa, diss_s_e, diss_s_nr, diss_s_pt, diss_s_q, &
157	diss_s_qr, diss_s_sa, e, e_p, flux_s_e, flux_s_nr, flux_s_pt, &
158	flux_s_q, flux_s_qr, flux_s_sa, flux_l_e, flux_l_nr, &
159	flux_l_pt, flux_l_q, flux_l_qr, flux_l_sa, nr, nr_p, nrsws, &
160	nrswst, pt, ptdf_x, ptdf_y, pt_init, pt_p, prho, q, q_init, &
161	q_p, qsws, qswst, qr, qr_p, qrsws, qrswst, rdf, rdf_sc, &
162	ref_state, rho, sa, sa_init, sa_p, saswsb, saswst, shf, tend, &
163	te_m, tnr_m, tpt_m, tq_m, tqr_m, tsa_m, tswst, tu_m, tv_m, &
164	tw_m, u, ug, u_init, u_p, v, vg, vpt, v_init, v_p, w, w_p
165
166	USE control_parameters, &
167	ONLY: call_microphysics_at_all_substeps, cloud_physics, &
168	constant_diffusion, dp_external, &
169	dp_level_ind_b, dp_smooth_factor, dpdxy, dt_3d, humidity, &
170	icloud_scheme, inflow_l, intermediate_timestep_count, &
171	intermediate_timestep_count_max, large_scale_forcing, &
172	large_scale_subsidence, neutral, nudging, ocean, outflow_l, &
173	outflow_s, passive_scalar, precipitation, &
174	prho_reference, prho_reference, prho_reference, pt_reference, &
175	pt_reference, pt_reference, radiation, scalar_advec, &
176	scalar_advec, simulated_time, sloping_surface, timestep_scheme, &
177	tsc, use_subsidence_tendencies, use_upstream_for_tke, &
178	wall_heatflux, &
179	wall_nrflux, wall_qflux, wall_qflux, wall_qflux, wall_qrflux, &
180	wall_salinityflux, ws_scheme_mom, ws_scheme_sca
181
182	USE cpulog, &
183	ONLY: cpu_log, log_point
184
185	USE eqn_state_seawater_mod, &
186	ONLY: eqn_state_seawater
187
188	USE indices, &
189	ONLY: i_left, i_right, j_north, j_south, nxl, nxlu, nxr, nyn, nys, &
190	nysv, nzb_s_inner, nzb_u_inner, nzb_v_inner, nzb_w_inner, nzt
191
192	USE advec_ws, &
193	ONLY: advec_s_ws, advec_s_ws_acc, advec_u_ws, advec_u_ws_acc, &
194	advec_v_ws, advec_v_ws_acc, advec_w_ws, advec_w_ws_acc
195
196	USE advec_s_pw_mod, &
197	ONLY: advec_s_pw
198
199	USE advec_s_up_mod, &
200	ONLY: advec_s_up
201
202	USE advec_u_pw_mod, &
203	ONLY: advec_u_pw
204
205	USE advec_u_up_mod, &
206	ONLY: advec_u_up
207
208	USE advec_v_pw_mod, &
209	ONLY: advec_v_pw
210
211	USE advec_v_up_mod, &
212	ONLY: advec_v_up
213
214	USE advec_w_pw_mod, &
215	ONLY: advec_w_pw
216
217	USE advec_w_up_mod, &
218	ONLY: advec_w_up
219
220	USE buoyancy_mod, &
221	ONLY: buoyancy, buoyancy_acc
222
223	USE calc_precipitation_mod, &
224	ONLY: calc_precipitation
225
226	USE calc_radiation_mod, &
227	ONLY: calc_radiation
228
229	USE coriolis_mod, &
230	ONLY: coriolis, coriolis_acc
231
232	USE diffusion_e_mod, &
233	ONLY: diffusion_e, diffusion_e_acc
234
235	USE diffusion_s_mod, &
236	ONLY: diffusion_s, diffusion_s_acc
237
238	USE diffusion_u_mod, &
239	ONLY: diffusion_u, diffusion_u_acc
240
241	USE diffusion_v_mod, &
242	ONLY: diffusion_v, diffusion_v_acc
243
244	USE diffusion_w_mod, &
245	ONLY: diffusion_w, diffusion_w_acc
246
247	USE impact_of_latent_heat_mod, &
248	ONLY: impact_of_latent_heat
249
250	USE kinds
251
252	USE ls_forcing_mod, &
253	ONLY: ls_advec
254
255	USE microphysics_mod, &
256	ONLY: microphysics_control
257
258	USE nudge_mod, &
259	ONLY: nudge
260
261	USE plant_canopy_model_mod, &
262	ONLY: cthf, plant_canopy, plant_canopy_model
263
264	USE production_e_mod, &
265	ONLY: production_e, production_e_acc
266
267	USE statistics, &
268	ONLY: hom
269
270	USE subsidence_mod, &
271	ONLY: subsidence
272
273	USE user_actions_mod, &
274	ONLY: user_actions
275
276
277	PRIVATE
278	PUBLIC prognostic_equations_cache, prognostic_equations_vector, &
279	prognostic_equations_acc
280
281	INTERFACE prognostic_equations_cache
282	MODULE PROCEDURE prognostic_equations_cache
283	END INTERFACE prognostic_equations_cache
284
285	INTERFACE prognostic_equations_vector
286	MODULE PROCEDURE prognostic_equations_vector
287	END INTERFACE prognostic_equations_vector
288
289	INTERFACE prognostic_equations_acc
290	MODULE PROCEDURE prognostic_equations_acc
291	END INTERFACE prognostic_equations_acc
292
293
294	CONTAINS
295
296
297	SUBROUTINE prognostic_equations_cache
298
299	!------------------------------------------------------------------------------!
300	! Version with one optimized loop over all equations. It is only allowed to
301	! be called for the Wicker and Skamarock or Piascek-Williams advection scheme.
302	!
303	! Here the calls of most subroutines are embedded in two DO loops over i and j,
304	! so communication between CPUs is not allowed (does not make sense) within
305	! these loops.
306	!
307	! (Optimized to avoid cache missings, i.e. for Power4/5-architectures.)
308	!------------------------------------------------------------------------------!
309
310	IMPLICIT NONE
311
312	INTEGER(iwp) :: i !:
313	INTEGER(iwp) :: i_omp_start !:
314	INTEGER(iwp) :: j !:
315	INTEGER(iwp) :: k !:
316	INTEGER(iwp) :: omp_get_thread_num !:
317	INTEGER(iwp) :: tn = 0 !:
318
319	LOGICAL :: loop_start !:
320
321
322	!
323	!-- Time measurement can only be performed for the whole set of equations
324	CALL cpu_log( log_point(32), 'all progn.equations', 'start' )
325
326	!
327	!-- Loop over all prognostic equations
328	!$OMP PARALLEL private (i,i_omp_start,j,k,loop_start,tn)
329
330	!$ tn = omp_get_thread_num()
331	loop_start = .TRUE.
332	!$OMP DO
333	DO i = nxl, nxr
334
335	!
336	!-- Store the first loop index. It differs for each thread and is required
337	!-- later in advec_ws
338	IF ( loop_start ) THEN
339	loop_start = .FALSE.
340	i_omp_start = i
341	ENDIF
342
343	DO j = nys, nyn
344	!
345	!-- If required, calculate cloud microphysical impacts (two-moment scheme)
346	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. &
347	( intermediate_timestep_count == 1 .OR. &
348	call_microphysics_at_all_substeps ) &
349	) THEN
350	CALL microphysics_control( i, j )
351	ENDIF
352	!
353	!-- Tendency terms for u-velocity component
354	IF ( .NOT. outflow_l .OR. i > nxl ) THEN
355
356	tend(:,j,i) = 0.0_wp
357	IF ( timestep_scheme(1:5) == 'runge' ) THEN
358	IF ( ws_scheme_mom ) THEN
359	CALL advec_u_ws( i, j, i_omp_start, tn )
360	ELSE
361	CALL advec_u_pw( i, j )
362	ENDIF
363	ELSE
364	CALL advec_u_up( i, j )
365	ENDIF
366	CALL diffusion_u( i, j )
367	CALL coriolis( i, j, 1 )
368	IF ( sloping_surface .AND. .NOT. neutral ) THEN
369	CALL buoyancy( i, j, pt, 1 )
370	ENDIF
371
372	!
373	!-- Drag by plant canopy
374	IF ( plant_canopy ) CALL plant_canopy_model( i, j, 1 )
375
376	!
377	!-- External pressure gradient
378	IF ( dp_external ) THEN
379	DO k = dp_level_ind_b+1, nzt
380	tend(k,j,i) = tend(k,j,i) - dpdxy(1) * dp_smooth_factor(k)
381	ENDDO
382	ENDIF
383
384	!
385	!-- Nudging
386	IF ( nudging ) CALL nudge( i, j, simulated_time, 'u' )
387
388	CALL user_actions( i, j, 'u-tendency' )
389	!
390	!-- Prognostic equation for u-velocity component
391	DO k = nzb_u_inner(j,i)+1, nzt
392	u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
393	tsc(3) * tu_m(k,j,i) ) &
394	- tsc(5) * rdf(k) * ( u(k,j,i) - u_init(k) )
395	ENDDO
396
397	!
398	!-- Calculate tendencies for the next Runge-Kutta step
399	IF ( timestep_scheme(1:5) == 'runge' ) THEN
400	IF ( intermediate_timestep_count == 1 ) THEN
401	DO k = nzb_u_inner(j,i)+1, nzt
402	tu_m(k,j,i) = tend(k,j,i)
403	ENDDO
404	ELSEIF ( intermediate_timestep_count < &
405	intermediate_timestep_count_max ) THEN
406	DO k = nzb_u_inner(j,i)+1, nzt
407	tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tu_m(k,j,i)
408	ENDDO
409	ENDIF
410	ENDIF
411
412	ENDIF
413
414	!
415	!-- Tendency terms for v-velocity component
416	IF ( .NOT. outflow_s .OR. j > nys ) THEN
417
418	tend(:,j,i) = 0.0_wp
419	IF ( timestep_scheme(1:5) == 'runge' ) THEN
420	IF ( ws_scheme_mom ) THEN
421	CALL advec_v_ws( i, j, i_omp_start, tn )
422	ELSE
423	CALL advec_v_pw( i, j )
424	ENDIF
425	ELSE
426	CALL advec_v_up( i, j )
427	ENDIF
428	CALL diffusion_v( i, j )
429	CALL coriolis( i, j, 2 )
430
431	!
432	!-- Drag by plant canopy
433	IF ( plant_canopy ) CALL plant_canopy_model( i, j, 2 )
434
435	!
436	!-- External pressure gradient
437	IF ( dp_external ) THEN
438	DO k = dp_level_ind_b+1, nzt
439	tend(k,j,i) = tend(k,j,i) - dpdxy(2) * dp_smooth_factor(k)
440	ENDDO
441	ENDIF
442
443	!
444	!-- Nudging
445	IF ( nudging ) CALL nudge( i, j, simulated_time, 'v' )
446
447	CALL user_actions( i, j, 'v-tendency' )
448	!
449	!-- Prognostic equation for v-velocity component
450	DO k = nzb_v_inner(j,i)+1, nzt
451	v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
452	tsc(3) * tv_m(k,j,i) ) &
453	- tsc(5) * rdf(k) * ( v(k,j,i) - v_init(k) )
454	ENDDO
455
456	!
457	!-- Calculate tendencies for the next Runge-Kutta step
458	IF ( timestep_scheme(1:5) == 'runge' ) THEN
459	IF ( intermediate_timestep_count == 1 ) THEN
460	DO k = nzb_v_inner(j,i)+1, nzt
461	tv_m(k,j,i) = tend(k,j,i)
462	ENDDO
463	ELSEIF ( intermediate_timestep_count < &
464	intermediate_timestep_count_max ) THEN
465	DO k = nzb_v_inner(j,i)+1, nzt
466	tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tv_m(k,j,i)
467	ENDDO
468	ENDIF
469	ENDIF
470
471	ENDIF
472
473	!
474	!-- Tendency terms for w-velocity component
475	tend(:,j,i) = 0.0_wp
476	IF ( timestep_scheme(1:5) == 'runge' ) THEN
477	IF ( ws_scheme_mom ) THEN
478	CALL advec_w_ws( i, j, i_omp_start, tn )
479	ELSE
480	CALL advec_w_pw( i, j )
481	END IF
482	ELSE
483	CALL advec_w_up( i, j )
484	ENDIF
485	CALL diffusion_w( i, j )
486	CALL coriolis( i, j, 3 )
487
488	IF ( .NOT. neutral ) THEN
489	IF ( ocean ) THEN
490	CALL buoyancy( i, j, rho, 3 )
491	ELSE
492	IF ( .NOT. humidity ) THEN
493	CALL buoyancy( i, j, pt, 3 )
494	ELSE
495	CALL buoyancy( i, j, vpt, 3 )
496	ENDIF
497	ENDIF
498	ENDIF
499
500	!
501	!-- Drag by plant canopy
502	IF ( plant_canopy ) CALL plant_canopy_model( i, j, 3 )
503
504	CALL user_actions( i, j, 'w-tendency' )
505
506	!
507	!-- Prognostic equation for w-velocity component
508	DO k = nzb_w_inner(j,i)+1, nzt-1
509	w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
510	tsc(3) * tw_m(k,j,i) ) &
511	- tsc(5) * rdf(k) * w(k,j,i)
512	ENDDO
513
514	!
515	!-- Calculate tendencies for the next Runge-Kutta step
516	IF ( timestep_scheme(1:5) == 'runge' ) THEN
517	IF ( intermediate_timestep_count == 1 ) THEN
518	DO k = nzb_w_inner(j,i)+1, nzt-1
519	tw_m(k,j,i) = tend(k,j,i)
520	ENDDO
521	ELSEIF ( intermediate_timestep_count < &
522	intermediate_timestep_count_max ) THEN
523	DO k = nzb_w_inner(j,i)+1, nzt-1
524	tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tw_m(k,j,i)
525	ENDDO
526	ENDIF
527	ENDIF
528
529	!
530	!-- If required, compute prognostic equation for potential temperature
531	IF ( .NOT. neutral ) THEN
532	!
533	!-- Tendency terms for potential temperature
534	tend(:,j,i) = 0.0_wp
535	IF ( timestep_scheme(1:5) == 'runge' ) THEN
536	IF ( ws_scheme_sca ) THEN
537	CALL advec_s_ws( i, j, pt, 'pt', flux_s_pt, diss_s_pt, &
538	flux_l_pt, diss_l_pt, i_omp_start, tn )
539	ELSE
540	CALL advec_s_pw( i, j, pt )
541	ENDIF
542	ELSE
543	CALL advec_s_up( i, j, pt )
544	ENDIF
545	CALL diffusion_s( i, j, pt, shf, tswst, wall_heatflux )
546
547	!
548	!-- If required compute heating/cooling due to long wave radiation
549	!-- processes
550	IF ( radiation ) THEN
551	CALL calc_radiation( i, j )
552	ENDIF
553
554	!
555	!-- If required compute impact of latent heat due to precipitation
556	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. &
557	precipitation ) THEN
558	CALL impact_of_latent_heat( i, j )
559	ENDIF
560
561	!
562	!-- Consideration of heat sources within the plant canopy
563	IF ( plant_canopy .AND. cthf /= 0.0_wp ) THEN
564	CALL plant_canopy_model( i, j, 4 )
565	ENDIF
566
567	!
568	!-- Large scale advection
569	IF ( large_scale_forcing ) THEN
570	CALL ls_advec( i, j, simulated_time, 'pt' )
571	ENDIF
572
573	!
574	!-- Nudging
575	IF ( nudging ) CALL nudge( i, j, simulated_time, 'pt' )
576
577	!
578	!-- If required, compute effect of large-scale subsidence/ascent
579	IF ( large_scale_subsidence .AND. &
580	.NOT. use_subsidence_tendencies ) THEN
581	CALL subsidence( i, j, tend, pt, pt_init, 2 )
582	ENDIF
583
584	CALL user_actions( i, j, 'pt-tendency' )
585
586	!
587	!-- Prognostic equation for potential temperature
588	DO k = nzb_s_inner(j,i)+1, nzt
589	pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
590	tsc(3) * tpt_m(k,j,i) ) &
591	- tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
592	( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
593	ENDDO
594
595	!
596	!-- Calculate tendencies for the next Runge-Kutta step
597	IF ( timestep_scheme(1:5) == 'runge' ) THEN
598	IF ( intermediate_timestep_count == 1 ) THEN
599	DO k = nzb_s_inner(j,i)+1, nzt
600	tpt_m(k,j,i) = tend(k,j,i)
601	ENDDO
602	ELSEIF ( intermediate_timestep_count < &
603	intermediate_timestep_count_max ) THEN
604	DO k = nzb_s_inner(j,i)+1, nzt
605	tpt_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
606	5.3125_wp * tpt_m(k,j,i)
607	ENDDO
608	ENDIF
609	ENDIF
610
611	ENDIF
612
613	!
614	!-- If required, compute prognostic equation for salinity
615	IF ( ocean ) THEN
616
617	!
618	!-- Tendency-terms for salinity
619	tend(:,j,i) = 0.0_wp
620	IF ( timestep_scheme(1:5) == 'runge' ) &
621	THEN
622	IF ( ws_scheme_sca ) THEN
623	CALL advec_s_ws( i, j, sa, 'sa', flux_s_sa, &
624	diss_s_sa, flux_l_sa, diss_l_sa, i_omp_start, tn )
625	ELSE
626	CALL advec_s_pw( i, j, sa )
627	ENDIF
628	ELSE
629	CALL advec_s_up( i, j, sa )
630	ENDIF
631	CALL diffusion_s( i, j, sa, saswsb, saswst, wall_salinityflux )
632
633	CALL user_actions( i, j, 'sa-tendency' )
634
635	!
636	!-- Prognostic equation for salinity
637	DO k = nzb_s_inner(j,i)+1, nzt
638	sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
639	tsc(3) * tsa_m(k,j,i) ) &
640	- tsc(5) * rdf_sc(k) * &
641	( sa(k,j,i) - sa_init(k) )
642	IF ( sa_p(k,j,i) < 0.0_wp ) sa_p(k,j,i) = 0.1_wp * sa(k,j,i)
643	ENDDO
644
645	!
646	!-- Calculate tendencies for the next Runge-Kutta step
647	IF ( timestep_scheme(1:5) == 'runge' ) THEN
648	IF ( intermediate_timestep_count == 1 ) THEN
649	DO k = nzb_s_inner(j,i)+1, nzt
650	tsa_m(k,j,i) = tend(k,j,i)
651	ENDDO
652	ELSEIF ( intermediate_timestep_count < &
653	intermediate_timestep_count_max ) THEN
654	DO k = nzb_s_inner(j,i)+1, nzt
655	tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
656	5.3125_wp * tsa_m(k,j,i)
657	ENDDO
658	ENDIF
659	ENDIF
660
661	!
662	!-- Calculate density by the equation of state for seawater
663	CALL eqn_state_seawater( i, j )
664
665	ENDIF
666
667	!
668	!-- If required, compute prognostic equation for total water content /
669	!-- scalar
670	IF ( humidity .OR. passive_scalar ) THEN
671
672	!
673	!-- Tendency-terms for total water content / scalar
674	tend(:,j,i) = 0.0_wp
675	IF ( timestep_scheme(1:5) == 'runge' ) &
676	THEN
677	IF ( ws_scheme_sca ) THEN
678	CALL advec_s_ws( i, j, q, 'q', flux_s_q, &
679	diss_s_q, flux_l_q, diss_l_q, i_omp_start, tn )
680	ELSE
681	CALL advec_s_pw( i, j, q )
682	ENDIF
683	ELSE
684	CALL advec_s_up( i, j, q )
685	ENDIF
686	CALL diffusion_s( i, j, q, qsws, qswst, wall_qflux )
687
688	!
689	!-- If required compute decrease of total water content due to
690	!-- precipitation
691	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. &
692	precipitation ) THEN
693	CALL calc_precipitation( i, j )
694	ENDIF
695	!
696	!-- Sink or source of scalar concentration due to canopy elements
697	IF ( plant_canopy ) CALL plant_canopy_model( i, j, 5 )
698
699	!
700	!-- Large scale advection
701	IF ( large_scale_forcing ) THEN
702	CALL ls_advec( i, j, simulated_time, 'q' )
703	ENDIF
704
705	!
706	!-- Nudging
707	IF ( nudging ) CALL nudge( i, j, simulated_time, 'q' )
708
709	!
710	!-- If required compute influence of large-scale subsidence/ascent
711	IF ( large_scale_subsidence .AND. &
712	.NOT. use_subsidence_tendencies ) THEN
713	CALL subsidence( i, j, tend, q, q_init, 3 )
714	ENDIF
715
716	CALL user_actions( i, j, 'q-tendency' )
717
718	!
719	!-- Prognostic equation for total water content / scalar
720	DO k = nzb_s_inner(j,i)+1, nzt
721	q_p(k,j,i) = q(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
722	tsc(3) * tq_m(k,j,i) ) &
723	- tsc(5) * rdf_sc(k) * &
724	( q(k,j,i) - q_init(k) )
725	IF ( q_p(k,j,i) < 0.0_wp ) q_p(k,j,i) = 0.1_wp * q(k,j,i)
726	ENDDO
727
728	!
729	!-- Calculate tendencies for the next Runge-Kutta step
730	IF ( timestep_scheme(1:5) == 'runge' ) THEN
731	IF ( intermediate_timestep_count == 1 ) THEN
732	DO k = nzb_s_inner(j,i)+1, nzt
733	tq_m(k,j,i) = tend(k,j,i)
734	ENDDO
735	ELSEIF ( intermediate_timestep_count < &
736	intermediate_timestep_count_max ) THEN
737	DO k = nzb_s_inner(j,i)+1, nzt
738	tq_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
739	5.3125_wp * tq_m(k,j,i)
740	ENDDO
741	ENDIF
742	ENDIF
743
744	!
745	!-- If required, calculate prognostic equations for rain water content
746	!-- and rain drop concentration
747	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. &
748	precipitation ) THEN
749	!
750	!-- Calculate prognostic equation for rain water content
751	tend(:,j,i) = 0.0_wp
752	IF ( timestep_scheme(1:5) == 'runge' ) &
753	THEN
754	IF ( ws_scheme_sca ) THEN
755	CALL advec_s_ws( i, j, qr, 'qr', flux_s_qr, &
756	diss_s_qr, flux_l_qr, diss_l_qr, &
757	i_omp_start, tn )
758	ELSE
759	CALL advec_s_pw( i, j, qr )
760	ENDIF
761	ELSE
762	CALL advec_s_up( i, j, qr )
763	ENDIF
764	CALL diffusion_s( i, j, qr, qrsws, qrswst, wall_qrflux )
765
766	CALL user_actions( i, j, 'qr-tendency' )
767	!
768	!-- Prognostic equation for rain water content
769	DO k = nzb_s_inner(j,i)+1, nzt
770	qr_p(k,j,i) = qr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
771	tsc(3) * tqr_m(k,j,i) ) &
772	- tsc(5) * rdf_sc(k) * qr(k,j,i)
773	IF ( qr_p(k,j,i) < 0.0_wp ) qr_p(k,j,i) = 0.0_wp
774	ENDDO
775	!
776	!-- Calculate tendencies for the next Runge-Kutta step
777	IF ( timestep_scheme(1:5) == 'runge' ) THEN
778	IF ( intermediate_timestep_count == 1 ) THEN
779	DO k = nzb_s_inner(j,i)+1, nzt
780	tqr_m(k,j,i) = tend(k,j,i)
781	ENDDO
782	ELSEIF ( intermediate_timestep_count < &
783	intermediate_timestep_count_max ) THEN
784	DO k = nzb_s_inner(j,i)+1, nzt
785	tqr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
786	5.3125_wp * tqr_m(k,j,i)
787	ENDDO
788	ENDIF
789	ENDIF
790
791	!
792	!-- Calculate prognostic equation for rain drop concentration.
793	tend(:,j,i) = 0.0_wp
794	IF ( timestep_scheme(1:5) == 'runge' ) THEN
795	IF ( ws_scheme_sca ) THEN
796	CALL advec_s_ws( i, j, nr, 'nr', flux_s_nr, &
797	diss_s_nr, flux_l_nr, diss_l_nr, &
798	i_omp_start, tn )
799	ELSE
800	CALL advec_s_pw( i, j, nr )
801	ENDIF
802	ELSE
803	CALL advec_s_up( i, j, nr )
804	ENDIF
805	CALL diffusion_s( i, j, nr, nrsws, nrswst, wall_nrflux )
806
807	CALL user_actions( i, j, 'nr-tendency' )
808	!
809	!-- Prognostic equation for rain drop concentration
810	DO k = nzb_s_inner(j,i)+1, nzt
811	nr_p(k,j,i) = nr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
812	tsc(3) * tnr_m(k,j,i) ) &
813	- tsc(5) * rdf_sc(k) * nr(k,j,i)
814	IF ( nr_p(k,j,i) < 0.0_wp ) nr_p(k,j,i) = 0.0_wp
815	ENDDO
816	!
817	!-- Calculate tendencies for the next Runge-Kutta step
818	IF ( timestep_scheme(1:5) == 'runge' ) THEN
819	IF ( intermediate_timestep_count == 1 ) THEN
820	DO k = nzb_s_inner(j,i)+1, nzt
821	tnr_m(k,j,i) = tend(k,j,i)
822	ENDDO
823	ELSEIF ( intermediate_timestep_count < &
824	intermediate_timestep_count_max ) THEN
825	DO k = nzb_s_inner(j,i)+1, nzt
826	tnr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
827	5.3125_wp * tnr_m(k,j,i)
828	ENDDO
829	ENDIF
830	ENDIF
831
832	ENDIF
833
834	ENDIF
835
836	!
837	!-- If required, compute prognostic equation for turbulent kinetic
838	!-- energy (TKE)
839	IF ( .NOT. constant_diffusion ) THEN
840
841	!
842	!-- Tendency-terms for TKE
843	tend(:,j,i) = 0.0_wp
844	IF ( timestep_scheme(1:5) == 'runge' &
845	.AND. .NOT. use_upstream_for_tke ) THEN
846	IF ( ws_scheme_sca ) THEN
847	CALL advec_s_ws( i, j, e, 'e', flux_s_e, diss_s_e, &
848	flux_l_e, diss_l_e , i_omp_start, tn )
849	ELSE
850	CALL advec_s_pw( i, j, e )
851	ENDIF
852	ELSE
853	CALL advec_s_up( i, j, e )
854	ENDIF
855	IF ( .NOT. humidity ) THEN
856	IF ( ocean ) THEN
857	CALL diffusion_e( i, j, prho, prho_reference )
858	ELSE
859	CALL diffusion_e( i, j, pt, pt_reference )
860	ENDIF
861	ELSE
862	CALL diffusion_e( i, j, vpt, pt_reference )
863	ENDIF
864	CALL production_e( i, j )
865
866	!
867	!-- Additional sink term for flows through plant canopies
868	IF ( plant_canopy ) CALL plant_canopy_model( i, j, 6 )
869
870	CALL user_actions( i, j, 'e-tendency' )
871
872	!
873	!-- Prognostic equation for TKE.
874	!-- Eliminate negative TKE values, which can occur due to numerical
875	!-- reasons in the course of the integration. In such cases the old
876	!-- TKE value is reduced by 90%.
877	DO k = nzb_s_inner(j,i)+1, nzt
878	e_p(k,j,i) = e(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
879	tsc(3) * te_m(k,j,i) )
880	IF ( e_p(k,j,i) < 0.0_wp ) e_p(k,j,i) = 0.1_wp * e(k,j,i)
881	ENDDO
882
883	!
884	!-- Calculate tendencies for the next Runge-Kutta step
885	IF ( timestep_scheme(1:5) == 'runge' ) THEN
886	IF ( intermediate_timestep_count == 1 ) THEN
887	DO k = nzb_s_inner(j,i)+1, nzt
888	te_m(k,j,i) = tend(k,j,i)
889	ENDDO
890	ELSEIF ( intermediate_timestep_count < &
891	intermediate_timestep_count_max ) THEN
892	DO k = nzb_s_inner(j,i)+1, nzt
893	te_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
894	5.3125_wp * te_m(k,j,i)
895	ENDDO
896	ENDIF
897	ENDIF
898
899	ENDIF ! TKE equation
900
901	ENDDO
902	ENDDO
903	!$OMP END PARALLEL
904
905	CALL cpu_log( log_point(32), 'all progn.equations', 'stop' )
906
907
908	END SUBROUTINE prognostic_equations_cache
909
910
911	SUBROUTINE prognostic_equations_vector
912
913	!------------------------------------------------------------------------------!
914	! Version for vector machines
915	!------------------------------------------------------------------------------!
916
917	IMPLICIT NONE
918
919	INTEGER(iwp) :: i !:
920	INTEGER(iwp) :: j !:
921	INTEGER(iwp) :: k !:
922
923	REAL(wp) :: sbt !:
924
925
926	!
927	!-- If required, calculate cloud microphysical impacts (two-moment scheme)
928	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. &
929	( intermediate_timestep_count == 1 .OR. &
930	call_microphysics_at_all_substeps ) &
931	) THEN
932	CALL cpu_log( log_point(51), 'microphysics', 'start' )
933	CALL microphysics_control
934	CALL cpu_log( log_point(51), 'microphysics', 'stop' )
935	ENDIF
936
937	!
938	!-- u-velocity component
939	CALL cpu_log( log_point(5), 'u-equation', 'start' )
940
941	tend = 0.0_wp
942	IF ( timestep_scheme(1:5) == 'runge' ) THEN
943	IF ( ws_scheme_mom ) THEN
944	CALL advec_u_ws
945	ELSE
946	CALL advec_u_pw
947	ENDIF
948	ELSE
949	CALL advec_u_up
950	ENDIF
951	CALL diffusion_u
952	CALL coriolis( 1 )
953	IF ( sloping_surface .AND. .NOT. neutral ) THEN
954	CALL buoyancy( pt, 1 )
955	ENDIF
956
957	!
958	!-- Drag by plant canopy
959	IF ( plant_canopy ) CALL plant_canopy_model( 1 )
960
961	!
962	!-- External pressure gradient
963	IF ( dp_external ) THEN
964	DO i = nxlu, nxr
965	DO j = nys, nyn
966	DO k = dp_level_ind_b+1, nzt
967	tend(k,j,i) = tend(k,j,i) - dpdxy(1) * dp_smooth_factor(k)
968	ENDDO
969	ENDDO
970	ENDDO
971	ENDIF
972
973	!
974	!-- Nudging
975	IF ( nudging ) CALL nudge( simulated_time, 'u' )
976
977	CALL user_actions( 'u-tendency' )
978
979	!
980	!-- Prognostic equation for u-velocity component
981	DO i = nxlu, nxr
982	DO j = nys, nyn
983	DO k = nzb_u_inner(j,i)+1, nzt
984	u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
985	tsc(3) * tu_m(k,j,i) ) &
986	- tsc(5) * rdf(k) * ( u(k,j,i) - u_init(k) )
987	ENDDO
988	ENDDO
989	ENDDO
990
991	!
992	!-- Calculate tendencies for the next Runge-Kutta step
993	IF ( timestep_scheme(1:5) == 'runge' ) THEN
994	IF ( intermediate_timestep_count == 1 ) THEN
995	DO i = nxlu, nxr
996	DO j = nys, nyn
997	DO k = nzb_u_inner(j,i)+1, nzt
998	tu_m(k,j,i) = tend(k,j,i)
999	ENDDO
1000	ENDDO
1001	ENDDO
1002	ELSEIF ( intermediate_timestep_count < &
1003	intermediate_timestep_count_max ) THEN
1004	DO i = nxlu, nxr
1005	DO j = nys, nyn
1006	DO k = nzb_u_inner(j,i)+1, nzt
1007	tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tu_m(k,j,i)
1008	ENDDO
1009	ENDDO
1010	ENDDO
1011	ENDIF
1012	ENDIF
1013
1014	CALL cpu_log( log_point(5), 'u-equation', 'stop' )
1015
1016	!
1017	!-- v-velocity component
1018	CALL cpu_log( log_point(6), 'v-equation', 'start' )
1019
1020	tend = 0.0_wp
1021	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1022	IF ( ws_scheme_mom ) THEN
1023	CALL advec_v_ws
1024	ELSE
1025	CALL advec_v_pw
1026	END IF
1027	ELSE
1028	CALL advec_v_up
1029	ENDIF
1030	CALL diffusion_v
1031	CALL coriolis( 2 )
1032
1033	!
1034	!-- Drag by plant canopy
1035	IF ( plant_canopy ) CALL plant_canopy_model( 2 )
1036
1037	!
1038	!-- External pressure gradient
1039	IF ( dp_external ) THEN
1040	DO i = nxl, nxr
1041	DO j = nysv, nyn
1042	DO k = dp_level_ind_b+1, nzt
1043	tend(k,j,i) = tend(k,j,i) - dpdxy(2) * dp_smooth_factor(k)
1044	ENDDO
1045	ENDDO
1046	ENDDO
1047	ENDIF
1048
1049	!
1050	!-- Nudging
1051	IF ( nudging ) CALL nudge( simulated_time, 'v' )
1052
1053	CALL user_actions( 'v-tendency' )
1054
1055	!
1056	!-- Prognostic equation for v-velocity component
1057	DO i = nxl, nxr
1058	DO j = nysv, nyn
1059	DO k = nzb_v_inner(j,i)+1, nzt
1060	v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
1061	tsc(3) * tv_m(k,j,i) ) &
1062	- tsc(5) * rdf(k) * ( v(k,j,i) - v_init(k) )
1063	ENDDO
1064	ENDDO
1065	ENDDO
1066
1067	!
1068	!-- Calculate tendencies for the next Runge-Kutta step
1069	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1070	IF ( intermediate_timestep_count == 1 ) THEN
1071	DO i = nxl, nxr
1072	DO j = nysv, nyn
1073	DO k = nzb_v_inner(j,i)+1, nzt
1074	tv_m(k,j,i) = tend(k,j,i)
1075	ENDDO
1076	ENDDO
1077	ENDDO
1078	ELSEIF ( intermediate_timestep_count < &
1079	intermediate_timestep_count_max ) THEN
1080	DO i = nxl, nxr
1081	DO j = nysv, nyn
1082	DO k = nzb_v_inner(j,i)+1, nzt
1083	tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tv_m(k,j,i)
1084	ENDDO
1085	ENDDO
1086	ENDDO
1087	ENDIF
1088	ENDIF
1089
1090	CALL cpu_log( log_point(6), 'v-equation', 'stop' )
1091
1092	!
1093	!-- w-velocity component
1094	CALL cpu_log( log_point(7), 'w-equation', 'start' )
1095
1096	tend = 0.0_wp
1097	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1098	IF ( ws_scheme_mom ) THEN
1099	CALL advec_w_ws
1100	ELSE
1101	CALL advec_w_pw
1102	ENDIF
1103	ELSE
1104	CALL advec_w_up
1105	ENDIF
1106	CALL diffusion_w
1107	CALL coriolis( 3 )
1108
1109	IF ( .NOT. neutral ) THEN
1110	IF ( ocean ) THEN
1111	CALL buoyancy( rho, 3 )
1112	ELSE
1113	IF ( .NOT. humidity ) THEN
1114	CALL buoyancy( pt, 3 )
1115	ELSE
1116	CALL buoyancy( vpt, 3 )
1117	ENDIF
1118	ENDIF
1119	ENDIF
1120
1121	!
1122	!-- Drag by plant canopy
1123	IF ( plant_canopy ) CALL plant_canopy_model( 3 )
1124
1125	CALL user_actions( 'w-tendency' )
1126
1127	!
1128	!-- Prognostic equation for w-velocity component
1129	DO i = nxl, nxr
1130	DO j = nys, nyn
1131	DO k = nzb_w_inner(j,i)+1, nzt-1
1132	w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
1133	tsc(3) * tw_m(k,j,i) ) &
1134	- tsc(5) * rdf(k) * w(k,j,i)
1135	ENDDO
1136	ENDDO
1137	ENDDO
1138
1139	!
1140	!-- Calculate tendencies for the next Runge-Kutta step
1141	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1142	IF ( intermediate_timestep_count == 1 ) THEN
1143	DO i = nxl, nxr
1144	DO j = nys, nyn
1145	DO k = nzb_w_inner(j,i)+1, nzt-1
1146	tw_m(k,j,i) = tend(k,j,i)
1147	ENDDO
1148	ENDDO
1149	ENDDO
1150	ELSEIF ( intermediate_timestep_count < &
1151	intermediate_timestep_count_max ) THEN
1152	DO i = nxl, nxr
1153	DO j = nys, nyn
1154	DO k = nzb_w_inner(j,i)+1, nzt-1
1155	tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tw_m(k,j,i)
1156	ENDDO
1157	ENDDO
1158	ENDDO
1159	ENDIF
1160	ENDIF
1161
1162	CALL cpu_log( log_point(7), 'w-equation', 'stop' )
1163
1164
1165	!
1166	!-- If required, compute prognostic equation for potential temperature
1167	IF ( .NOT. neutral ) THEN
1168
1169	CALL cpu_log( log_point(13), 'pt-equation', 'start' )
1170
1171	!
1172	!-- pt-tendency terms with communication
1173	sbt = tsc(2)
1174	IF ( scalar_advec == 'bc-scheme' ) THEN
1175
1176	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1177	!
1178	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1179	sbt = 1.0_wp
1180	ENDIF
1181	tend = 0.0_wp
1182	CALL advec_s_bc( pt, 'pt' )
1183
1184	ENDIF
1185
1186	!
1187	!-- pt-tendency terms with no communication
1188	IF ( scalar_advec /= 'bc-scheme' ) THEN
1189	tend = 0.0_wp
1190	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1191	IF ( ws_scheme_sca ) THEN
1192	CALL advec_s_ws( pt, 'pt' )
1193	ELSE
1194	CALL advec_s_pw( pt )
1195	ENDIF
1196	ELSE
1197	CALL advec_s_up( pt )
1198	ENDIF
1199	ENDIF
1200
1201	CALL diffusion_s( pt, shf, tswst, wall_heatflux )
1202
1203	!
1204	!-- If required compute heating/cooling due to long wave radiation processes
1205	IF ( radiation ) THEN
1206	CALL calc_radiation
1207	ENDIF
1208
1209	!
1210	!-- If required compute impact of latent heat due to precipitation
1211	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. precipitation ) THEN
1212	CALL impact_of_latent_heat
1213	ENDIF
1214
1215	!
1216	!-- Consideration of heat sources within the plant canopy
1217	IF ( plant_canopy .AND. ( cthf /= 0.0_wp ) ) THEN
1218	CALL plant_canopy_model( 4 )
1219	ENDIF
1220
1221	!
1222	!-- Large scale advection
1223	IF ( large_scale_forcing ) THEN
1224	CALL ls_advec( simulated_time, 'pt' )
1225	ENDIF
1226
1227	!
1228	!-- Nudging
1229	IF ( nudging ) CALL nudge( simulated_time, 'pt' )
1230
1231	!
1232	!-- If required compute influence of large-scale subsidence/ascent
1233	IF ( large_scale_subsidence .AND. &
1234	.NOT. use_subsidence_tendencies ) THEN
1235	CALL subsidence( tend, pt, pt_init, 2 )
1236	ENDIF
1237
1238	CALL user_actions( 'pt-tendency' )
1239
1240	!
1241	!-- Prognostic equation for potential temperature
1242	DO i = nxl, nxr
1243	DO j = nys, nyn
1244	DO k = nzb_s_inner(j,i)+1, nzt
1245	pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1246	tsc(3) * tpt_m(k,j,i) ) &
1247	- tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
1248	( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
1249	ENDDO
1250	ENDDO
1251	ENDDO
1252
1253	!
1254	!-- Calculate tendencies for the next Runge-Kutta step
1255	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1256	IF ( intermediate_timestep_count == 1 ) THEN
1257	DO i = nxl, nxr
1258	DO j = nys, nyn
1259	DO k = nzb_s_inner(j,i)+1, nzt
1260	tpt_m(k,j,i) = tend(k,j,i)
1261	ENDDO
1262	ENDDO
1263	ENDDO
1264	ELSEIF ( intermediate_timestep_count < &
1265	intermediate_timestep_count_max ) THEN
1266	DO i = nxl, nxr
1267	DO j = nys, nyn
1268	DO k = nzb_s_inner(j,i)+1, nzt
1269	tpt_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
1270	5.3125_wp * tpt_m(k,j,i)
1271	ENDDO
1272	ENDDO
1273	ENDDO
1274	ENDIF
1275	ENDIF
1276
1277	CALL cpu_log( log_point(13), 'pt-equation', 'stop' )
1278
1279	ENDIF
1280
1281	!
1282	!-- If required, compute prognostic equation for salinity
1283	IF ( ocean ) THEN
1284
1285	CALL cpu_log( log_point(37), 'sa-equation', 'start' )
1286
1287	!
1288	!-- sa-tendency terms with communication
1289	sbt = tsc(2)
1290	IF ( scalar_advec == 'bc-scheme' ) THEN
1291
1292	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1293	!
1294	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1295	sbt = 1.0_wp
1296	ENDIF
1297	tend = 0.0_wp
1298	CALL advec_s_bc( sa, 'sa' )
1299
1300	ENDIF
1301
1302	!
1303	!-- sa-tendency terms with no communication
1304	IF ( scalar_advec /= 'bc-scheme' ) THEN
1305	tend = 0.0_wp
1306	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1307	IF ( ws_scheme_sca ) THEN
1308	CALL advec_s_ws( sa, 'sa' )
1309	ELSE
1310	CALL advec_s_pw( sa )
1311	ENDIF
1312	ELSE
1313	CALL advec_s_up( sa )
1314	ENDIF
1315	ENDIF
1316
1317	CALL diffusion_s( sa, saswsb, saswst, wall_salinityflux )
1318
1319	CALL user_actions( 'sa-tendency' )
1320
1321	!
1322	!-- Prognostic equation for salinity
1323	DO i = nxl, nxr
1324	DO j = nys, nyn
1325	DO k = nzb_s_inner(j,i)+1, nzt
1326	sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1327	tsc(3) * tsa_m(k,j,i) ) &
1328	- tsc(5) * rdf_sc(k) * &
1329	( sa(k,j,i) - sa_init(k) )
1330	IF ( sa_p(k,j,i) < 0.0_wp ) sa_p(k,j,i) = 0.1_wp * sa(k,j,i)
1331	ENDDO
1332	ENDDO
1333	ENDDO
1334
1335	!
1336	!-- Calculate tendencies for the next Runge-Kutta step
1337	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1338	IF ( intermediate_timestep_count == 1 ) THEN
1339	DO i = nxl, nxr
1340	DO j = nys, nyn
1341	DO k = nzb_s_inner(j,i)+1, nzt
1342	tsa_m(k,j,i) = tend(k,j,i)
1343	ENDDO
1344	ENDDO
1345	ENDDO
1346	ELSEIF ( intermediate_timestep_count < &
1347	intermediate_timestep_count_max ) THEN
1348	DO i = nxl, nxr
1349	DO j = nys, nyn
1350	DO k = nzb_s_inner(j,i)+1, nzt
1351	tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
1352	5.3125_wp * tsa_m(k,j,i)
1353	ENDDO
1354	ENDDO
1355	ENDDO
1356	ENDIF
1357	ENDIF
1358
1359	CALL cpu_log( log_point(37), 'sa-equation', 'stop' )
1360
1361	!
1362	!-- Calculate density by the equation of state for seawater
1363	CALL cpu_log( log_point(38), 'eqns-seawater', 'start' )
1364	CALL eqn_state_seawater
1365	CALL cpu_log( log_point(38), 'eqns-seawater', 'stop' )
1366
1367	ENDIF
1368
1369	!
1370	!-- If required, compute prognostic equation for total water content / scalar
1371	IF ( humidity .OR. passive_scalar ) THEN
1372
1373	CALL cpu_log( log_point(29), 'q/s-equation', 'start' )
1374
1375	!
1376	!-- Scalar/q-tendency terms with communication
1377	sbt = tsc(2)
1378	IF ( scalar_advec == 'bc-scheme' ) THEN
1379
1380	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1381	!
1382	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1383	sbt = 1.0_wp
1384	ENDIF
1385	tend = 0.0_wp
1386	CALL advec_s_bc( q, 'q' )
1387
1388	ENDIF
1389
1390	!
1391	!-- Scalar/q-tendency terms with no communication
1392	IF ( scalar_advec /= 'bc-scheme' ) THEN
1393	tend = 0.0_wp
1394	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1395	IF ( ws_scheme_sca ) THEN
1396	CALL advec_s_ws( q, 'q' )
1397	ELSE
1398	CALL advec_s_pw( q )
1399	ENDIF
1400	ELSE
1401	CALL advec_s_up( q )
1402	ENDIF
1403	ENDIF
1404
1405	CALL diffusion_s( q, qsws, qswst, wall_qflux )
1406
1407	!
1408	!-- If required compute decrease of total water content due to
1409	!-- precipitation
1410	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. precipitation ) THEN
1411	CALL calc_precipitation
1412	ENDIF
1413
1414	!
1415	!-- Sink or source of scalar concentration due to canopy elements
1416	IF ( plant_canopy ) CALL plant_canopy_model( 5 )
1417
1418	!
1419	!-- Large scale advection
1420	IF ( large_scale_forcing ) THEN
1421	CALL ls_advec( simulated_time, 'q' )
1422	ENDIF
1423
1424	!
1425	!-- Nudging
1426	IF ( nudging ) CALL nudge( simulated_time, 'q' )
1427
1428	!
1429	!-- If required compute influence of large-scale subsidence/ascent
1430	IF ( large_scale_subsidence .AND. &
1431	.NOT. use_subsidence_tendencies ) THEN
1432	CALL subsidence( tend, q, q_init, 3 )
1433	ENDIF
1434
1435	CALL user_actions( 'q-tendency' )
1436
1437	!
1438	!-- Prognostic equation for total water content / scalar
1439	DO i = nxl, nxr
1440	DO j = nys, nyn
1441	DO k = nzb_s_inner(j,i)+1, nzt
1442	q_p(k,j,i) = q(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1443	tsc(3) * tq_m(k,j,i) ) &
1444	- tsc(5) * rdf_sc(k) * &
1445	( q(k,j,i) - q_init(k) )
1446	IF ( q_p(k,j,i) < 0.0_wp ) q_p(k,j,i) = 0.1_wp * q(k,j,i)
1447	ENDDO
1448	ENDDO
1449	ENDDO
1450
1451	!
1452	!-- Calculate tendencies for the next Runge-Kutta step
1453	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1454	IF ( intermediate_timestep_count == 1 ) THEN
1455	DO i = nxl, nxr
1456	DO j = nys, nyn
1457	DO k = nzb_s_inner(j,i)+1, nzt
1458	tq_m(k,j,i) = tend(k,j,i)
1459	ENDDO
1460	ENDDO
1461	ENDDO
1462	ELSEIF ( intermediate_timestep_count < &
1463	intermediate_timestep_count_max ) THEN
1464	DO i = nxl, nxr
1465	DO j = nys, nyn
1466	DO k = nzb_s_inner(j,i)+1, nzt
1467	tq_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tq_m(k,j,i)
1468	ENDDO
1469	ENDDO
1470	ENDDO
1471	ENDIF
1472	ENDIF
1473
1474	CALL cpu_log( log_point(29), 'q/s-equation', 'stop' )
1475
1476	!
1477	!-- If required, calculate prognostic equations for rain water content
1478	!-- and rain drop concentration
1479	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. precipitation ) THEN
1480
1481	CALL cpu_log( log_point(52), 'qr-equation', 'start' )
1482
1483	!
1484	!-- Calculate prognostic equation for rain water content
1485	sbt = tsc(2)
1486	IF ( scalar_advec == 'bc-scheme' ) THEN
1487
1488	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1489	!
1490	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1491	sbt = 1.0_wp
1492	ENDIF
1493	tend = 0.0_wp
1494	CALL advec_s_bc( qr, 'qr' )
1495
1496	ENDIF
1497
1498	!
1499	!-- qr-tendency terms with no communication
1500	IF ( scalar_advec /= 'bc-scheme' ) THEN
1501	tend = 0.0_wp
1502	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1503	IF ( ws_scheme_sca ) THEN
1504	CALL advec_s_ws( qr, 'qr' )
1505	ELSE
1506	CALL advec_s_pw( qr )
1507	ENDIF
1508	ELSE
1509	CALL advec_s_up( qr )
1510	ENDIF
1511	ENDIF
1512
1513	CALL diffusion_s( qr, qrsws, qrswst, wall_qrflux )
1514
1515	CALL user_actions( 'qr-tendency' )
1516
1517	!
1518	!-- Prognostic equation for rain water content
1519	DO i = nxl, nxr
1520	DO j = nys, nyn
1521	DO k = nzb_s_inner(j,i)+1, nzt
1522	qr_p(k,j,i) = qr(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1523	tsc(3) * tqr_m(k,j,i) ) &
1524	- tsc(5) * rdf_sc(k) * qr(k,j,i)
1525	IF ( qr_p(k,j,i) < 0.0_wp ) qr_p(k,j,i) = 0.0_wp
1526	ENDDO
1527	ENDDO
1528	ENDDO
1529
1530	!
1531	!-- Calculate tendencies for the next Runge-Kutta step
1532	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1533	IF ( intermediate_timestep_count == 1 ) THEN
1534	DO i = nxl, nxr
1535	DO j = nys, nyn
1536	DO k = nzb_s_inner(j,i)+1, nzt
1537	tqr_m(k,j,i) = tend(k,j,i)
1538	ENDDO
1539	ENDDO
1540	ENDDO
1541	ELSEIF ( intermediate_timestep_count < &
1542	intermediate_timestep_count_max ) THEN
1543	DO i = nxl, nxr
1544	DO j = nys, nyn
1545	DO k = nzb_s_inner(j,i)+1, nzt
1546	tqr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * &
1547	tqr_m(k,j,i)
1548	ENDDO
1549	ENDDO
1550	ENDDO
1551	ENDIF
1552	ENDIF
1553
1554	CALL cpu_log( log_point(52), 'qr-equation', 'stop' )
1555	CALL cpu_log( log_point(53), 'nr-equation', 'start' )
1556
1557	!
1558	!-- Calculate prognostic equation for rain drop concentration
1559	sbt = tsc(2)
1560	IF ( scalar_advec == 'bc-scheme' ) THEN
1561
1562	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1563	!
1564	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1565	sbt = 1.0_wp
1566	ENDIF
1567	tend = 0.0_wp
1568	CALL advec_s_bc( nr, 'nr' )
1569
1570	ENDIF
1571
1572	!
1573	!-- nr-tendency terms with no communication
1574	IF ( scalar_advec /= 'bc-scheme' ) THEN
1575	tend = 0.0_wp
1576	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1577	IF ( ws_scheme_sca ) THEN
1578	CALL advec_s_ws( nr, 'nr' )
1579	ELSE
1580	CALL advec_s_pw( nr )
1581	ENDIF
1582	ELSE
1583	CALL advec_s_up( nr )
1584	ENDIF
1585	ENDIF
1586
1587	CALL diffusion_s( nr, nrsws, nrswst, wall_nrflux )
1588
1589	CALL user_actions( 'nr-tendency' )
1590
1591	!
1592	!-- Prognostic equation for rain drop concentration
1593	DO i = nxl, nxr
1594	DO j = nys, nyn
1595	DO k = nzb_s_inner(j,i)+1, nzt
1596	nr_p(k,j,i) = nr(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1597	tsc(3) * tnr_m(k,j,i) ) &
1598	- tsc(5) * rdf_sc(k) * nr(k,j,i)
1599	IF ( nr_p(k,j,i) < 0.0_wp ) nr_p(k,j,i) = 0.0_wp
1600	ENDDO
1601	ENDDO
1602	ENDDO
1603
1604	!
1605	!-- Calculate tendencies for the next Runge-Kutta step
1606	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1607	IF ( intermediate_timestep_count == 1 ) THEN
1608	DO i = nxl, nxr
1609	DO j = nys, nyn
1610	DO k = nzb_s_inner(j,i)+1, nzt
1611	tnr_m(k,j,i) = tend(k,j,i)
1612	ENDDO
1613	ENDDO
1614	ENDDO
1615	ELSEIF ( intermediate_timestep_count < &
1616	intermediate_timestep_count_max ) THEN
1617	DO i = nxl, nxr
1618	DO j = nys, nyn
1619	DO k = nzb_s_inner(j,i)+1, nzt
1620	tnr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * &
1621	tnr_m(k,j,i)
1622	ENDDO
1623	ENDDO
1624	ENDDO
1625	ENDIF
1626	ENDIF
1627
1628	CALL cpu_log( log_point(53), 'nr-equation', 'stop' )
1629
1630	ENDIF
1631
1632	ENDIF
1633
1634	!
1635	!-- If required, compute prognostic equation for turbulent kinetic
1636	!-- energy (TKE)
1637	IF ( .NOT. constant_diffusion ) THEN
1638
1639	CALL cpu_log( log_point(16), 'tke-equation', 'start' )
1640
1641	sbt = tsc(2)
1642	IF ( .NOT. use_upstream_for_tke ) THEN
1643	IF ( scalar_advec == 'bc-scheme' ) THEN
1644
1645	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1646	!
1647	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1648	sbt = 1.0_wp
1649	ENDIF
1650	tend = 0.0_wp
1651	CALL advec_s_bc( e, 'e' )
1652
1653	ENDIF
1654	ENDIF
1655
1656	!
1657	!-- TKE-tendency terms with no communication
1658	IF ( scalar_advec /= 'bc-scheme' .OR. use_upstream_for_tke ) THEN
1659	IF ( use_upstream_for_tke ) THEN
1660	tend = 0.0_wp
1661	CALL advec_s_up( e )
1662	ELSE
1663	tend = 0.0_wp
1664	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1665	IF ( ws_scheme_sca ) THEN
1666	CALL advec_s_ws( e, 'e' )
1667	ELSE
1668	CALL advec_s_pw( e )
1669	ENDIF
1670	ELSE
1671	CALL advec_s_up( e )
1672	ENDIF
1673	ENDIF
1674	ENDIF
1675
1676	IF ( .NOT. humidity ) THEN
1677	IF ( ocean ) THEN
1678	CALL diffusion_e( prho, prho_reference )
1679	ELSE
1680	CALL diffusion_e( pt, pt_reference )
1681	ENDIF
1682	ELSE
1683	CALL diffusion_e( vpt, pt_reference )
1684	ENDIF
1685
1686	CALL production_e
1687
1688	!
1689	!-- Additional sink term for flows through plant canopies
1690	IF ( plant_canopy ) CALL plant_canopy_model( 6 )
1691	CALL user_actions( 'e-tendency' )
1692
1693	!
1694	!-- Prognostic equation for TKE.
1695	!-- Eliminate negative TKE values, which can occur due to numerical
1696	!-- reasons in the course of the integration. In such cases the old TKE
1697	!-- value is reduced by 90%.
1698	DO i = nxl, nxr
1699	DO j = nys, nyn
1700	DO k = nzb_s_inner(j,i)+1, nzt
1701	e_p(k,j,i) = e(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
1702	tsc(3) * te_m(k,j,i) )
1703	IF ( e_p(k,j,i) < 0.0_wp ) e_p(k,j,i) = 0.1_wp * e(k,j,i)
1704	ENDDO
1705	ENDDO
1706	ENDDO
1707
1708	!
1709	!-- Calculate tendencies for the next Runge-Kutta step
1710	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1711	IF ( intermediate_timestep_count == 1 ) THEN
1712	DO i = nxl, nxr
1713	DO j = nys, nyn
1714	DO k = nzb_s_inner(j,i)+1, nzt
1715	te_m(k,j,i) = tend(k,j,i)
1716	ENDDO
1717	ENDDO
1718	ENDDO
1719	ELSEIF ( intermediate_timestep_count < &
1720	intermediate_timestep_count_max ) THEN
1721	DO i = nxl, nxr
1722	DO j = nys, nyn
1723	DO k = nzb_s_inner(j,i)+1, nzt
1724	te_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * te_m(k,j,i)
1725	ENDDO
1726	ENDDO
1727	ENDDO
1728	ENDIF
1729	ENDIF
1730
1731	CALL cpu_log( log_point(16), 'tke-equation', 'stop' )
1732
1733	ENDIF
1734
1735	END SUBROUTINE prognostic_equations_vector
1736
1737
1738	SUBROUTINE prognostic_equations_acc
1739
1740	!------------------------------------------------------------------------------!
1741	! Version for accelerator boards
1742	!------------------------------------------------------------------------------!
1743
1744	IMPLICIT NONE
1745
1746	INTEGER(iwp) :: i !:
1747	INTEGER(iwp) :: j !:
1748	INTEGER(iwp) :: k !:
1749	INTEGER(iwp) :: runge_step !:
1750
1751	REAL(wp) :: sbt !:
1752
1753	!
1754	!-- Set switch for intermediate Runge-Kutta step
1755	runge_step = 0
1756	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1757	IF ( intermediate_timestep_count == 1 ) THEN
1758	runge_step = 1
1759	ELSEIF ( intermediate_timestep_count < &
1760	intermediate_timestep_count_max ) THEN
1761	runge_step = 2
1762	ENDIF
1763	ENDIF
1764
1765	!
1766	!-- If required, calculate cloud microphysical impacts (two-moment scheme)
1767	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. &
1768	( intermediate_timestep_count == 1 .OR. &
1769	call_microphysics_at_all_substeps ) &
1770	) THEN
1771	CALL cpu_log( log_point(51), 'microphysics', 'start' )
1772	CALL microphysics_control
1773	CALL cpu_log( log_point(51), 'microphysics', 'stop' )
1774	ENDIF
1775
1776	!
1777	!-- u-velocity component
1778	!++ Statistics still not completely ported to accelerators
1779	!$acc update device( hom, ref_state )
1780	CALL cpu_log( log_point(5), 'u-equation', 'start' )
1781
1782	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1783	IF ( ws_scheme_mom ) THEN
1784	CALL advec_u_ws_acc
1785	ELSE
1786	tend = 0.0_wp ! to be removed later??
1787	CALL advec_u_pw
1788	ENDIF
1789	ELSE
1790	CALL advec_u_up
1791	ENDIF
1792	CALL diffusion_u_acc
1793	CALL coriolis_acc( 1 )
1794	IF ( sloping_surface .AND. .NOT. neutral ) THEN
1795	CALL buoyancy( pt, 1 )
1796	ENDIF
1797
1798	!
1799	!-- Drag by plant canopy
1800	IF ( plant_canopy ) CALL plant_canopy_model( 1 )
1801
1802	!
1803	!-- External pressure gradient
1804	IF ( dp_external ) THEN
1805	DO i = i_left, i_right
1806	DO j = j_south, j_north
1807	DO k = dp_level_ind_b+1, nzt
1808	tend(k,j,i) = tend(k,j,i) - dpdxy(1) * dp_smooth_factor(k)
1809	ENDDO
1810	ENDDO
1811	ENDDO
1812	ENDIF
1813
1814	!
1815	!-- Nudging
1816	IF ( nudging ) CALL nudge( simulated_time, 'u' )
1817
1818	CALL user_actions( 'u-tendency' )
1819
1820	!
1821	!-- Prognostic equation for u-velocity component
1822	!$acc kernels present( nzb_u_inner, rdf, tend, tu_m, u, u_init, u_p )
1823	!$acc loop independent
1824	DO i = i_left, i_right
1825	!$acc loop independent
1826	DO j = j_south, j_north
1827	!$acc loop independent
1828	DO k = 1, nzt
1829	IF ( k > nzb_u_inner(j,i) ) THEN
1830	u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
1831	tsc(3) * tu_m(k,j,i) ) &
1832	- tsc(5) * rdf(k) * ( u(k,j,i) - u_init(k) )
1833	!
1834	!-- Tendencies for the next Runge-Kutta step
1835	IF ( runge_step == 1 ) THEN
1836	tu_m(k,j,i) = tend(k,j,i)
1837	ELSEIF ( runge_step == 2 ) THEN
1838	tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tu_m(k,j,i)
1839	ENDIF
1840	ENDIF
1841	ENDDO
1842	ENDDO
1843	ENDDO
1844	!$acc end kernels
1845
1846	CALL cpu_log( log_point(5), 'u-equation', 'stop' )
1847
1848	!
1849	!-- v-velocity component
1850	CALL cpu_log( log_point(6), 'v-equation', 'start' )
1851
1852	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1853	IF ( ws_scheme_mom ) THEN
1854	CALL advec_v_ws_acc
1855	ELSE
1856	tend = 0.0_wp ! to be removed later??
1857	CALL advec_v_pw
1858	END IF
1859	ELSE
1860	CALL advec_v_up
1861	ENDIF
1862	CALL diffusion_v_acc
1863	CALL coriolis_acc( 2 )
1864
1865	!
1866	!-- Drag by plant canopy
1867	IF ( plant_canopy ) CALL plant_canopy_model( 2 )
1868
1869	!
1870	!-- External pressure gradient
1871	IF ( dp_external ) THEN
1872	DO i = i_left, i_right
1873	DO j = j_south, j_north
1874	DO k = dp_level_ind_b+1, nzt
1875	tend(k,j,i) = tend(k,j,i) - dpdxy(2) * dp_smooth_factor(k)
1876	ENDDO
1877	ENDDO
1878	ENDDO
1879	ENDIF
1880
1881	!
1882	!-- Nudging
1883	IF ( nudging ) CALL nudge( simulated_time, 'v' )
1884
1885	CALL user_actions( 'v-tendency' )
1886
1887	!
1888	!-- Prognostic equation for v-velocity component
1889	!$acc kernels present( nzb_v_inner, rdf, tend, tv_m, v, v_init, v_p )
1890	!$acc loop independent
1891	DO i = i_left, i_right
1892	!$acc loop independent
1893	DO j = j_south, j_north
1894	!$acc loop independent
1895	DO k = 1, nzt
1896	IF ( k > nzb_v_inner(j,i) ) THEN
1897	v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
1898	tsc(3) * tv_m(k,j,i) ) &
1899	- tsc(5) * rdf(k) * ( v(k,j,i) - v_init(k) )
1900	!
1901	!-- Tendencies for the next Runge-Kutta step
1902	IF ( runge_step == 1 ) THEN
1903	tv_m(k,j,i) = tend(k,j,i)
1904	ELSEIF ( runge_step == 2 ) THEN
1905	tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tv_m(k,j,i)
1906	ENDIF
1907	ENDIF
1908	ENDDO
1909	ENDDO
1910	ENDDO
1911	!$acc end kernels
1912
1913	CALL cpu_log( log_point(6), 'v-equation', 'stop' )
1914
1915	!
1916	!-- w-velocity component
1917	CALL cpu_log( log_point(7), 'w-equation', 'start' )
1918
1919	IF ( timestep_scheme(1:5) == 'runge' ) THEN
1920	IF ( ws_scheme_mom ) THEN
1921	CALL advec_w_ws_acc
1922	ELSE
1923	tend = 0.0_wp ! to be removed later??
1924	CALL advec_w_pw
1925	ENDIF
1926	ELSE
1927	CALL advec_w_up
1928	ENDIF
1929	CALL diffusion_w_acc
1930	CALL coriolis_acc( 3 )
1931
1932	IF ( .NOT. neutral ) THEN
1933	IF ( ocean ) THEN
1934	CALL buoyancy( rho, 3 )
1935	ELSE
1936	IF ( .NOT. humidity ) THEN
1937	CALL buoyancy_acc( pt, 3 )
1938	ELSE
1939	CALL buoyancy( vpt, 3 )
1940	ENDIF
1941	ENDIF
1942	ENDIF
1943
1944	!
1945	!-- Drag by plant canopy
1946	IF ( plant_canopy ) CALL plant_canopy_model( 3 )
1947
1948	CALL user_actions( 'w-tendency' )
1949
1950	!
1951	!-- Prognostic equation for w-velocity component
1952	!$acc kernels present( nzb_w_inner, rdf, tend, tw_m, w, w_p )
1953	!$acc loop independent
1954	DO i = i_left, i_right
1955	!$acc loop independent
1956	DO j = j_south, j_north
1957	!$acc loop independent
1958	DO k = 1, nzt-1
1959	IF ( k > nzb_w_inner(j,i) ) THEN
1960	w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
1961	tsc(3) * tw_m(k,j,i) ) &
1962	- tsc(5) * rdf(k) * w(k,j,i)
1963	!
1964	!-- Tendencies for the next Runge-Kutta step
1965	IF ( runge_step == 1 ) THEN
1966	tw_m(k,j,i) = tend(k,j,i)
1967	ELSEIF ( runge_step == 2 ) THEN
1968	tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tw_m(k,j,i)
1969	ENDIF
1970	ENDIF
1971	ENDDO
1972	ENDDO
1973	ENDDO
1974	!$acc end kernels
1975
1976	CALL cpu_log( log_point(7), 'w-equation', 'stop' )
1977
1978
1979	!
1980	!-- If required, compute prognostic equation for potential temperature
1981	IF ( .NOT. neutral ) THEN
1982
1983	CALL cpu_log( log_point(13), 'pt-equation', 'start' )
1984
1985	!
1986	!-- pt-tendency terms with communication
1987	sbt = tsc(2)
1988	IF ( scalar_advec == 'bc-scheme' ) THEN
1989
1990	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
1991	!
1992	!-- Bott-Chlond scheme always uses Euler time step. Thus:
1993	sbt = 1.0_wp
1994	ENDIF
1995	tend = 0.0_wp
1996	CALL advec_s_bc( pt, 'pt' )
1997
1998	ENDIF
1999
2000	!
2001	!-- pt-tendency terms with no communication
2002	IF ( scalar_advec /= 'bc-scheme' ) THEN
2003	tend = 0.0_wp
2004	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2005	IF ( ws_scheme_sca ) THEN
2006	CALL advec_s_ws_acc( pt, 'pt' )
2007	ELSE
2008	tend = 0.0_wp ! to be removed later??
2009	CALL advec_s_pw( pt )
2010	ENDIF
2011	ELSE
2012	CALL advec_s_up( pt )
2013	ENDIF
2014	ENDIF
2015
2016	CALL diffusion_s_acc( pt, shf, tswst, wall_heatflux )
2017
2018	!
2019	!-- If required compute heating/cooling due to long wave radiation processes
2020	IF ( radiation ) THEN
2021	CALL calc_radiation
2022	ENDIF
2023
2024	!
2025	!-- If required compute impact of latent heat due to precipitation
2026	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. precipitation ) THEN
2027	CALL impact_of_latent_heat
2028	ENDIF
2029
2030	!
2031	!-- Consideration of heat sources within the plant canopy
2032	IF ( plant_canopy .AND. ( cthf /= 0.0_wp ) ) THEN
2033	CALL plant_canopy_model( 4 )
2034	ENDIF
2035
2036	!
2037	!-- Large scale advection
2038	IF ( large_scale_forcing ) THEN
2039	CALL ls_advec( simulated_time, 'pt' )
2040	ENDIF
2041
2042	!
2043	!-- Nudging
2044	IF ( nudging ) CALL nudge( simulated_time, 'pt' )
2045
2046	!
2047	!-- If required compute influence of large-scale subsidence/ascent
2048	IF ( large_scale_subsidence .AND. &
2049	.NOT. use_subsidence_tendencies ) THEN
2050	CALL subsidence( tend, pt, pt_init, 2 )
2051	ENDIF
2052
2053	CALL user_actions( 'pt-tendency' )
2054
2055	!
2056	!-- Prognostic equation for potential temperature
2057	!$acc kernels present( nzb_s_inner, rdf_sc, ptdf_x, ptdf_y, pt_init ) &
2058	!$acc present( tend, tpt_m, pt, pt_p )
2059	!$acc loop independent
2060	DO i = i_left, i_right
2061	!$acc loop independent
2062	DO j = j_south, j_north
2063	!$acc loop independent
2064	DO k = 1, nzt
2065	IF ( k > nzb_s_inner(j,i) ) THEN
2066	pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2067	tsc(3) * tpt_m(k,j,i) ) &
2068	- tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
2069	( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
2070	!
2071	!-- Tendencies for the next Runge-Kutta step
2072	IF ( runge_step == 1 ) THEN
2073	tpt_m(k,j,i) = tend(k,j,i)
2074	ELSEIF ( runge_step == 2 ) THEN
2075	tpt_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tpt_m(k,j,i)
2076	ENDIF
2077	ENDIF
2078	ENDDO
2079	ENDDO
2080	ENDDO
2081	!$acc end kernels
2082
2083	CALL cpu_log( log_point(13), 'pt-equation', 'stop' )
2084
2085	ENDIF
2086
2087	!
2088	!-- If required, compute prognostic equation for salinity
2089	IF ( ocean ) THEN
2090
2091	CALL cpu_log( log_point(37), 'sa-equation', 'start' )
2092
2093	!
2094	!-- sa-tendency terms with communication
2095	sbt = tsc(2)
2096	IF ( scalar_advec == 'bc-scheme' ) THEN
2097
2098	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
2099	!
2100	!-- Bott-Chlond scheme always uses Euler time step. Thus:
2101	sbt = 1.0_wp
2102	ENDIF
2103	tend = 0.0_wp
2104	CALL advec_s_bc( sa, 'sa' )
2105
2106	ENDIF
2107
2108	!
2109	!-- sa-tendency terms with no communication
2110	IF ( scalar_advec /= 'bc-scheme' ) THEN
2111	tend = 0.0_wp
2112	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2113	IF ( ws_scheme_sca ) THEN
2114	CALL advec_s_ws( sa, 'sa' )
2115	ELSE
2116	CALL advec_s_pw( sa )
2117	ENDIF
2118	ELSE
2119	CALL advec_s_up( sa )
2120	ENDIF
2121	ENDIF
2122
2123	CALL diffusion_s( sa, saswsb, saswst, wall_salinityflux )
2124
2125	CALL user_actions( 'sa-tendency' )
2126
2127	!
2128	!-- Prognostic equation for salinity
2129	DO i = i_left, i_right
2130	DO j = j_south, j_north
2131	DO k = nzb_s_inner(j,i)+1, nzt
2132	sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2133	tsc(3) * tsa_m(k,j,i) ) &
2134	- tsc(5) * rdf_sc(k) * &
2135	( sa(k,j,i) - sa_init(k) )
2136	IF ( sa_p(k,j,i) < 0.0_wp ) sa_p(k,j,i) = 0.1_wp * sa(k,j,i)
2137	!
2138	!-- Tendencies for the next Runge-Kutta step
2139	IF ( runge_step == 1 ) THEN
2140	tsa_m(k,j,i) = tend(k,j,i)
2141	ELSEIF ( runge_step == 2 ) THEN
2142	tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tsa_m(k,j,i)
2143	ENDIF
2144	ENDDO
2145	ENDDO
2146	ENDDO
2147
2148	CALL cpu_log( log_point(37), 'sa-equation', 'stop' )
2149
2150	!
2151	!-- Calculate density by the equation of state for seawater
2152	CALL cpu_log( log_point(38), 'eqns-seawater', 'start' )
2153	CALL eqn_state_seawater
2154	CALL cpu_log( log_point(38), 'eqns-seawater', 'stop' )
2155
2156	ENDIF
2157
2158	!
2159	!-- If required, compute prognostic equation for total water content / scalar
2160	IF ( humidity .OR. passive_scalar ) THEN
2161
2162	CALL cpu_log( log_point(29), 'q/s-equation', 'start' )
2163
2164	!
2165	!-- Scalar/q-tendency terms with communication
2166	sbt = tsc(2)
2167	IF ( scalar_advec == 'bc-scheme' ) THEN
2168
2169	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
2170	!
2171	!-- Bott-Chlond scheme always uses Euler time step. Thus:
2172	sbt = 1.0_wp
2173	ENDIF
2174	tend = 0.0_wp
2175	CALL advec_s_bc( q, 'q' )
2176
2177	ENDIF
2178
2179	!
2180	!-- Scalar/q-tendency terms with no communication
2181	IF ( scalar_advec /= 'bc-scheme' ) THEN
2182	tend = 0.0_wp
2183	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2184	IF ( ws_scheme_sca ) THEN
2185	CALL advec_s_ws( q, 'q' )
2186	ELSE
2187	CALL advec_s_pw( q )
2188	ENDIF
2189	ELSE
2190	CALL advec_s_up( q )
2191	ENDIF
2192	ENDIF
2193
2194	CALL diffusion_s( q, qsws, qswst, wall_qflux )
2195
2196	!
2197	!-- If required compute decrease of total water content due to
2198	!-- precipitation
2199	IF ( cloud_physics .AND. icloud_scheme == 1 .AND. precipitation ) THEN
2200	CALL calc_precipitation
2201	ENDIF
2202
2203	!
2204	!-- Sink or source of scalar concentration due to canopy elements
2205	IF ( plant_canopy ) CALL plant_canopy_model( 5 )
2206
2207	!
2208	!-- Large scale advection
2209	IF ( large_scale_forcing ) THEN
2210	CALL ls_advec( simulated_time, 'q' )
2211	ENDIF
2212
2213	!
2214	!-- Nudging
2215	IF ( nudging ) CALL nudge( simulated_time, 'q' )
2216
2217	!
2218	!-- If required compute influence of large-scale subsidence/ascent
2219	IF ( large_scale_subsidence .AND. &
2220	.NOT. use_subsidence_tendencies ) THEN
2221	CALL subsidence( tend, q, q_init, 3 )
2222	ENDIF
2223
2224	CALL user_actions( 'q-tendency' )
2225
2226	!
2227	!-- Prognostic equation for total water content / scalar
2228	DO i = i_left, i_right
2229	DO j = j_south, j_north
2230	DO k = nzb_s_inner(j,i)+1, nzt
2231	q_p(k,j,i) = q(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2232	tsc(3) * tq_m(k,j,i) ) &
2233	- tsc(5) * rdf_sc(k) * &
2234	( q(k,j,i) - q_init(k) )
2235	IF ( q_p(k,j,i) < 0.0_wp ) q_p(k,j,i) = 0.1_wp * q(k,j,i)
2236	!
2237	!-- Tendencies for the next Runge-Kutta step
2238	IF ( runge_step == 1 ) THEN
2239	tq_m(k,j,i) = tend(k,j,i)
2240	ELSEIF ( runge_step == 2 ) THEN
2241	tq_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tq_m(k,j,i)
2242	ENDIF
2243	ENDDO
2244	ENDDO
2245	ENDDO
2246
2247	CALL cpu_log( log_point(29), 'q/s-equation', 'stop' )
2248
2249	!
2250	!-- If required, calculate prognostic equations for rain water content
2251	!-- and rain drop concentration
2252	IF ( cloud_physics .AND. icloud_scheme == 0 .AND. precipitation ) THEN
2253
2254	CALL cpu_log( log_point(52), 'qr-equation', 'start' )
2255	!
2256	!-- qr-tendency terms with communication
2257	sbt = tsc(2)
2258	IF ( scalar_advec == 'bc-scheme' ) THEN
2259
2260	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
2261	!
2262	!-- Bott-Chlond scheme always uses Euler time step. Thus:
2263	sbt = 1.0_wp
2264	ENDIF
2265	tend = 0.0_wp
2266	CALL advec_s_bc( qr, 'qr' )
2267
2268	ENDIF
2269
2270	!
2271	!-- qr-tendency terms with no communication
2272	IF ( scalar_advec /= 'bc-scheme' ) THEN
2273	tend = 0.0_wp
2274	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2275	IF ( ws_scheme_sca ) THEN
2276	CALL advec_s_ws( qr, 'qr' )
2277	ELSE
2278	CALL advec_s_pw( qr )
2279	ENDIF
2280	ELSE
2281	CALL advec_s_up( qr )
2282	ENDIF
2283	ENDIF
2284
2285	CALL diffusion_s( qr, qrsws, qrswst, wall_qrflux )
2286
2287	CALL user_actions( 'qr-tendency' )
2288
2289	!
2290	!-- Prognostic equation for rain water content
2291	DO i = i_left, i_right
2292	DO j = j_south, j_north
2293	DO k = nzb_s_inner(j,i)+1, nzt
2294	qr_p(k,j,i) = qr(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2295	tsc(3) * tqr_m(k,j,i) ) &
2296	- tsc(5) * rdf_sc(k) * qr(k,j,i)
2297	IF ( qr_p(k,j,i) < 0.0_wp ) qr_p(k,j,i) = 0.0_wp
2298	!
2299	!-- Tendencies for the next Runge-Kutta step
2300	IF ( runge_step == 1 ) THEN
2301	tqr_m(k,j,i) = tend(k,j,i)
2302	ELSEIF ( runge_step == 2 ) THEN
2303	tqr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * &
2304	tqr_m(k,j,i)
2305	ENDIF
2306	ENDDO
2307	ENDDO
2308	ENDDO
2309
2310	CALL cpu_log( log_point(52), 'qr-equation', 'stop' )
2311	CALL cpu_log( log_point(53), 'nr-equation', 'start' )
2312
2313	!
2314	!-- nr-tendency terms with communication
2315	sbt = tsc(2)
2316	IF ( scalar_advec == 'bc-scheme' ) THEN
2317
2318	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
2319	!
2320	!-- Bott-Chlond scheme always uses Euler time step. Thus:
2321	sbt = 1.0_wp
2322	ENDIF
2323	tend = 0.0_wp
2324	CALL advec_s_bc( nr, 'nr' )
2325
2326	ENDIF
2327
2328	!
2329	!-- nr-tendency terms with no communication
2330	IF ( scalar_advec /= 'bc-scheme' ) THEN
2331	tend = 0.0_wp
2332	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2333	IF ( ws_scheme_sca ) THEN
2334	CALL advec_s_ws( nr, 'nr' )
2335	ELSE
2336	CALL advec_s_pw( nr )
2337	ENDIF
2338	ELSE
2339	CALL advec_s_up( nr )
2340	ENDIF
2341	ENDIF
2342
2343	CALL diffusion_s( nr, nrsws, nrswst, wall_nrflux )
2344
2345	CALL user_actions( 'nr-tendency' )
2346
2347	!
2348	!-- Prognostic equation for rain drop concentration
2349	DO i = i_left, i_right
2350	DO j = j_south, j_north
2351	DO k = nzb_s_inner(j,i)+1, nzt
2352	nr_p(k,j,i) = nr(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2353	tsc(3) * tnr_m(k,j,i) ) &
2354	- tsc(5) * rdf_sc(k) * nr(k,j,i)
2355	IF ( nr_p(k,j,i) < 0.0_wp ) nr_p(k,j,i) = 0.0_wp
2356	!
2357	!-- Tendencies for the next Runge-Kutta step
2358	IF ( runge_step == 1 ) THEN
2359	tnr_m(k,j,i) = tend(k,j,i)
2360	ELSEIF ( runge_step == 2 ) THEN
2361	tnr_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * &
2362	tnr_m(k,j,i)
2363	ENDIF
2364	ENDDO
2365	ENDDO
2366	ENDDO
2367
2368	CALL cpu_log( log_point(53), 'nr-equation', 'stop' )
2369
2370	ENDIF
2371
2372	ENDIF
2373
2374	!
2375	!-- If required, compute prognostic equation for turbulent kinetic
2376	!-- energy (TKE)
2377	IF ( .NOT. constant_diffusion ) THEN
2378
2379	CALL cpu_log( log_point(16), 'tke-equation', 'start' )
2380
2381	sbt = tsc(2)
2382	IF ( .NOT. use_upstream_for_tke ) THEN
2383	IF ( scalar_advec == 'bc-scheme' ) THEN
2384
2385	IF ( timestep_scheme(1:5) /= 'runge' ) THEN
2386	!
2387	!-- Bott-Chlond scheme always uses Euler time step. Thus:
2388	sbt = 1.0_wp
2389	ENDIF
2390	tend = 0.0_wp
2391	CALL advec_s_bc( e, 'e' )
2392
2393	ENDIF
2394	ENDIF
2395
2396	!
2397	!-- TKE-tendency terms with no communication
2398	IF ( scalar_advec /= 'bc-scheme' .OR. use_upstream_for_tke ) THEN
2399	IF ( use_upstream_for_tke ) THEN
2400	tend = 0.0_wp
2401	CALL advec_s_up( e )
2402	ELSE
2403	IF ( timestep_scheme(1:5) == 'runge' ) THEN
2404	IF ( ws_scheme_sca ) THEN
2405	CALL advec_s_ws_acc( e, 'e' )
2406	ELSE
2407	tend = 0.0_wp ! to be removed later??
2408	CALL advec_s_pw( e )
2409	ENDIF
2410	ELSE
2411	tend = 0.0_wp ! to be removed later??
2412	CALL advec_s_up( e )
2413	ENDIF
2414	ENDIF
2415	ENDIF
2416
2417	IF ( .NOT. humidity ) THEN
2418	IF ( ocean ) THEN
2419	CALL diffusion_e( prho, prho_reference )
2420	ELSE
2421	CALL diffusion_e_acc( pt, pt_reference )
2422	ENDIF
2423	ELSE
2424	CALL diffusion_e( vpt, pt_reference )
2425	ENDIF
2426
2427	CALL production_e_acc
2428
2429	!
2430	!-- Additional sink term for flows through plant canopies
2431	IF ( plant_canopy ) CALL plant_canopy_model( 6 )
2432	CALL user_actions( 'e-tendency' )
2433
2434	!
2435	!-- Prognostic equation for TKE.
2436	!-- Eliminate negative TKE values, which can occur due to numerical
2437	!-- reasons in the course of the integration. In such cases the old TKE
2438	!-- value is reduced by 90%.
2439	!$acc kernels present( e, e_p, nzb_s_inner, tend, te_m )
2440	!$acc loop independent
2441	DO i = i_left, i_right
2442	!$acc loop independent
2443	DO j = j_south, j_north
2444	!$acc loop independent
2445	DO k = 1, nzt
2446	IF ( k > nzb_s_inner(j,i) ) THEN
2447	e_p(k,j,i) = e(k,j,i) + dt_3d * ( sbt * tend(k,j,i) + &
2448	tsc(3) * te_m(k,j,i) )
2449	IF ( e_p(k,j,i) < 0.0_wp ) e_p(k,j,i) = 0.1_wp * e(k,j,i)
2450	!
2451	!-- Tendencies for the next Runge-Kutta step
2452	IF ( runge_step == 1 ) THEN
2453	te_m(k,j,i) = tend(k,j,i)
2454	ELSEIF ( runge_step == 2 ) THEN
2455	te_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * te_m(k,j,i)
2456	ENDIF
2457	ENDIF
2458	ENDDO
2459	ENDDO
2460	ENDDO
2461	!$acc end kernels
2462
2463	CALL cpu_log( log_point(16), 'tke-equation', 'stop' )
2464
2465	ENDIF
2466
2467	END SUBROUTINE prognostic_equations_acc
2468
2469
2470	END MODULE prognostic_equations_mod

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: palm/trunk/SOURCE/prognostic_equations.f90 @ 1484

Download in other formats: