Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |