Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |