Home

Context Navigation

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

Last change on this file since 709 was 709, checked in by raasch, 13 years ago
formatting adjustments
Property svn:keywords set to `Id`
File size: 75.3 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |