Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |