Home

Context Navigation

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

Last change on this file since 1691 was 1691, checked in by maronga, 8 years ago
various bugfixes and modifications of the atmosphere-land-surface-radiation interaction. Completely re-written routine to calculate surface fluxes (surface_layer_fluxes.f90) that replaces prandtl_fluxes. Minor formatting corrections and renamings
Property svn:keywords set to `Id`
File size: 85.4 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |