Home

Context Navigation

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

Last change on this file since 3826 was 3820, checked in by forkel, 5 years ago
renaming of get_mechanismname, do_emiss and do_depo, sorting in chem_modules
Property svn:keywords set to `Id`
File size: 104.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |