Home

Context Navigation

source: palm/trunk/SOURCE/surface_mod.f90 @ 3118

Last change on this file since 3118 was 3055, checked in by suehring, 6 years ago
Bugfix, initialization of surface elements also in case of restart runs and cyclic fill
Property svn:keywords set to `Id`
File size: 211.0 KB

Line
1	!> @file surface_mod.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
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.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 1997-2018 Leibniz Universitaet Hannover
18	!
19	!------------------------------------------------------------------------------!
20	!
21	! Current revisions:
22	! ------------------
23	!
24	!
25	! Former revisions:
26	! -----------------
27	! $Id: surface_mod.f90 3055 2018-06-01 17:15:08Z knoop $
28	! Initialize Obukhov length, friction velocity and momentum fluxes also
29	! in case of restarts and cyclic fill
30	!
31	! 3026 2018-05-22 10:30:53Z schwenkel
32	! Changed the name specific humidity to mixing ratio, since we are computing
33	! mixing ratios.
34	!
35	! 2977 2018-04-17 10:27:57Z kanani
36	! Implement changes from branch radiation (r2948-2971) with minor modifications,
37	! plus some formatting.
38	! (moh.hefny)
39	! Added flag to check the existence of vertical urban/land surfaces, required
40	! to activate RTM
41	!
42	! 2970 2018-04-13 15:09:23Z suehring
43	! Remove un-necessary initialization of surface elements in old large-scale
44	! forcing mode
45	!
46	! 2963 2018-04-12 14:47:44Z suehring
47	! Introduce index for vegetation/wall, pavement/green-wall and water/window
48	! surfaces, for clearer access of surface fraction, albedo, emissivity, etc. .
49	!
50	! 2942 2018-04-03 13:51:09Z suehring
51	! Bugfix in assigning surface element data after restart
52	!
53	! 2940 2018-04-03 11:22:42Z suehring
54	! Bugfix in reading restart data of vertical surface elements
55	!
56	! 2920 2018-03-22 11:22:01Z kanani
57	! Correct comment for surface directions
58	!
59	! 2894 2018-03-15 09:17:58Z Giersch
60	! Calculations of the index range of the subdomain on file which overlaps with
61	! the current subdomain are already done in read_restart_data_mod,
62	! surface_read/write_restart_data was renamed to surface_r/wrd_local, variable
63	! named found has been introduced for checking if restart data was found,
64	! reading of restart strings has been moved completely to
65	! read_restart_data_mod, surface_rrd_local is already inside the overlap loop
66	! programmed in read_restart_data_mod, SAVE attribute added where necessary,
67	! deallocation and allocation of some arrays have been changed to take care of
68	! different restart files that can be opened (index i), the marker *** end
69	! surf *** is not necessary anymore, strings and their respective lengths are
70	! written out and read now in case of restart runs to get rid of prescribed
71	! character lengths (Giersch)
72	!
73	! 2813 2018-02-16 16:28:14Z suehring
74	! Some more bugfixes concerning restart runs
75	!
76	! 2812 2018-02-16 13:40:25Z hellstea
77	! Entries 'u_out', 'v_out' and 'w_out' removed from the functions
78	! get_topography_top_index and get_topography_top_index_ji
79	!
80	! 2805 2018-02-14 17:00:09Z suehring
81	! Bugfix in re-mapping surface elements in case of restart runs
82	!
83	! 2766 2018-01-22 17:17:47Z kanani
84	! Removed preprocessor directive __chem
85	!
86	! 2759 2018-01-17 16:24:59Z suehring
87	! Bugfix, consider density in vertical fluxes of passive scalar as well as
88	! chemical species.
89	!
90	! 2753 2018-01-16 14:16:49Z suehring
91	! +r_a_green, r_a_window
92	!
93	! 2718 2018-01-02 08:49:38Z maronga
94	! Changes from last commit documented
95	!
96	! 2706 2017-12-18 18:33:49Z suehring
97	! In case of restarts read and write pt_surface
98	!
99	! 2698 2017-12-14 18:46:24Z suehring
100	!
101	! 2696 2017-12-14 17:12:51Z kanani
102	! - Change in file header (GPL part)
103	! - Missing code for css added to surf_*, handling of surface_csflux updated (FK)
104	! - Bugfixes in reading/writing restart data in case several surface types are
105	! present at the same time (MS)
106	! - Implementation of chemistry module (FK)
107	! - Allocation of pt1 and qv1 now done for all surface types (MS)
108	! - Revised classification of surface types
109	! - Introduce offset values to simplify index determination of surface elements
110	! - Dimensions of function get_topo_top_index (MS)
111	! - added variables for USM
112	! - adapted to changes in USM (RvT)
113	!
114	! 2688 2017-12-12 17:27:04Z Giersch
115	! Allocation and initialization of the latent heat flux (qsws) at the top of
116	! the ocean domain in case of coupled runs. In addtion, a double if-query has
117	! been removed.
118	!
119	! 2638 2017-11-23 12:44:23Z raasch
120	! bugfix for constant top momentumflux
121	!
122	! 2575 2017-10-24 09:57:58Z maronga
123	! Pavement parameterization revised
124	!
125	! 2547 2017-10-16 12:41:56Z schwenkel
126	! extended by cloud_droplets option
127	!
128	! 2508 2017-10-02 08:57:09Z suehring
129	! Minor formatting adjustment
130	!
131	! 2478 2017-09-18 13:37:24Z suehring
132	! Bugfixes in initializing model top
133	!
134	! 2378 2017-08-31 13:57:27Z suehring
135	! Bugfix in write restart data
136	!
137	! 2339 2017-08-07 13:55:26Z gronemeier
138	! corrected timestamp in header
139	!
140	! 2338 2017-08-07 12:15:38Z gronemeier
141	! Modularize 1D model
142	!
143	! 2318 2017-07-20 17:27:44Z suehring
144	! New function to obtain topography top index.
145	!
146	! 2317 2017-07-20 17:27:19Z suehring
147	! Implementation of new microphysic scheme: cloud_scheme = 'morrison'
148	! includes two more prognostic equations for cloud drop concentration (nc)
149	! and cloud water content (qc).
150	!
151	! 2270 2017-06-09 12:18:47Z maronga
152	! Parameters removed/added due to changes in the LSM
153	!
154	! 2269 2017-06-09 11:57:32Z suehring
155	! Formatting and description adjustments
156	!
157	! 2256 2017-06-07 13:58:08Z suehring
158	! Enable heating at downward-facing surfaces
159	!
160	! 2233 2017-05-30 18:08:54Z suehring
161	! Initial revision
162	!
163	!
164	! Description:
165	! ------------
166	!> Surface module defines derived data structures to treat surface-
167	!> bounded grid cells. Three different types of surfaces are defined:
168	!> default surfaces, natural surfaces, and urban surfaces. The module
169	!> encompasses the allocation and initialization of surface arrays, and handles
170	!> reading and writing restart data.
171	!> In addition, a further derived data structure is defined, in order to set
172	!> boundary conditions at surfaces.
173	!> @todo For the moment, downward-facing surfaces are only classified as
174	!> default type
175	!> @todo Clean up urban-surface variables (some of them are not used any more)
176	!> @todo Revise chemistry surface flux part (reduce loops?!)
177	!------------------------------------------------------------------------------!
178	MODULE surface_mod
179
180	USE arrays_3d, &
181	ONLY: heatflux_input_conversion, momentumflux_input_conversion, &
182	rho_air, rho_air_zw, zu, zw, waterflux_input_conversion
183
184	USE chem_modules
185
186	USE control_parameters
187
188	USE indices, &
189	ONLY: nxl, nxlg, nxr, nxrg, nys, nysg, nyn, nyng, nzb, nzt, wall_flags_0
190
191	USE grid_variables, &
192	ONLY: dx, dy
193
194	USE kinds
195
196	USE model_1d_mod, &
197	ONLY: rif1d, us1d, usws1d, vsws1d
198
199
200	IMPLICIT NONE
201
202	!
203	!-- Data type used to identify grid-points where horizontal boundary conditions
204	!-- are applied
205	TYPE bc_type
206
207	INTEGER(iwp) :: ns !< number of surface elements on the PE
208
209	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: i !< x-index linking to the PALM 3D-grid
210	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: j !< y-index linking to the PALM 3D-grid
211	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: k !< z-index linking to the PALM 3D-grid
212
213	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: start_index !< start index within surface data type for given (j,i)
214	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: end_index !< end index within surface data type for given (j,i)
215
216	END TYPE bc_type
217	!
218	!-- Data type used to identify and treat surface-bounded grid points
219	TYPE surf_type
220
221	INTEGER(iwp) :: ioff !< offset value in x-direction, used to determine index of surface element
222	INTEGER(iwp) :: joff !< offset value in y-direction, used to determine index of surface element
223	INTEGER(iwp) :: koff !< offset value in z-direction, used to determine index of surface element
224	INTEGER(iwp) :: ns !< number of surface elements on the PE
225
226	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: i !< x-index linking to the PALM 3D-grid
227	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: j !< y-index linking to the PALM 3D-grid
228	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: k !< z-index linking to the PALM 3D-grid
229
230	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: facing !< Bit indicating surface orientation
231
232	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: start_index !< Start index within surface data type for given (j,i)
233	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: end_index !< End index within surface data type for given (j,i)
234
235	REAL(wp), DIMENSION(:), ALLOCATABLE :: z_mo !< surface-layer height
236	REAL(wp), DIMENSION(:), ALLOCATABLE :: uvw_abs !< absolute surface-parallel velocity
237	REAL(wp), DIMENSION(:), ALLOCATABLE :: us !< friction velocity
238	REAL(wp), DIMENSION(:), ALLOCATABLE :: ts !< scaling parameter temerature
239	REAL(wp), DIMENSION(:), ALLOCATABLE :: qs !< scaling parameter humidity
240	REAL(wp), DIMENSION(:), ALLOCATABLE :: ss !< scaling parameter passive scalar
241	REAL(wp), DIMENSION(:), ALLOCATABLE :: qcs !< scaling parameter qc
242	REAL(wp), DIMENSION(:), ALLOCATABLE :: ncs !< scaling parameter nc
243	REAL(wp), DIMENSION(:), ALLOCATABLE :: qrs !< scaling parameter qr
244	REAL(wp), DIMENSION(:), ALLOCATABLE :: nrs !< scaling parameter nr
245
246	REAL(wp), DIMENSION(:), ALLOCATABLE :: ol !< Obukhov length
247	REAL(wp), DIMENSION(:), ALLOCATABLE :: rib !< Richardson bulk number
248
249	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0 !< roughness length for momentum
250	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0h !< roughness length for heat
251	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0q !< roughness length for humidity
252
253	REAL(wp), DIMENSION(:), ALLOCATABLE :: pt1 !< Potential temperature at first grid level
254	REAL(wp), DIMENSION(:), ALLOCATABLE :: qv1 !< mixing ratio at first grid level
255	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: css !< scaling parameter chemical species
256	!
257	!-- Define arrays for surface fluxes
258	REAL(wp), DIMENSION(:), ALLOCATABLE :: usws !< vertical momentum flux for u-component at horizontal surfaces
259	REAL(wp), DIMENSION(:), ALLOCATABLE :: vsws !< vertical momentum flux for v-component at horizontal surfaces
260
261	REAL(wp), DIMENSION(:), ALLOCATABLE :: shf !< surface flux sensible heat
262	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws !< surface flux latent heat
263	REAL(wp), DIMENSION(:), ALLOCATABLE :: ssws !< surface flux passive scalar
264	REAL(wp), DIMENSION(:), ALLOCATABLE :: qcsws !< surface flux qc
265	REAL(wp), DIMENSION(:), ALLOCATABLE :: ncsws !< surface flux nc
266	REAL(wp), DIMENSION(:), ALLOCATABLE :: qrsws !< surface flux qr
267	REAL(wp), DIMENSION(:), ALLOCATABLE :: nrsws !< surface flux nr
268	REAL(wp), DIMENSION(:), ALLOCATABLE :: sasws !< surface flux salinity
269
270	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: cssws !< surface flux chemical species
271	!
272	!-- Required for horizontal walls in production_e
273	REAL(wp), DIMENSION(:), ALLOCATABLE :: u_0 !< virtual velocity component (see production_e_init for further explanation)
274	REAL(wp), DIMENSION(:), ALLOCATABLE :: v_0 !< virtual velocity component (see production_e_init for further explanation)
275
276	REAL(wp), DIMENSION(:), ALLOCATABLE :: mom_flux_uv !< momentum flux usvs and vsus at vertical surfaces (used in diffusion_u and diffusion_v)
277	REAL(wp), DIMENSION(:), ALLOCATABLE :: mom_flux_w !< momentum flux wsus and wsvs at vertical surfaces (used in diffusion_w)
278	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: mom_flux_tke !< momentum flux usvs, vsus, wsus, wsvs at vertical surfaces at grid center (used in production_e)
279	!
280	!-- Variables required for LSM as well as for USM
281	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: nzt_pavement !< top index for pavement in soil
282	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: albedo_type !< albedo type, for each fraction (wall,green,window or vegetation,pavement water)
283
284	LOGICAL, DIMENSION(:), ALLOCATABLE :: building_surface !< flag parameter indicating that the surface element is covered by buildings (no LSM actions, not implemented yet)
285	LOGICAL, DIMENSION(:), ALLOCATABLE :: building_covered !< flag indicating that buildings are on top of orography, only used for vertical surfaces in LSM
286	LOGICAL, DIMENSION(:), ALLOCATABLE :: pavement_surface !< flag parameter for pavements
287	LOGICAL, DIMENSION(:), ALLOCATABLE :: water_surface !< flag parameter for water surfaces
288	LOGICAL, DIMENSION(:), ALLOCATABLE :: vegetation_surface !< flag parameter for natural land surfaces
289
290	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: albedo !< broadband albedo for each surface fraction (LSM: vegetation, water, pavement; USM: wall, green, window)
291	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: emissivity !< emissivity of the surface, for each fraction (LSM: vegetation, water, pavement; USM: wall, green, window)
292	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: frac !< relative surface fraction (LSM: vegetation, water, pavement; USM: wall, green, window)
293
294	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: aldif !< albedo for longwave diffusive radiation, solar angle of 60Â°
295	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: aldir !< albedo for longwave direct radiation, solar angle of 60Â°
296	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: asdif !< albedo for shortwave diffusive radiation, solar angle of 60Â°
297	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: asdir !< albedo for shortwave direct radiation, solar angle of 60Â°
298	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rrtm_aldif !< albedo for longwave diffusive radiation, solar angle of 60Â°
299	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rrtm_aldir !< albedo for longwave direct radiation, solar angle of 60Â°
300	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rrtm_asdif !< albedo for shortwave diffusive radiation, solar angle of 60Â°
301	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rrtm_asdir !< albedo for shortwave direct radiation, solar angle of 60Â°
302
303	REAL(wp), DIMENSION(:), ALLOCATABLE :: pt_surface !< skin-surface temperature
304	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net !< net radiation
305	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net_l !< net radiation, used in USM
306	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_h !< heat conductivity of soil/ wall (W/m/K)
307	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_h_green !< heat conductivity of green soil (W/m/K)
308	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_h_window !< heat conductivity of windows (W/m/K)
309	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_h_def !< default heat conductivity of soil (W/m/K)
310
311	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_lw_in !< incoming longwave radiation
312	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_lw_out !< emitted longwave radiation
313	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_sw_in !< incoming shortwave radiation
314	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_sw_out !< emitted shortwave radiation
315
316
317
318	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_liq !< liquid water coverage (of vegetated area)
319	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_veg !< vegetation coverage
320	REAL(wp), DIMENSION(:), ALLOCATABLE :: f_sw_in !< fraction of absorbed shortwave radiation by the surface layer (not implemented yet)
321	REAL(wp), DIMENSION(:), ALLOCATABLE :: ghf !< ground heat flux
322	REAL(wp), DIMENSION(:), ALLOCATABLE :: g_d !< coefficient for dependence of r_canopy on water vapour pressure deficit
323	REAL(wp), DIMENSION(:), ALLOCATABLE :: lai !< leaf area index
324	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surface_u !< coupling between surface and soil (depends on vegetation type) (W/m2/K)
325	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surface_s !< coupling between surface and soil (depends on vegetation type) (W/m2/K)
326	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_liq !< surface flux of latent heat (liquid water portion)
327	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_soil !< surface flux of latent heat (soil portion)
328	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_veg !< surface flux of latent heat (vegetation portion)
329	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_a !< aerodynamic resistance
330	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_a_green !< aerodynamic resistance at green fraction
331	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_a_window !< aerodynamic resistance at window fraction
332	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_canopy !< canopy resistance
333	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_soil !< soil resistance
334	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_soil_min !< minimum soil resistance
335	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_s !< total surface resistance (combination of r_soil and r_canopy)
336	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_canopy_min !< minimum canopy (stomatal) resistance
337
338	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: alpha_vg !< coef. of Van Genuchten
339	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_w !< hydraulic diffusivity of soil (?)
340	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: gamma_w !< hydraulic conductivity of soil (W/m/K)
341	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: gamma_w_sat !< hydraulic conductivity at saturation
342	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: l_vg !< coef. of Van Genuchten
343	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_fc !< soil moisture at field capacity (m3/m3)
344	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_res !< residual soil moisture
345	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_sat !< saturation soil moisture (m3/m3)
346	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_wilt !< soil moisture at permanent wilting point (m3/m3)
347	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: n_vg !< coef. Van Genuchten
348	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_total_def !< default volumetric heat capacity of the (soil) layer (J/m3/K)
349	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_total !< volumetric heat capacity of the actual soil matrix (J/m3/K)
350	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: root_fr !< root fraction within the soil layers
351	!
352	!-- Urban surface variables
353	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: surface_types !< array of types of wall parameters
354
355	LOGICAL, DIMENSION(:), ALLOCATABLE :: isroof_surf !< flag indicating roof surfaces
356	LOGICAL, DIMENSION(:), ALLOCATABLE :: ground_level !< flag indicating ground floor level surfaces
357
358	REAL(wp), DIMENSION(:), ALLOCATABLE :: target_temp_summer !< indoor target temperature summer
359	REAL(wp), DIMENSION(:), ALLOCATABLE :: target_temp_winter !< indoor target temperature summer
360
361	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_surface !< heat capacity of the wall surface skin (J/m2/K)
362	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_surface_green !< heat capacity of the green surface skin (J/m2/K)
363	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_surface_window !< heat capacity of the window surface skin (J/m2/K)
364	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surf !< heat conductivity between air and surface (W/m2/K)
365	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surf_green !< heat conductivity between air and green surface (W/m2/K)
366	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surf_window !< heat conductivity between air and window surface (W/m2/K)
367	REAL(wp), DIMENSION(:), ALLOCATABLE :: thickness_wall !< thickness of the wall, roof and soil layers
368	REAL(wp), DIMENSION(:), ALLOCATABLE :: thickness_green !< thickness of the green wall, roof and soil layers
369	REAL(wp), DIMENSION(:), ALLOCATABLE :: thickness_window !< thickness of the window wall, roof and soil layers
370	REAL(wp), DIMENSION(:), ALLOCATABLE :: transmissivity !< transmissivity of windows
371
372	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsl !< reflected shortwave radiation for local surface in i-th reflection
373	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutll !< reflected + emitted longwave radiation for local surface in i-th reflection
374	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf !< total radiation flux incoming to minus outgoing from local surface
375
376	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_m !< surface temperature tendency (K)
377	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_window_m !< window surface temperature tendency (K)
378	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_green_m !< green surface temperature tendency (K)
379	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf !< kinematic wall heat flux of sensible heat (actually no longer needed)
380	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf_eb !< wall heat flux of sensible heat in wall normal direction
381
382	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb !< wall ground heat flux
383	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_window !< window ground heat flux
384	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_green !< green ground heat flux
385	REAL(wp), DIMENSION(:), ALLOCATABLE :: iwghf_eb !< indoor wall ground heat flux
386	REAL(wp), DIMENSION(:), ALLOCATABLE :: iwghf_eb_window !< indoor window ground heat flux
387
388	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_lw_out_change_0
389
390	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinsw !< shortwave radiation falling to local surface including radiation from reflections
391	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsw !< total shortwave radiation outgoing from nonvirtual surfaces surfaces after all reflection
392	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlw !< longwave radiation falling to local surface including radiation from reflections
393	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw !< total longwave radiation outgoing from nonvirtual surfaces surfaces after all reflection
394
395
396	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_wall !< volumetric heat capacity of the material ( J m-3 K-1 ) (= 2.19E6)
397	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_wall !< wall grid spacing (center-center)
398	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_wall !< 1/dz_wall
399	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_wall_stag !< wall grid spacing (edge-edge)
400	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_wall_stag !< 1/dz_wall_stag
401	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_wall_m !< t_wall prognostic array
402	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zw !< wall layer depths (m)
403	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_window !< volumetric heat capacity of the window material ( J m-3 K-1 ) (= 2.19E6)
404	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_window !< window grid spacing (center-center)
405	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_window !< 1/dz_window
406	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_window_stag !< window grid spacing (edge-edge)
407	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_window_stag !< 1/dz_window_stag
408	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_window_m !< t_window prognostic array
409	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zw_window !< window layer depths (m)
410	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_green !< volumetric heat capacity of the green material ( J m-3 K-1 ) (= 2.19E6)
411	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_green !< green grid spacing (center-center)
412	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_green !< 1/dz_green
413	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_green_stag !< green grid spacing (edge-edge)
414	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_green_stag !< 1/dz_green_stag
415	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_green_m !< t_green prognostic array
416	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zw_green !< green layer depths (m)
417
418
419	!-- arrays for time averages
420	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net_av !< average of rad_net_l
421	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinsw_av !< average of sw radiation falling to local surface including radiation from reflections
422	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlw_av !< average of lw radiation falling to local surface including radiation from reflections
423	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdir_av !< average of direct sw radiation falling to local surface
424	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdif_av !< average of diffuse sw radiation from sky and model boundary falling to local surface
425	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwdif_av !< average of diffuse lw radiation from sky and model boundary falling to local surface
426	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswref_av !< average of sw radiation falling to surface from reflections
427	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwref_av !< average of lw radiation falling to surface from reflections
428	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsw_av !< average of total sw radiation outgoing from nonvirtual surfaces surfaces after all reflection
429	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw_av !< average of total lw radiation outgoing from nonvirtual surfaces surfaces after all reflection
430	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfins_av !< average of array of residua of sw radiation absorbed in surface after last reflection
431	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinl_av !< average of array of residua of lw radiation absorbed in surface after last reflection
432	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf_av !< average of total radiation flux incoming to minus outgoing from local surface
433	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_av !< average of wghf_eb
434	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_window_av !< average of wghf_eb window
435	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_green_av !< average of wghf_eb window
436	REAL(wp), DIMENSION(:), ALLOCATABLE :: iwghf_eb_av !< indoor average of wghf_eb
437	REAL(wp), DIMENSION(:), ALLOCATABLE :: iwghf_eb_window_av !< indoor average of wghf_eb window
438	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf_eb_av !< average of wshf_eb
439	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_av !< average of wall surface temperature (K)
440	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_window_av !< average of window surface temperature (K)
441	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_green_av !< average of green wall surface temperature (K)
442	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_whole_av !< average of whole wall surface temperature (K)
443	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_10cm_av !< average of the near surface temperature (K)
444
445	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_wall_av !< Average of t_wall
446	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_window_av !< Average of t_window
447	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_green_av !< Average of t_green
448
449	END TYPE surf_type
450
451	TYPE (bc_type), DIMENSION(0:1) :: bc_h !< boundary condition data type, horizontal upward- and downward facing surfaces
452
453	TYPE (surf_type), DIMENSION(0:2), TARGET :: surf_def_h !< horizontal default surfaces (Up, Down, and Top)
454	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_def_v !< vertical default surfaces (North, South, East, West)
455	TYPE (surf_type) , TARGET :: surf_lsm_h !< horizontal natural land surfaces, so far only upward-facing
456	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_lsm_v !< vertical land surfaces (North, South, East, West)
457	TYPE (surf_type) , TARGET :: surf_usm_h !< horizontal urban surfaces, so far only upward-facing
458	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_usm_v !< vertical urban surfaces (North, South, East, West)
459
460	INTEGER(iwp), PARAMETER :: ind_veg_wall = 0 !< index for vegetation / wall-surface fraction, used for access of albedo, emissivity, etc., for each surface type
461	INTEGER(iwp), PARAMETER :: ind_pav_green = 1 !< index for pavement / green-wall surface fraction, used for access of albedo, emissivity, etc., for each surface type
462	INTEGER(iwp), PARAMETER :: ind_wat_win = 2 !< index for water / window-surface fraction, used for access of albedo, emissivity, etc., for each surface type
463
464	INTEGER(iwp) :: ns_h_on_file(0:2) !< total number of horizontal surfaces with the same facing, required for writing restart data
465	INTEGER(iwp) :: ns_v_on_file(0:3) !< total number of vertical surfaces with the same facing, required for writing restart data
466
467	LOGICAL :: vertical_surfaces_exist = .FALSE. !< flag indicating that there are vertical urban/land surfaces
468	!< in the domain (required to activiate RTM)
469
470
471	SAVE
472
473	PRIVATE
474
475	INTERFACE init_bc
476	MODULE PROCEDURE init_bc
477	END INTERFACE init_bc
478
479	INTERFACE init_surfaces
480	MODULE PROCEDURE init_surfaces
481	END INTERFACE init_surfaces
482
483	INTERFACE init_surface_arrays
484	MODULE PROCEDURE init_surface_arrays
485	END INTERFACE init_surface_arrays
486
487	INTERFACE surface_rrd_local
488	MODULE PROCEDURE surface_rrd_local
489	END INTERFACE surface_rrd_local
490
491	INTERFACE surface_wrd_local
492	MODULE PROCEDURE surface_wrd_local
493	END INTERFACE surface_wrd_local
494
495	INTERFACE surface_last_actions
496	MODULE PROCEDURE surface_last_actions
497	END INTERFACE surface_last_actions
498
499	INTERFACE surface_restore_elements
500	MODULE PROCEDURE surface_restore_elements_1d
501	MODULE PROCEDURE surface_restore_elements_2d
502	END INTERFACE surface_restore_elements
503
504	!
505	!-- Public variables
506	PUBLIC bc_h, ind_pav_green, ind_veg_wall, ind_wat_win, ns_h_on_file, &
507	ns_v_on_file, surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, &
508	surf_usm_h, surf_usm_v, surf_type, vertical_surfaces_exist
509	!
510	!-- Public subroutines and functions
511	PUBLIC get_topography_top_index, get_topography_top_index_ji, init_bc, &
512	init_surfaces, &
513	init_surface_arrays, surface_rrd_local, &
514	surface_restore_elements, surface_wrd_local, &
515	surface_last_actions
516
517
518	CONTAINS
519
520	!------------------------------------------------------------------------------!
521	! Description:
522	! ------------
523	!> Initialize data type for setting boundary conditions at horizontal surfaces.
524	!------------------------------------------------------------------------------!
525	SUBROUTINE init_bc
526
527	IMPLICIT NONE
528
529	INTEGER(iwp) :: i !<
530	INTEGER(iwp) :: j !<
531	INTEGER(iwp) :: k !<
532
533	INTEGER(iwp), DIMENSION(0:1) :: num_h !<
534	INTEGER(iwp), DIMENSION(0:1) :: num_h_kji !<
535	INTEGER(iwp), DIMENSION(0:1) :: start_index_h !<
536
537	!
538	!-- First of all, count the number of upward- and downward-facing surfaces
539	num_h = 0
540	DO i = nxlg, nxrg
541	DO j = nysg, nyng
542	DO k = nzb+1, nzt
543	!
544	!-- Check if current gridpoint belongs to the atmosphere
545	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
546	!
547	!-- Upward-facing
548	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) &
549	num_h(0) = num_h(0) + 1
550	!
551	!-- Downward-facing
552	IF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) &
553	num_h(1) = num_h(1) + 1
554	ENDIF
555	ENDDO
556	ENDDO
557	ENDDO
558	!
559	!-- Save the number of surface elements
560	bc_h(0)%ns = num_h(0)
561	bc_h(1)%ns = num_h(1)
562	!
563	!-- ALLOCATE data type variables
564	!-- Upward facing
565	ALLOCATE( bc_h(0)%i(1:bc_h(0)%ns) )
566	ALLOCATE( bc_h(0)%j(1:bc_h(0)%ns) )
567	ALLOCATE( bc_h(0)%k(1:bc_h(0)%ns) )
568	ALLOCATE( bc_h(0)%start_index(nysg:nyng,nxlg:nxrg) )
569	ALLOCATE( bc_h(0)%end_index(nysg:nyng,nxlg:nxrg) )
570	bc_h(0)%start_index = 1
571	bc_h(0)%end_index = 0
572	!
573	!-- Downward facing
574	ALLOCATE( bc_h(1)%i(1:bc_h(1)%ns) )
575	ALLOCATE( bc_h(1)%j(1:bc_h(1)%ns) )
576	ALLOCATE( bc_h(1)%k(1:bc_h(1)%ns) )
577	ALLOCATE( bc_h(1)%start_index(nysg:nyng,nxlg:nxrg) )
578	ALLOCATE( bc_h(1)%end_index(nysg:nyng,nxlg:nxrg) )
579	bc_h(1)%start_index = 1
580	bc_h(1)%end_index = 0
581	!
582	!-- Store the respective indices on data type
583	num_h(0:1) = 1
584	start_index_h(0:1) = 1
585	DO i = nxlg, nxrg
586	DO j = nysg, nyng
587
588	num_h_kji(0:1) = 0
589	DO k = nzb+1, nzt
590	!
591	!-- Check if current gridpoint belongs to the atmosphere
592	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
593	!
594	!-- Upward-facing
595	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
596	bc_h(0)%i(num_h(0)) = i
597	bc_h(0)%j(num_h(0)) = j
598	bc_h(0)%k(num_h(0)) = k
599	num_h_kji(0) = num_h_kji(0) + 1
600	num_h(0) = num_h(0) + 1
601	ENDIF
602	!
603	!-- Downward-facing
604	IF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
605	bc_h(1)%i(num_h(1)) = i
606	bc_h(1)%j(num_h(1)) = j
607	bc_h(1)%k(num_h(1)) = k
608	num_h_kji(1) = num_h_kji(1) + 1
609	num_h(1) = num_h(1) + 1
610	ENDIF
611	ENDIF
612	ENDDO
613	bc_h(0)%start_index(j,i) = start_index_h(0)
614	bc_h(0)%end_index(j,i) = bc_h(0)%start_index(j,i) + num_h_kji(0) - 1
615	start_index_h(0) = bc_h(0)%end_index(j,i) + 1
616
617	bc_h(1)%start_index(j,i) = start_index_h(1)
618	bc_h(1)%end_index(j,i) = bc_h(1)%start_index(j,i) + num_h_kji(1) - 1
619	start_index_h(1) = bc_h(1)%end_index(j,i) + 1
620	ENDDO
621	ENDDO
622
623
624	END SUBROUTINE init_bc
625
626
627	!------------------------------------------------------------------------------!
628	! Description:
629	! ------------
630	!> Initialize horizontal and vertical surfaces. Counts the number of default-,
631	!> natural and urban surfaces and allocates memory, respectively.
632	!------------------------------------------------------------------------------!
633	SUBROUTINE init_surface_arrays
634
635
636	USE pegrid
637
638
639	IMPLICIT NONE
640
641	INTEGER(iwp) :: i !< running index x-direction
642	INTEGER(iwp) :: j !< running index y-direction
643	INTEGER(iwp) :: k !< running index z-direction
644	INTEGER(iwp) :: l !< index variable for surface facing
645	INTEGER(iwp) :: num_lsm_h !< number of horizontally-aligned natural surfaces
646	INTEGER(iwp) :: num_usm_h !< number of horizontally-aligned urban surfaces
647
648	INTEGER(iwp), DIMENSION(0:2) :: num_def_h !< number of horizontally-aligned default surfaces
649	INTEGER(iwp), DIMENSION(0:3) :: num_def_v !< number of vertically-aligned default surfaces
650	INTEGER(iwp), DIMENSION(0:3) :: num_lsm_v !< number of vertically-aligned natural surfaces
651	INTEGER(iwp), DIMENSION(0:3) :: num_usm_v !< number of vertically-aligned urban surfaces
652
653	INTEGER(iwp) :: num_surf_v_l !< number of vertically-aligned local urban/land surfaces
654	INTEGER(iwp) :: num_surf_v !< number of vertically-aligned total urban/land surfaces
655
656	LOGICAL :: building !< flag indicating building grid point
657	LOGICAL :: terrain !< flag indicating natural terrain grid point
658
659	num_def_h = 0
660	num_def_v = 0
661	num_lsm_h = 0
662	num_lsm_v = 0
663	num_usm_h = 0
664	num_usm_v = 0
665	!
666	!-- Surfaces are classified according to the input data read from static
667	!-- input file. If no input file is present, all surfaces are classified
668	!-- either as natural, urban, or default, depending on the setting of
669	!-- land_surface and urban_surface. To control this, use the control
670	!-- flag topo_no_distinct
671	!
672	!-- Count number of horizontal surfaces on local domain
673	DO i = nxl, nxr
674	DO j = nys, nyn
675	DO k = nzb+1, nzt
676	!
677	!-- Check if current gridpoint belongs to the atmosphere
678	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
679	!
680	!-- Check if grid point adjoins to any upward-facing horizontal
681	!-- surface, e.g. the Earth surface, plane roofs, or ceilings.
682
683	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
684	!
685	!-- Determine flags indicating terrain or building.
686	terrain = BTEST( wall_flags_0(k-1,j,i), 5 ) .OR. &
687	topo_no_distinct
688	building = BTEST( wall_flags_0(k-1,j,i), 6 ) .OR. &
689	topo_no_distinct
690	!
691	!-- Land-surface type
692	IF ( land_surface .AND. terrain ) THEN
693	num_lsm_h = num_lsm_h + 1
694	!
695	!-- Urban surface tpye
696	ELSEIF ( urban_surface .AND. building ) THEN
697	num_usm_h = num_usm_h + 1
698	!
699	!-- Default-surface type
700	ELSEIF ( .NOT. land_surface .AND. &
701	.NOT. urban_surface ) THEN
702
703	num_def_h(0) = num_def_h(0) + 1
704	!
705	!-- Unclassifified surface-grid point. Give error message.
706	ELSE
707	WRITE( message_string, * ) &
708	'Unclassified upward-facing ' // &
709	'surface element at '// &
710	'grid point (k,j,i) = ', k, j, i
711	CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 )
712	ENDIF
713
714	ENDIF
715	!
716	!-- Check for top-fluxes
717	IF ( k == nzt .AND. use_top_fluxes ) THEN
718	num_def_h(2) = num_def_h(2) + 1
719	!
720	!-- Check for any other downward-facing surface. So far only for
721	!-- default surface type.
722	ELSEIF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
723	num_def_h(1) = num_def_h(1) + 1
724	ENDIF
725
726	ENDIF
727	ENDDO
728	ENDDO
729	ENDDO
730	!
731	!-- Count number of vertical surfaces on local domain
732	DO i = nxl, nxr
733	DO j = nys, nyn
734	DO k = nzb+1, nzt
735	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
736	!
737	!-- Northward-facing
738	IF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 0 ) ) THEN
739	!
740	!-- Determine flags indicating terrain or building
741
742	terrain = BTEST( wall_flags_0(k,j-1,i), 5 ) .OR. &
743	topo_no_distinct
744	building = BTEST( wall_flags_0(k,j-1,i), 6 ) .OR. &
745	topo_no_distinct
746	IF ( land_surface .AND. terrain ) THEN
747	num_lsm_v(0) = num_lsm_v(0) + 1
748	ELSEIF ( urban_surface .AND. building ) THEN
749	num_usm_v(0) = num_usm_v(0) + 1
750	!
751	!-- Default-surface type
752	ELSEIF ( .NOT. land_surface .AND. &
753	.NOT. urban_surface ) THEN
754	num_def_v(0) = num_def_v(0) + 1
755	!
756	!-- Unclassifified surface-grid point. Give error message.
757	ELSE
758	WRITE( message_string, * ) &
759	'Unclassified northward-facing ' // &
760	'surface element at '// &
761	'grid point (k,j,i) = ', k, j, i
762	CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 )
763
764	ENDIF
765	ENDIF
766	!
767	!-- Southward-facing
768	IF ( .NOT. BTEST( wall_flags_0(k,j+1,i), 0 ) ) THEN
769	!
770	!-- Determine flags indicating terrain or building
771	terrain = BTEST( wall_flags_0(k,j+1,i), 5 ) .OR. &
772	topo_no_distinct
773	building = BTEST( wall_flags_0(k,j+1,i), 6 ) .OR. &
774	topo_no_distinct
775	IF ( land_surface .AND. terrain ) THEN
776	num_lsm_v(1) = num_lsm_v(1) + 1
777	ELSEIF ( urban_surface .AND. building ) THEN
778	num_usm_v(1) = num_usm_v(1) + 1
779	!
780	!-- Default-surface type
781	ELSEIF ( .NOT. land_surface .AND. &
782	.NOT. urban_surface ) THEN
783	num_def_v(1) = num_def_v(1) + 1
784	!
785	!-- Unclassifified surface-grid point. Give error message.
786	ELSE
787	WRITE( message_string, * ) &
788	'Unclassified southward-facing ' // &
789	'surface element at '// &
790	'grid point (k,j,i) = ', k, j, i
791	CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 )
792
793	ENDIF
794	ENDIF
795	!
796	!-- Eastward-facing
797	IF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 0 ) ) THEN
798	!
799	!-- Determine flags indicating terrain or building
800	terrain = BTEST( wall_flags_0(k,j,i-1), 5 ) .OR. &
801	topo_no_distinct
802	building = BTEST( wall_flags_0(k,j,i-1), 6 ) .OR. &
803	topo_no_distinct
804	IF ( land_surface .AND. terrain ) THEN
805	num_lsm_v(2) = num_lsm_v(2) + 1
806	ELSEIF ( urban_surface .AND. building ) THEN
807	num_usm_v(2) = num_usm_v(2) + 1
808	!
809	!-- Default-surface type
810	ELSEIF ( .NOT. land_surface .AND. &
811	.NOT. urban_surface ) THEN
812	num_def_v(2) = num_def_v(2) + 1
813	!
814	!-- Unclassifified surface-grid point. Give error message.
815	ELSE
816	WRITE( message_string, * ) &
817	'Unclassified eastward-facing ' // &
818	'surface element at '// &
819	'grid point (k,j,i) = ', k, j, i
820	CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 )
821
822	ENDIF
823	ENDIF
824	!
825	!-- Westward-facing
826	IF ( .NOT. BTEST( wall_flags_0(k,j,i+1), 0 ) ) THEN
827	!
828	!-- Determine flags indicating terrain or building
829	terrain = BTEST( wall_flags_0(k,j,i+1), 5 ) .OR. &
830	topo_no_distinct
831	building = BTEST( wall_flags_0(k,j,i+1), 6 ) .OR. &
832	topo_no_distinct
833	IF ( land_surface .AND. terrain ) THEN
834	num_lsm_v(3) = num_lsm_v(3) + 1
835	ELSEIF ( urban_surface .AND. building ) THEN
836	num_usm_v(3) = num_usm_v(3) + 1
837	!
838	!-- Default-surface type
839	ELSEIF ( .NOT. land_surface .AND. &
840	.NOT. urban_surface ) THEN
841	num_def_v(3) = num_def_v(3) + 1
842	!
843	!-- Unclassifified surface-grid point. Give error message.
844	ELSE
845	WRITE( message_string, * ) &
846	'Unclassified westward-facing ' // &
847	'surface element at '// &
848	'grid point (k,j,i) = ', k, j, i
849	CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 )
850
851	ENDIF
852	ENDIF
853	ENDIF
854	ENDDO
855	ENDDO
856	ENDDO
857
858	!
859	!-- Store number of surfaces per core.
860	!-- Horizontal surface, default type, upward facing
861	surf_def_h(0)%ns = num_def_h(0)
862	!
863	!-- Horizontal surface, default type, downward facing
864	surf_def_h(1)%ns = num_def_h(1)
865	!
866	!-- Horizontal surface, default type, top downward facing
867	surf_def_h(2)%ns = num_def_h(2)
868	!
869	!-- Horizontal surface, natural type, so far only upward-facing
870	surf_lsm_h%ns = num_lsm_h
871	!
872	!-- Horizontal surface, urban type, so far only upward-facing
873	surf_usm_h%ns = num_usm_h
874	!
875	!-- Vertical surface, default type, northward facing
876	surf_def_v(0)%ns = num_def_v(0)
877	!
878	!-- Vertical surface, default type, southward facing
879	surf_def_v(1)%ns = num_def_v(1)
880	!
881	!-- Vertical surface, default type, eastward facing
882	surf_def_v(2)%ns = num_def_v(2)
883	!
884	!-- Vertical surface, default type, westward facing
885	surf_def_v(3)%ns = num_def_v(3)
886	!
887	!-- Vertical surface, natural type, northward facing
888	surf_lsm_v(0)%ns = num_lsm_v(0)
889	!
890	!-- Vertical surface, natural type, southward facing
891	surf_lsm_v(1)%ns = num_lsm_v(1)
892	!
893	!-- Vertical surface, natural type, eastward facing
894	surf_lsm_v(2)%ns = num_lsm_v(2)
895	!
896	!-- Vertical surface, natural type, westward facing
897	surf_lsm_v(3)%ns = num_lsm_v(3)
898	!
899	!-- Vertical surface, urban type, northward facing
900	surf_usm_v(0)%ns = num_usm_v(0)
901	!
902	!-- Vertical surface, urban type, southward facing
903	surf_usm_v(1)%ns = num_usm_v(1)
904	!
905	!-- Vertical surface, urban type, eastward facing
906	surf_usm_v(2)%ns = num_usm_v(2)
907	!
908	!-- Vertical surface, urban type, westward facing
909	surf_usm_v(3)%ns = num_usm_v(3)
910	!
911	!-- Allocate required attributes for horizontal surfaces - default type.
912	!-- Upward-facing (l=0) and downward-facing (l=1).
913	DO l = 0, 1
914	CALL allocate_surface_attributes_h ( surf_def_h(l), nys, nyn, nxl, nxr )
915	ENDDO
916	!
917	!-- Allocate required attributes for model top
918	CALL allocate_surface_attributes_h_top ( surf_def_h(2), nys, nyn, nxl, nxr )
919	!
920	!-- Allocate required attributes for horizontal surfaces - natural type.
921	CALL allocate_surface_attributes_h ( surf_lsm_h, nys, nyn, nxl, nxr )
922	!
923	!-- Allocate required attributes for horizontal surfaces - urban type.
924	CALL allocate_surface_attributes_h ( surf_usm_h, nys, nyn, nxl, nxr )
925
926	!
927	!-- Allocate required attributes for vertical surfaces.
928	!-- Northward-facing (l=0), southward-facing (l=1), eastward-facing (l=2)
929	!-- and westward-facing (l=3).
930	!-- Default type.
931	DO l = 0, 3
932	CALL allocate_surface_attributes_v ( surf_def_v(l), &
933	nys, nyn, nxl, nxr )
934	ENDDO
935	!
936	!-- Natural type
937	DO l = 0, 3
938	CALL allocate_surface_attributes_v ( surf_lsm_v(l), &
939	nys, nyn, nxl, nxr )
940	ENDDO
941	!
942	!-- Urban type
943	DO l = 0, 3
944	CALL allocate_surface_attributes_v ( surf_usm_v(l), &
945	nys, nyn, nxl, nxr )
946	ENDDO
947	!
948	!-- Set the flag for the existence of vertical urban/land surfaces
949	num_surf_v_l = 0
950	DO l = 0, 3
951	num_surf_v_l = num_surf_v_l + surf_usm_v(l)%ns + surf_lsm_v(l)%ns
952	ENDDO
953
954	#if defined( __parallel )
955	CALL MPI_ALLREDUCE( num_surf_v_l, num_surf_v, 1, MPI_INTEGER, MPI_SUM, comm2d, ierr)
956	#else
957	num_surf_v = num_surf_v_l
958	#endif
959	IF ( num_surf_v > 0 ) vertical_surfaces_exist = .TRUE.
960
961
962	END SUBROUTINE init_surface_arrays
963
964
965	!------------------------------------------------------------------------------!
966	! Description:
967	! ------------
968	!> Deallocating memory for upward and downward-facing horizontal surface types,
969	!> except for top fluxes.
970	!------------------------------------------------------------------------------!
971	SUBROUTINE deallocate_surface_attributes_h( surfaces )
972
973	IMPLICIT NONE
974
975
976	TYPE(surf_type) :: surfaces !< respective surface type
977
978
979	DEALLOCATE ( surfaces%start_index )
980	DEALLOCATE ( surfaces%end_index )
981	!
982	!-- Indices to locate surface element
983	DEALLOCATE ( surfaces%i )
984	DEALLOCATE ( surfaces%j )
985	DEALLOCATE ( surfaces%k )
986	!
987	!-- Surface-layer height
988	DEALLOCATE ( surfaces%z_mo )
989	!
990	!-- Surface orientation
991	DEALLOCATE ( surfaces%facing )
992	!
993	!-- Surface-parallel wind velocity
994	DEALLOCATE ( surfaces%uvw_abs )
995	!
996	!-- Roughness
997	DEALLOCATE ( surfaces%z0 )
998	DEALLOCATE ( surfaces%z0h )
999	DEALLOCATE ( surfaces%z0q )
1000	!
1001	!-- Friction velocity
1002	DEALLOCATE ( surfaces%us )
1003	!
1004	!-- Stability parameter
1005	DEALLOCATE ( surfaces%ol )
1006	!
1007	!-- Bulk Richardson number
1008	DEALLOCATE ( surfaces%rib )
1009	!
1010	!-- Vertical momentum fluxes of u and v
1011	DEALLOCATE ( surfaces%usws )
1012	DEALLOCATE ( surfaces%vsws )
1013	!
1014	!-- Required in production_e
1015	IF ( .NOT. constant_diffusion ) THEN
1016	DEALLOCATE ( surfaces%u_0 )
1017	DEALLOCATE ( surfaces%v_0 )
1018	ENDIF
1019	!
1020	!-- Characteristic temperature and surface flux of sensible heat
1021	DEALLOCATE ( surfaces%ts )
1022	DEALLOCATE ( surfaces%shf )
1023	!
1024	!-- surface temperature
1025	DEALLOCATE ( surfaces%pt_surface )
1026	!
1027	!-- Characteristic humidity and surface flux of latent heat
1028	IF ( humidity ) THEN
1029	DEALLOCATE ( surfaces%qs )
1030	DEALLOCATE ( surfaces%qsws )
1031	ENDIF
1032	!
1033	!-- Characteristic scalar and surface flux of scalar
1034	IF ( passive_scalar ) THEN
1035	DEALLOCATE ( surfaces%ss )
1036	DEALLOCATE ( surfaces%ssws )
1037	ENDIF
1038	!
1039	!-- Scaling parameter (cs*) and surface flux of chemical species
1040	IF ( air_chemistry ) THEN
1041	DEALLOCATE ( surfaces%css )
1042	DEALLOCATE ( surfaces%cssws )
1043	ENDIF
1044	!
1045	!-- Arrays for storing potential temperature and
1046	!-- mixing ratio at first grid level
1047	DEALLOCATE ( surfaces%pt1 )
1048	DEALLOCATE ( surfaces%qv1 )
1049	!
1050	!--
1051	IF ( cloud_physics .AND. microphysics_morrison) THEN
1052	DEALLOCATE ( surfaces%qcs )
1053	DEALLOCATE ( surfaces%ncs )
1054	DEALLOCATE ( surfaces%qcsws )
1055	DEALLOCATE ( surfaces%ncsws )
1056	ENDIF
1057	!
1058	!--
1059	IF ( cloud_physics .AND. microphysics_seifert) THEN
1060	DEALLOCATE ( surfaces%qrs )
1061	DEALLOCATE ( surfaces%nrs )
1062	DEALLOCATE ( surfaces%qrsws )
1063	DEALLOCATE ( surfaces%nrsws )
1064	ENDIF
1065	!
1066	!-- Salinity surface flux
1067	IF ( ocean ) DEALLOCATE ( surfaces%sasws )
1068
1069	END SUBROUTINE deallocate_surface_attributes_h
1070
1071
1072	!------------------------------------------------------------------------------!
1073	! Description:
1074	! ------------
1075	!> Allocating memory for upward and downward-facing horizontal surface types,
1076	!> except for top fluxes.
1077	!------------------------------------------------------------------------------!
1078	SUBROUTINE allocate_surface_attributes_h( surfaces, &
1079	nys_l, nyn_l, nxl_l, nxr_l )
1080
1081	IMPLICIT NONE
1082
1083	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1084	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1085	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1086	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1087
1088	TYPE(surf_type) :: surfaces !< respective surface type
1089
1090	!
1091	!-- Allocate arrays for start and end index of horizontal surface type
1092	!-- for each (j,i)-grid point. This is required e.g. in diffion_x, which is
1093	!-- called for each (j,i). In order to find the location where the
1094	!-- respective flux is store within the surface-type, start- and end-
1095	!-- index are stored for each (j,i). For example, each (j,i) can have
1096	!-- several entries where fluxes for horizontal surfaces might be stored,
1097	!-- e.g. for overhanging structures where several upward-facing surfaces
1098	!-- might exist for given (j,i).
1099	!-- If no surface of respective type exist at current (j,i), set indicies
1100	!-- such that loop in diffusion routines will not be entered.
1101	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
1102	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
1103	surfaces%start_index = 0
1104	surfaces%end_index = -1
1105	!
1106	!-- Indices to locate surface element
1107	ALLOCATE ( surfaces%i(1:surfaces%ns) )
1108	ALLOCATE ( surfaces%j(1:surfaces%ns) )
1109	ALLOCATE ( surfaces%k(1:surfaces%ns) )
1110	!
1111	!-- Surface-layer height
1112	ALLOCATE ( surfaces%z_mo(1:surfaces%ns) )
1113	!
1114	!-- Surface orientation
1115	ALLOCATE ( surfaces%facing(1:surfaces%ns) )
1116	!
1117	!-- Surface-parallel wind velocity
1118	ALLOCATE ( surfaces%uvw_abs(1:surfaces%ns) )
1119	!
1120	!-- Roughness
1121	ALLOCATE ( surfaces%z0(1:surfaces%ns) )
1122	ALLOCATE ( surfaces%z0h(1:surfaces%ns) )
1123	ALLOCATE ( surfaces%z0q(1:surfaces%ns) )
1124	!
1125	!-- Friction velocity
1126	ALLOCATE ( surfaces%us(1:surfaces%ns) )
1127	!
1128	!-- Stability parameter
1129	ALLOCATE ( surfaces%ol(1:surfaces%ns) )
1130	!
1131	!-- Bulk Richardson number
1132	ALLOCATE ( surfaces%rib(1:surfaces%ns) )
1133	!
1134	!-- Vertical momentum fluxes of u and v
1135	ALLOCATE ( surfaces%usws(1:surfaces%ns) )
1136	ALLOCATE ( surfaces%vsws(1:surfaces%ns) )
1137	!
1138	!-- Required in production_e
1139	IF ( .NOT. constant_diffusion ) THEN
1140	ALLOCATE ( surfaces%u_0(1:surfaces%ns) )
1141	ALLOCATE ( surfaces%v_0(1:surfaces%ns) )
1142	ENDIF
1143	!
1144	!-- Characteristic temperature and surface flux of sensible heat
1145	ALLOCATE ( surfaces%ts(1:surfaces%ns) )
1146	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
1147	!
1148	!-- surface temperature
1149	ALLOCATE ( surfaces%pt_surface(1:surfaces%ns) )
1150	!
1151	!-- Characteristic humidity and surface flux of latent heat
1152	IF ( humidity ) THEN
1153	ALLOCATE ( surfaces%qs(1:surfaces%ns) )
1154	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
1155	ENDIF
1156	!
1157	!-- Characteristic scalar and surface flux of scalar
1158	IF ( passive_scalar ) THEN
1159	ALLOCATE ( surfaces%ss(1:surfaces%ns) )
1160	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
1161	ENDIF
1162	!
1163	!-- Scaling parameter (cs*) and surface flux of chemical species
1164	IF ( air_chemistry ) THEN
1165	ALLOCATE ( surfaces%css(1:nvar,1:surfaces%ns) )
1166	ALLOCATE ( surfaces%cssws(1:nvar,1:surfaces%ns) )
1167	ENDIF
1168	!
1169	!-- Arrays for storing potential temperature and
1170	!-- mixing ratio at first grid level
1171	ALLOCATE ( surfaces%pt1(1:surfaces%ns) )
1172	ALLOCATE ( surfaces%qv1(1:surfaces%ns) )
1173	!
1174	!--
1175	IF ( cloud_physics .AND. microphysics_morrison) THEN
1176	ALLOCATE ( surfaces%qcs(1:surfaces%ns) )
1177	ALLOCATE ( surfaces%ncs(1:surfaces%ns) )
1178	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
1179	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
1180	ENDIF
1181	!
1182	!--
1183	IF ( cloud_physics .AND. microphysics_seifert) THEN
1184	ALLOCATE ( surfaces%qrs(1:surfaces%ns) )
1185	ALLOCATE ( surfaces%nrs(1:surfaces%ns) )
1186	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
1187	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
1188	ENDIF
1189	!
1190	!-- Salinity surface flux
1191	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
1192
1193	END SUBROUTINE allocate_surface_attributes_h
1194
1195
1196	!------------------------------------------------------------------------------!
1197	! Description:
1198	! ------------
1199	!> Deallocating memory for model-top fluxes
1200	!------------------------------------------------------------------------------!
1201	SUBROUTINE deallocate_surface_attributes_h_top( surfaces )
1202
1203	IMPLICIT NONE
1204
1205
1206	TYPE(surf_type) :: surfaces !< respective surface type
1207
1208	DEALLOCATE ( surfaces%start_index )
1209	DEALLOCATE ( surfaces%end_index )
1210	!
1211	!-- Indices to locate surface (model-top) element
1212	DEALLOCATE ( surfaces%i )
1213	DEALLOCATE ( surfaces%j )
1214	DEALLOCATE ( surfaces%k )
1215
1216	IF ( .NOT. constant_diffusion ) THEN
1217	DEALLOCATE ( surfaces%u_0 )
1218	DEALLOCATE ( surfaces%v_0 )
1219	ENDIF
1220	!
1221	!-- Vertical momentum fluxes of u and v
1222	DEALLOCATE ( surfaces%usws )
1223	DEALLOCATE ( surfaces%vsws )
1224	!
1225	!-- Sensible heat flux
1226	DEALLOCATE ( surfaces%shf )
1227	!
1228	!-- Latent heat flux
1229	IF ( humidity .OR. coupling_mode == 'ocean_to_atmosphere') THEN
1230	DEALLOCATE ( surfaces%qsws )
1231	ENDIF
1232	!
1233	!-- Scalar flux
1234	IF ( passive_scalar ) THEN
1235	DEALLOCATE ( surfaces%ssws )
1236	ENDIF
1237	!
1238	!-- Chemical species flux
1239	IF ( air_chemistry ) THEN
1240	DEALLOCATE ( surfaces%cssws )
1241	ENDIF
1242	!
1243	!--
1244	IF ( cloud_physics .AND. microphysics_morrison) THEN
1245	DEALLOCATE ( surfaces%qcsws )
1246	DEALLOCATE ( surfaces%ncsws )
1247	ENDIF
1248	!
1249	!--
1250	IF ( cloud_physics .AND. microphysics_seifert) THEN
1251	DEALLOCATE ( surfaces%qrsws )
1252	DEALLOCATE ( surfaces%nrsws )
1253	ENDIF
1254	!
1255	!-- Salinity flux
1256	IF ( ocean ) DEALLOCATE ( surfaces%sasws )
1257
1258	END SUBROUTINE deallocate_surface_attributes_h_top
1259
1260
1261	!------------------------------------------------------------------------------!
1262	! Description:
1263	! ------------
1264	!> Allocating memory for model-top fluxes
1265	!------------------------------------------------------------------------------!
1266	SUBROUTINE allocate_surface_attributes_h_top( surfaces, &
1267	nys_l, nyn_l, nxl_l, nxr_l )
1268
1269	IMPLICIT NONE
1270
1271	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1272	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1273	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1274	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1275
1276	TYPE(surf_type) :: surfaces !< respective surface type
1277
1278	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
1279	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
1280	surfaces%start_index = 0
1281	surfaces%end_index = -1
1282	!
1283	!-- Indices to locate surface (model-top) element
1284	ALLOCATE ( surfaces%i(1:surfaces%ns) )
1285	ALLOCATE ( surfaces%j(1:surfaces%ns) )
1286	ALLOCATE ( surfaces%k(1:surfaces%ns) )
1287
1288	IF ( .NOT. constant_diffusion ) THEN
1289	ALLOCATE ( surfaces%u_0(1:surfaces%ns) )
1290	ALLOCATE ( surfaces%v_0(1:surfaces%ns) )
1291	ENDIF
1292	!
1293	!-- Vertical momentum fluxes of u and v
1294	ALLOCATE ( surfaces%usws(1:surfaces%ns) )
1295	ALLOCATE ( surfaces%vsws(1:surfaces%ns) )
1296	!
1297	!-- Sensible heat flux
1298	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
1299	!
1300	!-- Latent heat flux
1301	IF ( humidity .OR. coupling_mode == 'ocean_to_atmosphere') THEN
1302	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
1303	ENDIF
1304	!
1305	!-- Scalar flux
1306	IF ( passive_scalar ) THEN
1307	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
1308	ENDIF
1309	!
1310	!-- Chemical species flux
1311	IF ( air_chemistry ) THEN
1312	ALLOCATE ( surfaces%cssws(1:nvar,1:surfaces%ns) )
1313	ENDIF
1314	!
1315	!--
1316	IF ( cloud_physics .AND. microphysics_morrison) THEN
1317	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
1318	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
1319	ENDIF
1320	!
1321	!--
1322	IF ( cloud_physics .AND. microphysics_seifert) THEN
1323	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
1324	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
1325	ENDIF
1326	!
1327	!-- Salinity flux
1328	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
1329
1330	END SUBROUTINE allocate_surface_attributes_h_top
1331
1332
1333	!------------------------------------------------------------------------------!
1334	! Description:
1335	! ------------
1336	!> Deallocating memory for vertical surface types.
1337	!------------------------------------------------------------------------------!
1338	SUBROUTINE deallocate_surface_attributes_v( surfaces )
1339
1340	IMPLICIT NONE
1341
1342
1343	TYPE(surf_type) :: surfaces !< respective surface type
1344
1345	!
1346	!-- Allocate arrays for start and end index of vertical surface type
1347	!-- for each (j,i)-grid point. This is required in diffion_x, which is
1348	!-- called for each (j,i). In order to find the location where the
1349	!-- respective flux is store within the surface-type, start- and end-
1350	!-- index are stored for each (j,i). For example, each (j,i) can have
1351	!-- several entries where fluxes for vertical surfaces might be stored.
1352	!-- In the flat case, where no vertical walls exit, set indicies such
1353	!-- that loop in diffusion routines will not be entered.
1354	DEALLOCATE ( surfaces%start_index )
1355	DEALLOCATE ( surfaces%end_index )
1356	!
1357	!-- Indices to locate surface element.
1358	DEALLOCATE ( surfaces%i )
1359	DEALLOCATE ( surfaces%j )
1360	DEALLOCATE ( surfaces%k )
1361	!
1362	!-- Surface-layer height
1363	DEALLOCATE ( surfaces%z_mo )
1364	!
1365	!-- Surface orientation
1366	DEALLOCATE ( surfaces%facing )
1367	!
1368	!-- Surface parallel wind velocity
1369	DEALLOCATE ( surfaces%uvw_abs )
1370	!
1371	!-- Roughness
1372	DEALLOCATE ( surfaces%z0 )
1373	DEALLOCATE ( surfaces%z0h )
1374	DEALLOCATE ( surfaces%z0q )
1375
1376	!
1377	!-- Friction velocity
1378	DEALLOCATE ( surfaces%us )
1379	!
1380	!-- Allocate Obukhov length and bulk Richardson number. Actually, at
1381	!-- vertical surfaces these are only required for natural surfaces.
1382	!-- for natural land surfaces
1383	DEALLOCATE( surfaces%ol )
1384	DEALLOCATE( surfaces%rib )
1385	!
1386	!-- Allocate arrays for surface momentum fluxes for u and v. For u at north-
1387	!-- and south-facing surfaces, for v at east- and west-facing surfaces.
1388	DEALLOCATE ( surfaces%mom_flux_uv )
1389	!
1390	!-- Allocate array for surface momentum flux for w - wsus and wsvs
1391	DEALLOCATE ( surfaces%mom_flux_w )
1392	!
1393	!-- Allocate array for surface momentum flux for subgrid-scale tke wsus and
1394	!-- wsvs; first index usvs or vsws, second index for wsus or wsvs, depending
1395	!-- on surface.
1396	DEALLOCATE ( surfaces%mom_flux_tke )
1397	!
1398	!-- Characteristic temperature and surface flux of sensible heat
1399	DEALLOCATE ( surfaces%ts )
1400	DEALLOCATE ( surfaces%shf )
1401	!
1402	!-- surface temperature
1403	DEALLOCATE ( surfaces%pt_surface )
1404	!
1405	!-- Characteristic humidity and surface flux of latent heat
1406	IF ( humidity ) THEN
1407	DEALLOCATE ( surfaces%qs )
1408	DEALLOCATE ( surfaces%qsws )
1409	ENDIF
1410	!
1411	!-- Characteristic scalar and surface flux of scalar
1412	IF ( passive_scalar ) THEN
1413	DEALLOCATE ( surfaces%ss )
1414	DEALLOCATE ( surfaces%ssws )
1415	ENDIF
1416	!
1417	!-- Scaling parameter (cs*) and surface flux of chemical species
1418	IF ( air_chemistry ) THEN
1419	DEALLOCATE ( surfaces%css )
1420	DEALLOCATE ( surfaces%cssws )
1421	ENDIF
1422	!
1423	!-- Arrays for storing potential temperature and
1424	!-- mixing ratio at first grid level
1425	DEALLOCATE ( surfaces%pt1 )
1426	DEALLOCATE ( surfaces%qv1 )
1427
1428	IF ( cloud_physics .AND. microphysics_seifert) THEN
1429	DEALLOCATE ( surfaces%qcs )
1430	DEALLOCATE ( surfaces%ncs )
1431	DEALLOCATE ( surfaces%qcsws )
1432	DEALLOCATE ( surfaces%ncsws )
1433	ENDIF
1434
1435	IF ( cloud_physics .AND. microphysics_seifert) THEN
1436	DEALLOCATE ( surfaces%qrs )
1437	DEALLOCATE ( surfaces%nrs )
1438	DEALLOCATE ( surfaces%qrsws )
1439	DEALLOCATE ( surfaces%nrsws )
1440	ENDIF
1441	!
1442	!-- Salinity surface flux
1443	IF ( ocean ) DEALLOCATE ( surfaces%sasws )
1444
1445	END SUBROUTINE deallocate_surface_attributes_v
1446
1447
1448	!------------------------------------------------------------------------------!
1449	! Description:
1450	! ------------
1451	!> Allocating memory for vertical surface types.
1452	!------------------------------------------------------------------------------!
1453	SUBROUTINE allocate_surface_attributes_v( surfaces, &
1454	nys_l, nyn_l, nxl_l, nxr_l )
1455
1456	IMPLICIT NONE
1457
1458	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1459	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1460	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1461	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
1462
1463	TYPE(surf_type) :: surfaces !< respective surface type
1464
1465	!
1466	!-- Allocate arrays for start and end index of vertical surface type
1467	!-- for each (j,i)-grid point. This is required in diffion_x, which is
1468	!-- called for each (j,i). In order to find the location where the
1469	!-- respective flux is store within the surface-type, start- and end-
1470	!-- index are stored for each (j,i). For example, each (j,i) can have
1471	!-- several entries where fluxes for vertical surfaces might be stored.
1472	!-- In the flat case, where no vertical walls exit, set indicies such
1473	!-- that loop in diffusion routines will not be entered.
1474	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
1475	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
1476	surfaces%start_index = 0
1477	surfaces%end_index = -1
1478	!
1479	!-- Indices to locate surface element.
1480	ALLOCATE ( surfaces%i(1:surfaces%ns) )
1481	ALLOCATE ( surfaces%j(1:surfaces%ns) )
1482	ALLOCATE ( surfaces%k(1:surfaces%ns) )
1483	!
1484	!-- Surface-layer height
1485	ALLOCATE ( surfaces%z_mo(1:surfaces%ns) )
1486	!
1487	!-- Surface orientation
1488	ALLOCATE ( surfaces%facing(1:surfaces%ns) )
1489	!
1490	!-- Surface parallel wind velocity
1491	ALLOCATE ( surfaces%uvw_abs(1:surfaces%ns) )
1492	!
1493	!-- Roughness
1494	ALLOCATE ( surfaces%z0(1:surfaces%ns) )
1495	ALLOCATE ( surfaces%z0h(1:surfaces%ns) )
1496	ALLOCATE ( surfaces%z0q(1:surfaces%ns) )
1497
1498	!
1499	!-- Friction velocity
1500	ALLOCATE ( surfaces%us(1:surfaces%ns) )
1501	!
1502	!-- Allocate Obukhov length and bulk Richardson number. Actually, at
1503	!-- vertical surfaces these are only required for natural surfaces.
1504	!-- for natural land surfaces
1505	ALLOCATE( surfaces%ol(1:surfaces%ns) )
1506	ALLOCATE( surfaces%rib(1:surfaces%ns) )
1507	!
1508	!-- Allocate arrays for surface momentum fluxes for u and v. For u at north-
1509	!-- and south-facing surfaces, for v at east- and west-facing surfaces.
1510	ALLOCATE ( surfaces%mom_flux_uv(1:surfaces%ns) )
1511	!
1512	!-- Allocate array for surface momentum flux for w - wsus and wsvs
1513	ALLOCATE ( surfaces%mom_flux_w(1:surfaces%ns) )
1514	!
1515	!-- Allocate array for surface momentum flux for subgrid-scale tke wsus and
1516	!-- wsvs; first index usvs or vsws, second index for wsus or wsvs, depending
1517	!-- on surface.
1518	ALLOCATE ( surfaces%mom_flux_tke(0:1,1:surfaces%ns) )
1519	!
1520	!-- Characteristic temperature and surface flux of sensible heat
1521	ALLOCATE ( surfaces%ts(1:surfaces%ns) )
1522	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
1523	!
1524	!-- surface temperature
1525	ALLOCATE ( surfaces%pt_surface(1:surfaces%ns) )
1526	!
1527	!-- Characteristic humidity and surface flux of latent heat
1528	IF ( humidity ) THEN
1529	ALLOCATE ( surfaces%qs(1:surfaces%ns) )
1530	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
1531	ENDIF
1532	!
1533	!-- Characteristic scalar and surface flux of scalar
1534	IF ( passive_scalar ) THEN
1535	ALLOCATE ( surfaces%ss(1:surfaces%ns) )
1536	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
1537	ENDIF
1538	!
1539	!-- Scaling parameter (cs*) and surface flux of chemical species
1540	IF ( air_chemistry ) THEN
1541	ALLOCATE ( surfaces%css(1:nvar,1:surfaces%ns) )
1542	ALLOCATE ( surfaces%cssws(1:nvar,1:surfaces%ns) )
1543	ENDIF
1544	!
1545	!-- Arrays for storing potential temperature and
1546	!-- mixing ratio at first grid level
1547	ALLOCATE ( surfaces%pt1(1:surfaces%ns) )
1548	ALLOCATE ( surfaces%qv1(1:surfaces%ns) )
1549
1550	IF ( cloud_physics .AND. microphysics_seifert) THEN
1551	ALLOCATE ( surfaces%qcs(1:surfaces%ns) )
1552	ALLOCATE ( surfaces%ncs(1:surfaces%ns) )
1553	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
1554	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
1555	ENDIF
1556
1557	IF ( cloud_physics .AND. microphysics_seifert) THEN
1558	ALLOCATE ( surfaces%qrs(1:surfaces%ns) )
1559	ALLOCATE ( surfaces%nrs(1:surfaces%ns) )
1560	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
1561	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
1562	ENDIF
1563	!
1564	!-- Salinity surface flux
1565	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
1566
1567	END SUBROUTINE allocate_surface_attributes_v
1568
1569
1570	!------------------------------------------------------------------------------!
1571	! Description:
1572	! ------------
1573	!> Initialize surface elements, i.e. set initial values for surface fluxes,
1574	!> friction velocity, calcuation of start/end indices, etc. .
1575	!> Please note, further initialization concerning
1576	!> special surface characteristics, e.g. soil- and vegatation type,
1577	!> building type, etc., is done in the land-surface and urban-surface module,
1578	!> respectively.
1579	!------------------------------------------------------------------------------!
1580	SUBROUTINE init_surfaces
1581
1582	IMPLICIT NONE
1583
1584	INTEGER(iwp) :: i !< running index x-direction
1585	INTEGER(iwp) :: j !< running index y-direction
1586	INTEGER(iwp) :: k !< running index z-direction
1587	INTEGER(iwp) :: l !< index variable used to distinguish surface facing
1588	INTEGER(iwp) :: m !< running index surface elements
1589
1590	INTEGER(iwp) :: start_index_lsm_h !< dummy to determing local start index in surface type for given (j,i), for horizontal natural surfaces
1591	INTEGER(iwp) :: start_index_usm_h !< dummy to determing local start index in surface type for given (j,i), for horizontal urban surfaces
1592
1593	INTEGER(iwp) :: num_lsm_h !< current number of horizontal surface element, natural type
1594	INTEGER(iwp) :: num_lsm_h_kji !< dummy to determing local end index in surface type for given (j,i), for for horizonal natural surfaces
1595	INTEGER(iwp) :: num_usm_h !< current number of horizontal surface element, urban type
1596	INTEGER(iwp) :: num_usm_h_kji !< dummy to determing local end index in surface type for given (j,i), for for horizonal urban surfaces
1597
1598	INTEGER(iwp), DIMENSION(0:2) :: num_def_h !< current number of horizontal surface element, default type
1599	INTEGER(iwp), DIMENSION(0:2) :: num_def_h_kji !< dummy to determing local end index in surface type for given (j,i), for horizonal default surfaces
1600	INTEGER(iwp), DIMENSION(0:2) :: start_index_def_h !< dummy to determing local start index in surface type for given (j,i), for horizontal default surfaces
1601
1602	INTEGER(iwp), DIMENSION(0:3) :: num_def_v !< current number of vertical surface element, default type
1603	INTEGER(iwp), DIMENSION(0:3) :: num_def_v_kji !< dummy to determing local end index in surface type for given (j,i), for vertical default surfaces
1604	INTEGER(iwp), DIMENSION(0:3) :: num_lsm_v !< current number of vertical surface element, natural type
1605	INTEGER(iwp), DIMENSION(0:3) :: num_lsm_v_kji !< dummy to determing local end index in surface type for given (j,i), for vertical natural surfaces
1606	INTEGER(iwp), DIMENSION(0:3) :: num_usm_v !< current number of vertical surface element, urban type
1607	INTEGER(iwp), DIMENSION(0:3) :: num_usm_v_kji !< dummy to determing local end index in surface type for given (j,i), for vertical urban surfaces
1608
1609	INTEGER(iwp), DIMENSION(0:3) :: start_index_def_v !< dummy to determing local start index in surface type for given (j,i), for vertical default surfaces
1610	INTEGER(iwp), DIMENSION(0:3) :: start_index_lsm_v !< dummy to determing local start index in surface type for given (j,i), for vertical natural surfaces
1611	INTEGER(iwp), DIMENSION(0:3) :: start_index_usm_v !< dummy to determing local start index in surface type for given (j,i), for vertical urban surfaces
1612
1613	LOGICAL :: building !< flag indicating building grid point
1614	LOGICAL :: terrain !< flag indicating natural terrain grid point
1615
1616	!
1617	!-- Set offset indices, i.e. index difference between surface element and
1618	!-- surface-bounded grid point.
1619	!-- Upward facing - no horizontal offsets
1620	surf_def_h(0:2)%ioff = 0
1621	surf_def_h(0:2)%joff = 0
1622
1623	surf_lsm_h%ioff = 0
1624	surf_lsm_h%joff = 0
1625
1626	surf_usm_h%ioff = 0
1627	surf_usm_h%joff = 0
1628	!
1629	!-- Upward facing vertical offsets
1630	surf_def_h(0)%koff = -1
1631	surf_lsm_h%koff = -1
1632	surf_usm_h%koff = -1
1633	!
1634	!-- Downward facing vertical offset
1635	surf_def_h(1:2)%koff = 1
1636	!
1637	!-- Vertical surfaces - no vertical offset
1638	surf_def_v(0:3)%koff = 0
1639	surf_lsm_v(0:3)%koff = 0
1640	surf_usm_v(0:3)%koff = 0
1641	!
1642	!-- North- and southward facing - no offset in x
1643	surf_def_v(0:1)%ioff = 0
1644	surf_lsm_v(0:1)%ioff = 0
1645	surf_usm_v(0:1)%ioff = 0
1646	!
1647	!-- Northward facing offset in y
1648	surf_def_v(0)%joff = -1
1649	surf_lsm_v(0)%joff = -1
1650	surf_usm_v(0)%joff = -1
1651	!
1652	!-- Southward facing offset in y
1653	surf_def_v(1)%joff = 1
1654	surf_lsm_v(1)%joff = 1
1655	surf_usm_v(1)%joff = 1
1656
1657	!
1658	!-- East- and westward facing - no offset in y
1659	surf_def_v(2:3)%joff = 0
1660	surf_lsm_v(2:3)%joff = 0
1661	surf_usm_v(2:3)%joff = 0
1662	!
1663	!-- Eastward facing offset in x
1664	surf_def_v(2)%ioff = -1
1665	surf_lsm_v(2)%ioff = -1
1666	surf_usm_v(2)%ioff = -1
1667	!
1668	!-- Westward facing offset in y
1669	surf_def_v(3)%ioff = 1
1670	surf_lsm_v(3)%ioff = 1
1671	surf_usm_v(3)%ioff = 1
1672
1673	!
1674	!-- Initialize surface attributes, store indicies, surfaces orientation, etc.,
1675	num_def_h(0:2) = 1
1676	num_def_v(0:3) = 1
1677
1678	num_lsm_h = 1
1679	num_lsm_v(0:3) = 1
1680
1681	num_usm_h = 1
1682	num_usm_v(0:3) = 1
1683
1684	start_index_def_h(0:2) = 1
1685	start_index_def_v(0:3) = 1
1686
1687	start_index_lsm_h = 1
1688	start_index_lsm_v(0:3) = 1
1689
1690	start_index_usm_h = 1
1691	start_index_usm_v(0:3) = 1
1692
1693	DO i = nxl, nxr
1694	DO j = nys, nyn
1695
1696	num_def_h_kji = 0
1697	num_def_v_kji = 0
1698	num_lsm_h_kji = 0
1699	num_lsm_v_kji = 0
1700	num_usm_h_kji = 0
1701	num_usm_v_kji = 0
1702
1703	DO k = nzb+1, nzt
1704	!
1705	!-- Check if current gridpoint belongs to the atmosphere
1706	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
1707	!
1708	!-- Upward-facing surface. Distinguish between differet surface types.
1709	!-- To do, think about method to flag natural and non-natural
1710	!-- surfaces.
1711	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
1712	!
1713	!-- Determine flags indicating terrain or building
1714	terrain = BTEST( wall_flags_0(k-1,j,i), 5 ) .OR. &
1715	topo_no_distinct
1716	building = BTEST( wall_flags_0(k-1,j,i), 6 ) .OR. &
1717	topo_no_distinct
1718	!
1719	!-- Natural surface type
1720	IF ( land_surface .AND. terrain ) THEN
1721	CALL initialize_horizontal_surfaces( k, j, i, &
1722	surf_lsm_h, &
1723	num_lsm_h, &
1724	num_lsm_h_kji, &
1725	.TRUE., .FALSE. )
1726	!
1727	!-- Urban surface tpye
1728	ELSEIF ( urban_surface .AND. building ) THEN
1729	CALL initialize_horizontal_surfaces( k, j, i, &
1730	surf_usm_h, &
1731	num_usm_h, &
1732	num_usm_h_kji, &
1733	.TRUE., .FALSE. )
1734	!
1735	!-- Default surface type
1736	ELSE
1737	CALL initialize_horizontal_surfaces( k, j, i, &
1738	surf_def_h(0), &
1739	num_def_h(0), &
1740	num_def_h_kji(0),&
1741	.TRUE., .FALSE. )
1742	ENDIF
1743	ENDIF
1744	!
1745	!-- downward-facing surface, first, model top. Please note,
1746	!-- for the moment, downward-facing surfaces are always of
1747	!-- default type
1748	IF ( k == nzt .AND. use_top_fluxes ) THEN
1749	CALL initialize_top( k, j, i, surf_def_h(2), &
1750	num_def_h(2), num_def_h_kji(2) )
1751	!
1752	!-- Check for any other downward-facing surface. So far only for
1753	!-- default surface type.
1754	ELSEIF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
1755	CALL initialize_horizontal_surfaces( k, j, i, &
1756	surf_def_h(1), &
1757	num_def_h(1), &
1758	num_def_h_kji(1), &
1759	.FALSE., .TRUE. )
1760	ENDIF
1761	!
1762	!-- Check for vertical walls and, if required, initialize it.
1763	! Start with northward-facing surface.
1764	IF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 0 ) ) THEN
1765	!
1766	!-- Determine flags indicating terrain or building
1767	terrain = BTEST( wall_flags_0(k,j-1,i), 5 ) .OR. &
1768	topo_no_distinct
1769	building = BTEST( wall_flags_0(k,j-1,i), 6 ) .OR. &
1770	topo_no_distinct
1771	IF ( urban_surface .AND. building ) THEN
1772	CALL initialize_vertical_surfaces( 0, k, j, i, &
1773	surf_usm_v(0), &
1774	num_usm_v(0), &
1775	num_usm_v_kji(0), &
1776	.FALSE., .FALSE., &
1777	.FALSE., .TRUE. )
1778	ELSEIF ( land_surface .AND. terrain ) THEN
1779	CALL initialize_vertical_surfaces( 0, k, j, i, &
1780	surf_lsm_v(0), &
1781	num_lsm_v(0), &
1782	num_lsm_v_kji(0), &
1783	.FALSE., .FALSE., &
1784	.FALSE., .TRUE. )
1785	ELSE
1786	CALL initialize_vertical_surfaces( 0, k, j, i, &
1787	surf_def_v(0), &
1788	num_def_v(0), &
1789	num_def_v_kji(0), &
1790	.FALSE., .FALSE., &
1791	.FALSE., .TRUE. )
1792	ENDIF
1793	ENDIF
1794	!
1795	!-- southward-facing surface
1796	IF ( .NOT. BTEST( wall_flags_0(k,j+1,i), 0 ) ) THEN
1797	!
1798	!-- Determine flags indicating terrain or building
1799	terrain = BTEST( wall_flags_0(k,j+1,i), 5 ) .OR. &
1800	topo_no_distinct
1801	building = BTEST( wall_flags_0(k,j+1,i), 6 ) .OR. &
1802	topo_no_distinct
1803	IF ( urban_surface .AND. building ) THEN
1804	CALL initialize_vertical_surfaces( 1, k, j, i, &
1805	surf_usm_v(1), &
1806	num_usm_v(1), &
1807	num_usm_v_kji(1), &
1808	.FALSE., .FALSE., &
1809	.TRUE., .FALSE. )
1810	ELSEIF ( land_surface .AND. terrain ) THEN
1811	CALL initialize_vertical_surfaces( 1, k, j, i, &
1812	surf_lsm_v(1), &
1813	num_lsm_v(1), &
1814	num_lsm_v_kji(1), &
1815	.FALSE., .FALSE., &
1816	.TRUE., .FALSE. )
1817	ELSE
1818	CALL initialize_vertical_surfaces( 1, k, j, i, &
1819	surf_def_v(1), &
1820	num_def_v(1), &
1821	num_def_v_kji(1), &
1822	.FALSE., .FALSE., &
1823	.TRUE., .FALSE. )
1824	ENDIF
1825	ENDIF
1826	!
1827	!-- eastward-facing surface
1828	IF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 0 ) ) THEN
1829	!
1830	!-- Determine flags indicating terrain or building
1831	terrain = BTEST( wall_flags_0(k,j,i-1), 5 ) .OR. &
1832	topo_no_distinct
1833	building = BTEST( wall_flags_0(k,j,i-1), 6 ) .OR. &
1834	topo_no_distinct
1835	IF ( urban_surface .AND. building ) THEN
1836	CALL initialize_vertical_surfaces( 2, k, j, i, &
1837	surf_usm_v(2), &
1838	num_usm_v(2), &
1839	num_usm_v_kji(2), &
1840	.TRUE., .FALSE., &
1841	.FALSE., .FALSE. )
1842	ELSEIF ( land_surface .AND. terrain ) THEN
1843	CALL initialize_vertical_surfaces( 2, k, j, i, &
1844	surf_lsm_v(2), &
1845	num_lsm_v(2), &
1846	num_lsm_v_kji(2), &
1847	.TRUE., .FALSE., &
1848	.FALSE., .FALSE. )
1849	ELSE
1850	CALL initialize_vertical_surfaces( 2, k, j, i, &
1851	surf_def_v(2), &
1852	num_def_v(2), &
1853	num_def_v_kji(2), &
1854	.TRUE., .FALSE., &
1855	.FALSE., .FALSE. )
1856	ENDIF
1857	ENDIF
1858	!
1859	!-- westward-facing surface
1860	IF ( .NOT. BTEST( wall_flags_0(k,j,i+1), 0 ) ) THEN
1861	!
1862	!-- Determine flags indicating terrain or building
1863	terrain = BTEST( wall_flags_0(k,j,i+1), 5 ) .OR. &
1864	topo_no_distinct
1865	building = BTEST( wall_flags_0(k,j,i+1), 6 ) .OR. &
1866	topo_no_distinct
1867	IF ( urban_surface .AND. building ) THEN
1868	CALL initialize_vertical_surfaces( 3, k, j, i, &
1869	surf_usm_v(3), &
1870	num_usm_v(3), &
1871	num_usm_v_kji(3), &
1872	.FALSE., .TRUE., &
1873	.FALSE., .FALSE. )
1874	ELSEIF ( land_surface .AND. terrain ) THEN
1875	CALL initialize_vertical_surfaces( 3, k, j, i, &
1876	surf_lsm_v(3), &
1877	num_lsm_v(3), &
1878	num_lsm_v_kji(3), &
1879	.FALSE., .TRUE., &
1880	.FALSE., .FALSE. )
1881	ELSE
1882	CALL initialize_vertical_surfaces( 3, k, j, i, &
1883	surf_def_v(3), &
1884	num_def_v(3), &
1885	num_def_v_kji(3), &
1886	.FALSE., .TRUE., &
1887	.FALSE., .FALSE. )
1888	ENDIF
1889	ENDIF
1890	ENDIF
1891
1892
1893	ENDDO
1894	!
1895	!-- Determine start- and end-index at grid point (j,i). Also, for
1896	!-- horizontal surfaces more than 1 horizontal surface element can
1897	!-- exist at grid point (j,i) if overhanging structures are present.
1898	!-- Upward-facing surfaces
1899	surf_def_h(0)%start_index(j,i) = start_index_def_h(0)
1900	surf_def_h(0)%end_index(j,i) = surf_def_h(0)%start_index(j,i) + &
1901	num_def_h_kji(0) - 1
1902	start_index_def_h(0) = surf_def_h(0)%end_index(j,i) + 1
1903	!
1904	!-- Downward-facing surfaces, except model top
1905	surf_def_h(1)%start_index(j,i) = start_index_def_h(1)
1906	surf_def_h(1)%end_index(j,i) = surf_def_h(1)%start_index(j,i) + &
1907	num_def_h_kji(1) - 1
1908	start_index_def_h(1) = surf_def_h(1)%end_index(j,i) + 1
1909	!
1910	!-- Downward-facing surfaces -- model top fluxes
1911	surf_def_h(2)%start_index(j,i) = start_index_def_h(2)
1912	surf_def_h(2)%end_index(j,i) = surf_def_h(2)%start_index(j,i) + &
1913	num_def_h_kji(2) - 1
1914	start_index_def_h(2) = surf_def_h(2)%end_index(j,i) + 1
1915	!
1916	!-- Horizontal natural land surfaces
1917	surf_lsm_h%start_index(j,i) = start_index_lsm_h
1918	surf_lsm_h%end_index(j,i) = surf_lsm_h%start_index(j,i) + &
1919	num_lsm_h_kji - 1
1920	start_index_lsm_h = surf_lsm_h%end_index(j,i) + 1
1921	!
1922	!-- Horizontal urban surfaces
1923	surf_usm_h%start_index(j,i) = start_index_usm_h
1924	surf_usm_h%end_index(j,i) = surf_usm_h%start_index(j,i) + &
1925	num_usm_h_kji - 1
1926	start_index_usm_h = surf_usm_h%end_index(j,i) + 1
1927
1928	!
1929	!-- Vertical surfaces - Default type
1930	surf_def_v(0)%start_index(j,i) = start_index_def_v(0)
1931	surf_def_v(1)%start_index(j,i) = start_index_def_v(1)
1932	surf_def_v(2)%start_index(j,i) = start_index_def_v(2)
1933	surf_def_v(3)%start_index(j,i) = start_index_def_v(3)
1934	surf_def_v(0)%end_index(j,i) = start_index_def_v(0) + &
1935	num_def_v_kji(0) - 1
1936	surf_def_v(1)%end_index(j,i) = start_index_def_v(1) + &
1937	num_def_v_kji(1) - 1
1938	surf_def_v(2)%end_index(j,i) = start_index_def_v(2) + &
1939	num_def_v_kji(2) - 1
1940	surf_def_v(3)%end_index(j,i) = start_index_def_v(3) + &
1941	num_def_v_kji(3) - 1
1942	start_index_def_v(0) = surf_def_v(0)%end_index(j,i) + 1
1943	start_index_def_v(1) = surf_def_v(1)%end_index(j,i) + 1
1944	start_index_def_v(2) = surf_def_v(2)%end_index(j,i) + 1
1945	start_index_def_v(3) = surf_def_v(3)%end_index(j,i) + 1
1946	!
1947	!-- Natural type
1948	surf_lsm_v(0)%start_index(j,i) = start_index_lsm_v(0)
1949	surf_lsm_v(1)%start_index(j,i) = start_index_lsm_v(1)
1950	surf_lsm_v(2)%start_index(j,i) = start_index_lsm_v(2)
1951	surf_lsm_v(3)%start_index(j,i) = start_index_lsm_v(3)
1952	surf_lsm_v(0)%end_index(j,i) = start_index_lsm_v(0) + &
1953	num_lsm_v_kji(0) - 1
1954	surf_lsm_v(1)%end_index(j,i) = start_index_lsm_v(1) + &
1955	num_lsm_v_kji(1) - 1
1956	surf_lsm_v(2)%end_index(j,i) = start_index_lsm_v(2) + &
1957	num_lsm_v_kji(2) - 1
1958	surf_lsm_v(3)%end_index(j,i) = start_index_lsm_v(3) + &
1959	num_lsm_v_kji(3) - 1
1960	start_index_lsm_v(0) = surf_lsm_v(0)%end_index(j,i) + 1
1961	start_index_lsm_v(1) = surf_lsm_v(1)%end_index(j,i) + 1
1962	start_index_lsm_v(2) = surf_lsm_v(2)%end_index(j,i) + 1
1963	start_index_lsm_v(3) = surf_lsm_v(3)%end_index(j,i) + 1
1964	!
1965	!-- Urban type
1966	surf_usm_v(0)%start_index(j,i) = start_index_usm_v(0)
1967	surf_usm_v(1)%start_index(j,i) = start_index_usm_v(1)
1968	surf_usm_v(2)%start_index(j,i) = start_index_usm_v(2)
1969	surf_usm_v(3)%start_index(j,i) = start_index_usm_v(3)
1970	surf_usm_v(0)%end_index(j,i) = start_index_usm_v(0) + &
1971	num_usm_v_kji(0) - 1
1972	surf_usm_v(1)%end_index(j,i) = start_index_usm_v(1) + &
1973	num_usm_v_kji(1) - 1
1974	surf_usm_v(2)%end_index(j,i) = start_index_usm_v(2) + &
1975	num_usm_v_kji(2) - 1
1976	surf_usm_v(3)%end_index(j,i) = start_index_usm_v(3) + &
1977	num_usm_v_kji(3) - 1
1978	start_index_usm_v(0) = surf_usm_v(0)%end_index(j,i) + 1
1979	start_index_usm_v(1) = surf_usm_v(1)%end_index(j,i) + 1
1980	start_index_usm_v(2) = surf_usm_v(2)%end_index(j,i) + 1
1981	start_index_usm_v(3) = surf_usm_v(3)%end_index(j,i) + 1
1982
1983
1984	ENDDO
1985	ENDDO
1986
1987	CONTAINS
1988
1989	!------------------------------------------------------------------------------!
1990	! Description:
1991	! ------------
1992	!> Initialize horizontal surface elements, upward- and downward-facing.
1993	!> Note, horizontal surface type alsw comprises model-top fluxes, which are,
1994	!> initialized in a different routine.
1995	!------------------------------------------------------------------------------!
1996	SUBROUTINE initialize_horizontal_surfaces( k, j, i, surf, num_h, &
1997	num_h_kji, upward_facing, &
1998	downward_facing )
1999
2000	IMPLICIT NONE
2001
2002	INTEGER(iwp) :: i !< running index x-direction
2003	INTEGER(iwp) :: j !< running index y-direction
2004	INTEGER(iwp) :: k !< running index z-direction
2005	INTEGER(iwp) :: num_h !< current number of surface element
2006	INTEGER(iwp) :: num_h_kji !< dummy increment
2007	INTEGER(iwp) :: lsp !< running index chemical species
2008	INTEGER(iwp) :: lsp_pr !< running index chemical species??
2009
2010	LOGICAL :: upward_facing !< flag indicating upward-facing surface
2011	LOGICAL :: downward_facing !< flag indicating downward-facing surface
2012
2013	TYPE( surf_type ) :: surf !< respective surface type
2014	!
2015	!-- Store indices of respective surface element
2016	surf%i(num_h) = i
2017	surf%j(num_h) = j
2018	surf%k(num_h) = k
2019	!
2020	!-- Surface orientation, bit 0 is set to 1 for upward-facing surfaces,
2021	!-- bit 1 is for downward-facing surfaces.
2022	IF ( upward_facing ) surf%facing(num_h) = IBSET( surf%facing(num_h), 0 )
2023	IF ( downward_facing ) surf%facing(num_h) = IBSET( surf%facing(num_h), 1 )
2024	!
2025	!-- Initialize surface-layer height
2026	IF ( upward_facing ) THEN
2027	surf%z_mo(num_h) = zu(k) - zw(k-1)
2028	ELSE
2029	surf%z_mo(num_h) = zw(k) - zu(k)
2030	ENDIF
2031
2032	surf%z0(num_h) = roughness_length
2033	surf%z0h(num_h) = z0h_factor * roughness_length
2034	surf%z0q(num_h) = z0h_factor * roughness_length
2035	!
2036	!-- Initialization in case of 1D pre-cursor run
2037	IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )&
2038	THEN
2039	IF ( .NOT. constant_diffusion ) THEN
2040	IF ( constant_flux_layer ) THEN
2041	surf%ol(num_h) = surf%z_mo(num_h) / &
2042	( rif1d(nzb+1) + 1.0E-20_wp )
2043	surf%us(num_h) = us1d
2044	surf%usws(num_h) = usws1d
2045	surf%vsws(num_h) = vsws1d
2046	ELSE
2047	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
2048	surf%us(num_h) = 0.0_wp
2049	surf%usws(num_h) = 0.0_wp
2050	surf%vsws(num_h) = 0.0_wp
2051	ENDIF
2052	ELSE
2053	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
2054	surf%us(num_h) = 0.0_wp
2055	surf%usws(num_h) = 0.0_wp
2056	surf%vsws(num_h) = 0.0_wp
2057	ENDIF
2058	!
2059	!-- Initialization in all other cases
2060	ELSE
2061
2062	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
2063	!
2064	!-- Very small number is required for calculation of Obukhov length
2065	!-- at first timestep
2066	surf%us(num_h) = 1E-30_wp
2067	surf%usws(num_h) = 0.0_wp
2068	surf%vsws(num_h) = 0.0_wp
2069
2070	ENDIF
2071
2072	surf%rib(num_h) = 0.0_wp
2073	surf%uvw_abs(num_h) = 0.0_wp
2074
2075	IF ( .NOT. constant_diffusion ) THEN
2076	surf%u_0(num_h) = 0.0_wp
2077	surf%v_0(num_h) = 0.0_wp
2078	ENDIF
2079
2080	surf%ts(num_h) = 0.0_wp
2081
2082	IF ( humidity ) THEN
2083	surf%qs(num_h) = 0.0_wp
2084	IF ( cloud_physics .AND. microphysics_morrison) THEN
2085	surf%qcs(num_h) = 0.0_wp
2086	surf%ncs(num_h) = 0.0_wp
2087
2088	surf%qcsws(num_h) = 0.0_wp
2089	surf%ncsws(num_h) = 0.0_wp
2090
2091	ENDIF
2092	IF ( cloud_physics .AND. microphysics_seifert) THEN
2093	surf%qrs(num_h) = 0.0_wp
2094	surf%nrs(num_h) = 0.0_wp
2095
2096	surf%qrsws(num_h) = 0.0_wp
2097	surf%nrsws(num_h) = 0.0_wp
2098
2099	surf%pt1(num_h) = 0.0_wp
2100	surf%qv1(num_h) = 0.0_wp
2101
2102	ENDIF
2103	ENDIF
2104
2105	IF ( passive_scalar ) surf%ss(num_h) = 0.0_wp
2106
2107	DO lsp = 1, nvar
2108	IF ( air_chemistry ) surf%css(lsp,num_h) = 0.0_wp
2109	ENDDO
2110	!
2111	!-- Set initial value for surface temperature
2112	surf%pt_surface(num_h) = pt_surface
2113	!
2114	!-- Inititalize surface fluxes of sensible and latent heat, as well as
2115	!-- passive scalar
2116	IF ( use_surface_fluxes ) THEN
2117
2118	IF ( upward_facing ) THEN
2119	IF ( constant_heatflux ) THEN
2120	!
2121	!-- Initialize surface heatflux. However, skip this for now if
2122	!-- if random_heatflux is set. This case, shf is initialized later.
2123	IF ( .NOT. random_heatflux ) THEN
2124	surf%shf(num_h) = surface_heatflux * &
2125	heatflux_input_conversion(k-1)
2126	!
2127	!-- Check if surface heat flux might be replaced by
2128	!-- prescribed wall heatflux
2129	IF ( k-1 /= 0 ) THEN
2130	surf%shf(num_h) = wall_heatflux(0) * &
2131	heatflux_input_conversion(k-1)
2132	ENDIF
2133	ENDIF
2134	ELSE
2135	surf%shf(num_h) = 0.0_wp
2136	ENDIF
2137	!
2138	!-- Set heat-flux at downward-facing surfaces
2139	ELSE
2140	surf%shf(num_h) = wall_heatflux(5) * &
2141	heatflux_input_conversion(k)
2142	ENDIF
2143
2144	IF ( humidity ) THEN
2145	IF ( upward_facing ) THEN
2146	IF ( constant_waterflux ) THEN
2147	surf%qsws(num_h) = surface_waterflux * &
2148	waterflux_input_conversion(k-1)
2149	IF ( k-1 /= 0 ) THEN
2150	surf%qsws(num_h) = wall_humidityflux(0) * &
2151	waterflux_input_conversion(k-1)
2152	ENDIF
2153	ELSE
2154	surf%qsws(num_h) = 0.0_wp
2155	ENDIF
2156	ELSE
2157	surf%qsws(num_h) = wall_humidityflux(5) * &
2158	heatflux_input_conversion(k)
2159	ENDIF
2160	ENDIF
2161
2162	IF ( passive_scalar ) THEN
2163	IF ( upward_facing ) THEN
2164	IF ( constant_scalarflux ) THEN
2165	surf%ssws(num_h) = surface_scalarflux * rho_air_zw(k-1)
2166
2167	IF ( k-1 /= 0 ) &
2168	surf%ssws(num_h) = wall_scalarflux(0) * &
2169	rho_air_zw(k-1)
2170
2171	ELSE
2172	surf%ssws(num_h) = 0.0_wp
2173	ENDIF
2174	ELSE
2175	surf%ssws(num_h) = wall_scalarflux(5) * rho_air_zw(k)
2176	ENDIF
2177	ENDIF
2178
2179	IF ( air_chemistry ) THEN
2180	lsp_pr = 1
2181	DO WHILE ( TRIM( surface_csflux_name( lsp_pr ) ) /= 'novalue' ) !<'novalue' is the default
2182	DO lsp = 1, nvar
2183	!
2184	!-- Assign surface flux for each variable species
2185	IF ( TRIM( spc_names(lsp) ) == TRIM( surface_csflux_name(lsp_pr) ) ) THEN
2186	IF ( upward_facing ) THEN
2187	IF ( constant_csflux(lsp_pr) ) THEN
2188	surf%cssws(lsp,num_h) = &
2189	surface_csflux(lsp_pr) *&
2190	rho_air_zw(k-1)
2191
2192	IF ( k-1 /= 0 ) &
2193	surf%cssws(lsp,num_h) = &
2194	wall_csflux(lsp,0) * &
2195	rho_air_zw(k-1)
2196	ELSE
2197	surf%cssws(lsp,num_h) = 0.0_wp
2198	ENDIF
2199	ELSE
2200	surf%cssws(lsp,num_h) = wall_csflux(lsp,5) * &
2201	rho_air_zw(k)
2202	ENDIF
2203	ENDIF
2204	ENDDO
2205	lsp_pr = lsp_pr + 1
2206	ENDDO
2207	ENDIF
2208
2209	IF ( ocean ) THEN
2210	IF ( upward_facing ) THEN
2211	surf%sasws(num_h) = bottom_salinityflux * rho_air_zw(k-1)
2212	ELSE
2213	surf%sasws(num_h) = 0.0_wp
2214	ENDIF
2215	ENDIF
2216	ENDIF
2217	!
2218	!-- Increment surface indices
2219	num_h = num_h + 1
2220	num_h_kji = num_h_kji + 1
2221
2222
2223	END SUBROUTINE initialize_horizontal_surfaces
2224
2225
2226	!------------------------------------------------------------------------------!
2227	! Description:
2228	! ------------
2229	!> Initialize model-top fluxes. Currently, only the heatflux and salinity flux
2230	!> can be prescribed, latent flux is zero in this case!
2231	!------------------------------------------------------------------------------!
2232	SUBROUTINE initialize_top( k, j, i, surf, num_h, num_h_kji )
2233
2234	IMPLICIT NONE
2235
2236	INTEGER(iwp) :: i !< running index x-direction
2237	INTEGER(iwp) :: j !< running index y-direction
2238	INTEGER(iwp) :: k !< running index z-direction
2239	INTEGER(iwp) :: num_h !< current number of surface element
2240	INTEGER(iwp) :: num_h_kji !< dummy increment
2241	INTEGER(iwp) :: lsp !< running index for chemical species
2242
2243	TYPE( surf_type ) :: surf !< respective surface type
2244	!
2245	!-- Store indices of respective surface element
2246	surf%i(num_h) = i
2247	surf%j(num_h) = j
2248	surf%k(num_h) = k
2249	!
2250	!-- Initialize top heat flux
2251	IF ( constant_top_heatflux ) &
2252	surf%shf(num_h) = top_heatflux * heatflux_input_conversion(nzt+1)
2253	!
2254	!-- Initialization in case of a coupled model run
2255	IF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
2256	surf%shf(num_h) = 0.0_wp
2257	surf%qsws(num_h) = 0.0_wp
2258	ENDIF
2259	!
2260	!-- Prescribe latent heat flux at the top
2261	IF ( humidity ) THEN
2262	surf%qsws(num_h) = 0.0_wp
2263	IF ( cloud_physics .AND. microphysics_morrison ) THEN
2264	surf%ncsws(num_h) = 0.0_wp
2265	surf%qcsws(num_h) = 0.0_wp
2266	ENDIF
2267	IF ( cloud_physics .AND. microphysics_seifert ) THEN
2268	surf%nrsws(num_h) = 0.0_wp
2269	surf%qrsws(num_h) = 0.0_wp
2270	ENDIF
2271	ENDIF
2272	!
2273	!-- Prescribe top scalar flux
2274	IF ( passive_scalar .AND. constant_top_scalarflux ) &
2275	surf%ssws(num_h) = top_scalarflux * rho_air_zw(nzt+1)
2276	!
2277	!-- Prescribe top chemical species' flux
2278	DO lsp = 1, nvar
2279	IF ( air_chemistry .AND. constant_top_csflux(lsp) ) THEN
2280	surf%cssws(lsp,num_h) = top_csflux(lsp) * rho_air_zw(nzt+1)
2281	ENDIF
2282	ENDDO
2283	!
2284	!-- Prescribe top salinity flux
2285	IF ( ocean .AND. constant_top_salinityflux) &
2286	surf%sasws(num_h) = top_salinityflux * rho_air_zw(nzt+1)
2287	!
2288	!-- Top momentum fluxes
2289	IF ( constant_top_momentumflux ) THEN
2290	surf%usws(num_h) = top_momentumflux_u * &
2291	momentumflux_input_conversion(nzt+1)
2292	surf%vsws(num_h) = top_momentumflux_v * &
2293	momentumflux_input_conversion(nzt+1)
2294	ENDIF
2295	!
2296	!-- Increment surface indices
2297	num_h = num_h + 1
2298	num_h_kji = num_h_kji + 1
2299
2300
2301	END SUBROUTINE initialize_top
2302
2303
2304	!------------------------------------------------------------------------------!
2305	! Description:
2306	! ------------
2307	!> Initialize vertical surface elements.
2308	!------------------------------------------------------------------------------!
2309	SUBROUTINE initialize_vertical_surfaces( l, k, j, i, surf, num_v, &
2310	num_v_kji, east_facing, &
2311	west_facing, south_facing, &
2312	north_facing )
2313
2314	IMPLICIT NONE
2315
2316	INTEGER(iwp) :: component !< index of wall_fluxes_ array for respective orientation
2317	INTEGER(iwp) :: i !< running index x-direction
2318	INTEGER(iwp) :: j !< running index x-direction
2319	INTEGER(iwp) :: k !< running index x-direction
2320	INTEGER(iwp) :: l !< index variable for the surface type, indicating the facing
2321	INTEGER(iwp) :: num_v !< current number of surface element
2322	INTEGER(iwp) :: num_v_kji !< current number of surface element at (j,i)
2323	INTEGER(iwp) :: lsp !< running index for chemical species
2324
2325
2326	LOGICAL :: east_facing !< flag indicating east-facing surfaces
2327	LOGICAL :: north_facing !< flag indicating north-facing surfaces
2328	LOGICAL :: south_facing !< flag indicating south-facing surfaces
2329	LOGICAL :: west_facing !< flag indicating west-facing surfaces
2330
2331	TYPE( surf_type ) :: surf !< respective surface type
2332
2333	!
2334	!-- Store indices of respective wall element
2335	surf%i(num_v) = i
2336	surf%j(num_v) = j
2337	surf%k(num_v) = k
2338	!
2339	!-- Initialize surface-layer height, or more precisely, distance to surface
2340	IF ( north_facing .OR. south_facing ) THEN
2341	surf%z_mo(num_v) = 0.5_wp * dy
2342	ELSE
2343	surf%z_mo(num_v) = 0.5_wp * dx
2344	ENDIF
2345
2346	surf%facing(num_v) = 0
2347	!
2348	!-- Surface orientation. Moreover, set component id to map wall_heatflux,
2349	!-- etc., on surface type (further below)
2350	IF ( north_facing ) THEN
2351	surf%facing(num_v) = 5 !IBSET( surf%facing(num_v), 0 )
2352	component = 4
2353	ENDIF
2354
2355	IF ( south_facing ) THEN
2356	surf%facing(num_v) = 6 !IBSET( surf%facing(num_v), 1 )
2357	component = 3
2358	ENDIF
2359
2360	IF ( east_facing ) THEN
2361	surf%facing(num_v) = 7 !IBSET( surf%facing(num_v), 2 )
2362	component = 2
2363	ENDIF
2364
2365	IF ( west_facing ) THEN
2366	surf%facing(num_v) = 8 !IBSET( surf%facing(num_v), 3 )
2367	component = 1
2368	ENDIF
2369
2370
2371	surf%z0(num_v) = roughness_length
2372	surf%z0h(num_v) = z0h_factor * roughness_length
2373	surf%z0q(num_v) = z0h_factor * roughness_length
2374
2375	surf%us(num_v) = 0.0_wp
2376	!
2377	!-- If required, initialize Obukhov length
2378	IF ( ALLOCATED( surf%ol ) ) &
2379	surf%ol(num_v) = surf%z_mo(num_v) / zeta_min
2380
2381	surf%uvw_abs(num_v) = 0.0_wp
2382
2383	surf%mom_flux_uv(num_v) = 0.0_wp
2384	surf%mom_flux_w(num_v) = 0.0_wp
2385	surf%mom_flux_tke(0:1,num_v) = 0.0_wp
2386
2387	surf%ts(num_v) = 0.0_wp
2388	surf%shf(num_v) = wall_heatflux(component)
2389	!
2390	!-- Set initial value for surface temperature
2391	surf%pt_surface(num_v) = pt_surface
2392
2393	IF ( humidity ) THEN
2394	surf%qs(num_v) = 0.0_wp
2395	surf%qsws(num_v) = wall_humidityflux(component)
2396	!
2397	!-- Following wall fluxes are assumed to be zero
2398	IF ( cloud_physics .AND. microphysics_morrison) THEN
2399	surf%qcs(num_v) = 0.0_wp
2400	surf%ncs(num_v) = 0.0_wp
2401
2402	surf%qcsws(num_v) = 0.0_wp
2403	surf%ncsws(num_v) = 0.0_wp
2404	ENDIF
2405	IF ( cloud_physics .AND. microphysics_seifert) THEN
2406	surf%qrs(num_v) = 0.0_wp
2407	surf%nrs(num_v) = 0.0_wp
2408
2409	surf%qrsws(num_v) = 0.0_wp
2410	surf%nrsws(num_v) = 0.0_wp
2411	ENDIF
2412	ENDIF
2413
2414	IF ( passive_scalar ) THEN
2415	surf%ss(num_v) = 0.0_wp
2416	surf%ssws(num_v) = wall_scalarflux(component)
2417	ENDIF
2418
2419	IF ( air_chemistry ) THEN
2420	DO lsp = 1, nvar
2421	surf%css(lsp,num_v) = 0.0_wp
2422	surf%cssws(lsp,num_v) = wall_csflux(lsp,component)
2423	ENDDO
2424	ENDIF
2425
2426	!
2427	!-- So far, salinityflux at vertical surfaces is simply zero
2428	!-- at the moment
2429	IF ( ocean ) surf%sasws(num_v) = wall_salinityflux(component)
2430	!
2431	!-- Increment wall indices
2432	num_v = num_v + 1
2433	num_v_kji = num_v_kji + 1
2434
2435	END SUBROUTINE initialize_vertical_surfaces
2436
2437	END SUBROUTINE init_surfaces
2438
2439
2440	!------------------------------------------------------------------------------!
2441	! Description:
2442	! ------------
2443	!> Determines topography-top index at given (j,i)-position.
2444	!------------------------------------------------------------------------------!
2445	FUNCTION get_topography_top_index_ji( j, i, grid )
2446
2447	IMPLICIT NONE
2448
2449	CHARACTER(LEN=*) :: grid !< flag to distinquish between staggered grids
2450	INTEGER(iwp) :: i !< grid index in x-dimension
2451	INTEGER(iwp) :: ibit !< bit position where topography information is stored on respective grid
2452	INTEGER(iwp) :: j !< grid index in y-dimension
2453	INTEGER(iwp) :: get_topography_top_index_ji !< topography top index
2454
2455	SELECT CASE ( TRIM( grid ) )
2456
2457	CASE ( 's' )
2458	ibit = 12
2459	CASE ( 'u' )
2460	ibit = 14
2461	CASE ( 'v' )
2462	ibit = 16
2463	CASE ( 'w' )
2464	ibit = 18
2465	CASE ( 's_out' )
2466	ibit = 24
2467	CASE DEFAULT
2468	!
2469	!-- Set default to scalar grid
2470	ibit = 12
2471
2472	END SELECT
2473
2474	get_topography_top_index_ji = MAXLOC( &
2475	MERGE( 1, 0, &
2476	BTEST( wall_flags_0(:,j,i), ibit ) &
2477	), DIM = 1 &
2478	) - 1
2479
2480	RETURN
2481
2482	END FUNCTION get_topography_top_index_ji
2483
2484	!------------------------------------------------------------------------------!
2485	! Description:
2486	! ------------
2487	!> Determines topography-top index at each (j,i)-position.
2488	!------------------------------------------------------------------------------!
2489	FUNCTION get_topography_top_index( grid )
2490
2491	IMPLICIT NONE
2492
2493	CHARACTER(LEN=*) :: grid !< flag to distinquish between staggered grids
2494	INTEGER(iwp) :: ibit !< bit position where topography information is stored on respective grid
2495	INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) :: get_topography_top_index !< topography top index
2496
2497	SELECT CASE ( TRIM( grid ) )
2498
2499	CASE ( 's' )
2500	ibit = 12
2501	CASE ( 'u' )
2502	ibit = 14
2503	CASE ( 'v' )
2504	ibit = 16
2505	CASE ( 'w' )
2506	ibit = 18
2507	CASE ( 's_out' )
2508	ibit = 24
2509	CASE DEFAULT
2510	!
2511	!-- Set default to scalar grid
2512	ibit = 12
2513
2514	END SELECT
2515
2516	get_topography_top_index(nys:nyn,nxl:nxr) = MAXLOC( &
2517	MERGE( 1, 0, &
2518	BTEST( wall_flags_0(:,nys:nyn,nxl:nxr), ibit )&
2519	), DIM = 1 &
2520	) - 1
2521
2522	RETURN
2523
2524	END FUNCTION get_topography_top_index
2525
2526	!------------------------------------------------------------------------------!
2527	! Description:
2528	! ------------
2529	!> Gathers all surface elements with the same facing (but possibly different
2530	!> type) onto a surface type, and writes binary data into restart files.
2531	!------------------------------------------------------------------------------!
2532	SUBROUTINE surface_wrd_local
2533
2534
2535	IMPLICIT NONE
2536
2537	CHARACTER(LEN=1) :: dum !< dummy string to create output-variable name
2538
2539	INTEGER(iwp) :: i !< running index x-direction
2540	INTEGER(iwp) :: j !< running index y-direction
2541	INTEGER(iwp) :: l !< index surface type orientation
2542	INTEGER(iwp) :: lsp !< running index chemical species
2543	INTEGER(iwp) :: m !< running index for surface elements on individual surface array
2544	INTEGER(iwp), DIMENSION(0:2) :: start_index_h !< start index for horizontal surface elements on gathered surface array
2545	INTEGER(iwp), DIMENSION(0:3) :: mm !< running index for surface elements on gathered surface array
2546	INTEGER(iwp), DIMENSION(0:3) :: start_index_v !< start index for vertical surface elements on gathered surface array
2547
2548	TYPE(surf_type), DIMENSION(0:2) :: surf_h !< gathered horizontal surfaces, contains all surface types
2549	TYPE(surf_type), DIMENSION(0:3) :: surf_v !< gathered vertical surfaces, contains all surface types
2550
2551	!
2552	!-- Determine total number of horizontal and vertical surface elements before
2553	!-- writing var_list
2554	CALL surface_last_actions
2555	!
2556	!-- Count number of grid points with same facing and allocate attributes respectively
2557	!-- Horizontal upward facing
2558	surf_h(0)%ns = ns_h_on_file(0)
2559	CALL allocate_surface_attributes_h( surf_h(0), nys, nyn, nxl, nxr )
2560	!
2561	!-- Horizontal downward facing
2562	surf_h(1)%ns = ns_h_on_file(1)
2563	CALL allocate_surface_attributes_h( surf_h(1), nys, nyn, nxl, nxr )
2564	!
2565	!-- Model top
2566	surf_h(2)%ns = ns_h_on_file(2)
2567	CALL allocate_surface_attributes_h_top( surf_h(2), nys, nyn, nxl, nxr )
2568	!
2569	!-- Vertical surfaces
2570	DO l = 0, 3
2571	surf_v(l)%ns = ns_v_on_file(l)
2572	CALL allocate_surface_attributes_v( surf_v(l), &
2573	nys, nyn, nxl, nxr )
2574	ENDDO
2575	!
2576	!-- In the following, gather data from surfaces elements with the same
2577	!-- facing (but possibly differt type) on 1 data-type array.
2578	mm(0:2) = 1
2579	DO l = 0, 2
2580	DO i = nxl, nxr
2581	DO j = nys, nyn
2582	DO m = surf_def_h(l)%start_index(j,i), &
2583	surf_def_h(l)%end_index(j,i)
2584	IF ( ALLOCATED( surf_def_h(l)%us ) ) &
2585	surf_h(l)%us(mm(l)) = surf_def_h(l)%us(m)
2586	IF ( ALLOCATED( surf_def_h(l)%ts ) ) &
2587	surf_h(l)%ts(mm(l)) = surf_def_h(l)%ts(m)
2588	IF ( ALLOCATED( surf_def_h(l)%qs ) ) &
2589	surf_h(l)%qs(mm(l)) = surf_def_h(l)%qs(m)
2590	IF ( ALLOCATED( surf_def_h(l)%ss ) ) &
2591	surf_h(l)%ss(mm(l)) = surf_def_h(l)%ss(m)
2592	IF ( ALLOCATED( surf_def_h(l)%qcs ) ) &
2593	surf_h(l)%qcs(mm(l)) = surf_def_h(l)%qcs(m)
2594	IF ( ALLOCATED( surf_def_h(l)%ncs ) ) &
2595	surf_h(l)%ncs(mm(l)) = surf_def_h(l)%ncs(m)
2596	IF ( ALLOCATED( surf_def_h(l)%qrs ) ) &
2597	surf_h(l)%qrs(mm(l)) = surf_def_h(l)%qrs(m)
2598	IF ( ALLOCATED( surf_def_h(l)%nrs ) ) &
2599	surf_h(l)%nrs(mm(l)) = surf_def_h(l)%nrs(m)
2600	IF ( ALLOCATED( surf_def_h(l)%ol ) ) &
2601	surf_h(l)%ol(mm(l)) = surf_def_h(l)%ol(m)
2602	IF ( ALLOCATED( surf_def_h(l)%rib ) ) &
2603	surf_h(l)%rib(mm(l)) = surf_def_h(l)%rib(m)
2604	IF ( ALLOCATED( surf_def_h(l)%pt_surface ) ) &
2605	surf_h(l)%pt_surface(mm(l)) = surf_def_h(l)%pt_surface(m)
2606	IF ( ALLOCATED( surf_def_h(l)%usws ) ) &
2607	surf_h(l)%usws(mm(l)) = surf_def_h(l)%usws(m)
2608	IF ( ALLOCATED( surf_def_h(l)%vsws ) ) &
2609	surf_h(l)%vsws(mm(l)) = surf_def_h(l)%vsws(m)
2610	IF ( ALLOCATED( surf_def_h(l)%shf ) ) &
2611	surf_h(l)%shf(mm(l)) = surf_def_h(l)%shf(m)
2612	IF ( ALLOCATED( surf_def_h(l)%qsws ) ) &
2613	surf_h(l)%qsws(mm(l)) = surf_def_h(l)%qsws(m)
2614	IF ( ALLOCATED( surf_def_h(l)%ssws ) ) &
2615	surf_h(l)%ssws(mm(l)) = surf_def_h(l)%ssws(m)
2616	IF ( ALLOCATED( surf_def_h(l)%css ) ) THEN
2617	DO lsp = 1,nvar
2618	surf_h(l)%css(lsp,mm(l)) = surf_def_h(l)%css(lsp,m)
2619	ENDDO
2620	ENDIF
2621	IF ( ALLOCATED( surf_def_h(l)%cssws ) ) THEN
2622	DO lsp = 1,nvar
2623	surf_h(l)%cssws(lsp,mm(l)) = surf_def_h(l)%cssws(lsp,m)
2624	ENDDO
2625	ENDIF
2626	IF ( ALLOCATED( surf_def_h(l)%ncsws ) ) &
2627	surf_h(l)%ncsws(mm(l)) = surf_def_h(l)%ncsws(m)
2628	IF ( ALLOCATED( surf_def_h(l)%nrsws ) ) &
2629	surf_h(l)%nrsws(mm(l)) = surf_def_h(l)%nrsws(m)
2630	IF ( ALLOCATED( surf_def_h(l)%sasws ) ) &
2631	surf_h(l)%sasws(mm(l)) = surf_def_h(l)%sasws(m)
2632
2633	mm(l) = mm(l) + 1
2634	ENDDO
2635
2636	IF ( l == 0 ) THEN
2637	DO m = surf_lsm_h%start_index(j,i), &
2638	surf_lsm_h%end_index(j,i)
2639	IF ( ALLOCATED( surf_lsm_h%us ) ) &
2640	surf_h(0)%us(mm(0)) = surf_lsm_h%us(m)
2641	IF ( ALLOCATED( surf_lsm_h%ts ) ) &
2642	surf_h(0)%ts(mm(0)) = surf_lsm_h%ts(m)
2643	IF ( ALLOCATED( surf_lsm_h%qs ) ) &
2644	surf_h(0)%qs(mm(0)) = surf_lsm_h%qs(m)
2645	IF ( ALLOCATED( surf_lsm_h%ss ) ) &
2646	surf_h(0)%ss(mm(0)) = surf_lsm_h%ss(m)
2647	IF ( ALLOCATED( surf_lsm_h%qcs ) ) &
2648	surf_h(0)%qcs(mm(0)) = surf_lsm_h%qcs(m)
2649	IF ( ALLOCATED( surf_lsm_h%ncs ) ) &
2650	surf_h(0)%ncs(mm(0)) = surf_lsm_h%ncs(m)
2651	IF ( ALLOCATED( surf_lsm_h%qrs ) ) &
2652	surf_h(0)%qrs(mm(0)) = surf_lsm_h%qrs(m)
2653	IF ( ALLOCATED( surf_lsm_h%nrs ) ) &
2654	surf_h(0)%nrs(mm(0)) = surf_lsm_h%nrs(m)
2655	IF ( ALLOCATED( surf_lsm_h%ol ) ) &
2656	surf_h(0)%ol(mm(0)) = surf_lsm_h%ol(m)
2657	IF ( ALLOCATED( surf_lsm_h%rib ) ) &
2658	surf_h(0)%rib(mm(0)) = surf_lsm_h%rib(m)
2659	IF ( ALLOCATED( surf_lsm_h%pt_surface ) ) &
2660	surf_h(l)%pt_surface(mm(l)) = surf_lsm_h%pt_surface(m)
2661	IF ( ALLOCATED( surf_lsm_h%usws ) ) &
2662	surf_h(0)%usws(mm(0)) = surf_lsm_h%usws(m)
2663	IF ( ALLOCATED( surf_lsm_h%vsws ) ) &
2664	surf_h(0)%vsws(mm(0)) = surf_lsm_h%vsws(m)
2665	IF ( ALLOCATED( surf_lsm_h%shf ) ) &
2666	surf_h(0)%shf(mm(0)) = surf_lsm_h%shf(m)
2667	IF ( ALLOCATED( surf_lsm_h%qsws ) ) &
2668	surf_h(0)%qsws(mm(0)) = surf_lsm_h%qsws(m)
2669	IF ( ALLOCATED( surf_lsm_h%ssws ) ) &
2670	surf_h(0)%ssws(mm(0)) = surf_lsm_h%ssws(m)
2671	IF ( ALLOCATED( surf_lsm_h%css ) ) THEN
2672	DO lsp = 1, nvar
2673	surf_h(0)%css(lsp,mm(0)) = surf_lsm_h%css(lsp,m)
2674	ENDDO
2675	ENDIF
2676	IF ( ALLOCATED( surf_lsm_h%cssws ) ) THEN
2677	DO lsp = 1, nvar
2678	surf_h(0)%cssws(lsp,mm(0)) = surf_lsm_h%cssws(lsp,m)
2679	ENDDO
2680	ENDIF
2681	IF ( ALLOCATED( surf_lsm_h%ncsws ) ) &
2682	surf_h(0)%ncsws(mm(0)) = surf_lsm_h%ncsws(m)
2683	IF ( ALLOCATED( surf_lsm_h%nrsws ) ) &
2684	surf_h(0)%nrsws(mm(0)) = surf_lsm_h%nrsws(m)
2685	IF ( ALLOCATED( surf_lsm_h%sasws ) ) &
2686	surf_h(0)%sasws(mm(0)) = surf_lsm_h%sasws(m)
2687
2688	mm(0) = mm(0) + 1
2689
2690	ENDDO
2691
2692	DO m = surf_usm_h%start_index(j,i), &
2693	surf_usm_h%end_index(j,i)
2694	IF ( ALLOCATED( surf_usm_h%us ) ) &
2695	surf_h(0)%us(mm(0)) = surf_usm_h%us(m)
2696	IF ( ALLOCATED( surf_usm_h%ts ) ) &
2697	surf_h(0)%ts(mm(0)) = surf_usm_h%ts(m)
2698	IF ( ALLOCATED( surf_usm_h%qs ) ) &
2699	surf_h(0)%qs(mm(0)) = surf_usm_h%qs(m)
2700	IF ( ALLOCATED( surf_usm_h%ss ) ) &
2701	surf_h(0)%ss(mm(0)) = surf_usm_h%ss(m)
2702	IF ( ALLOCATED( surf_usm_h%qcs ) ) &
2703	surf_h(0)%qcs(mm(0)) = surf_usm_h%qcs(m)
2704	IF ( ALLOCATED( surf_usm_h%ncs ) ) &
2705	surf_h(0)%ncs(mm(0)) = surf_usm_h%ncs(m)
2706	IF ( ALLOCATED( surf_usm_h%qrs ) ) &
2707	surf_h(0)%qrs(mm(0)) = surf_usm_h%qrs(m)
2708	IF ( ALLOCATED( surf_usm_h%nrs ) ) &
2709	surf_h(0)%nrs(mm(0)) = surf_usm_h%nrs(m)
2710	IF ( ALLOCATED( surf_usm_h%ol ) ) &
2711	surf_h(0)%ol(mm(0)) = surf_usm_h%ol(m)
2712	IF ( ALLOCATED( surf_usm_h%rib ) ) &
2713	surf_h(0)%rib(mm(0)) = surf_usm_h%rib(m)
2714	IF ( ALLOCATED( surf_usm_h%pt_surface ) ) &
2715	surf_h(l)%pt_surface(mm(l)) = surf_usm_h%pt_surface(m)
2716	IF ( ALLOCATED( surf_usm_h%usws ) ) &
2717	surf_h(0)%usws(mm(0)) = surf_usm_h%usws(m)
2718	IF ( ALLOCATED( surf_usm_h%vsws ) ) &
2719	surf_h(0)%vsws(mm(0)) = surf_usm_h%vsws(m)
2720	IF ( ALLOCATED( surf_usm_h%shf ) ) &
2721	surf_h(0)%shf(mm(0)) = surf_usm_h%shf(m)
2722	IF ( ALLOCATED( surf_usm_h%qsws ) ) &
2723	surf_h(0)%qsws(mm(0)) = surf_usm_h%qsws(m)
2724	IF ( ALLOCATED( surf_usm_h%ssws ) ) &
2725	surf_h(0)%ssws(mm(0)) = surf_usm_h%ssws(m)
2726	IF ( ALLOCATED( surf_usm_h%css ) ) THEN
2727	DO lsp = 1, nvar
2728	surf_h(0)%css(lsp,mm(0)) = surf_usm_h%css(lsp,m)
2729	ENDDO
2730	ENDIF
2731	IF ( ALLOCATED( surf_usm_h%cssws ) ) THEN
2732	DO lsp = 1, nvar
2733	surf_h(0)%cssws(lsp,mm(0)) = surf_usm_h%cssws(lsp,m)
2734	ENDDO
2735	ENDIF
2736	IF ( ALLOCATED( surf_usm_h%ncsws ) ) &
2737	surf_h(0)%ncsws(mm(0)) = surf_usm_h%ncsws(m)
2738	IF ( ALLOCATED( surf_usm_h%nrsws ) ) &
2739	surf_h(0)%nrsws(mm(0)) = surf_usm_h%nrsws(m)
2740	IF ( ALLOCATED( surf_usm_h%sasws ) ) &
2741	surf_h(0)%sasws(mm(0)) = surf_usm_h%sasws(m)
2742
2743	mm(0) = mm(0) + 1
2744
2745	ENDDO
2746
2747
2748	ENDIF
2749
2750	ENDDO
2751
2752	ENDDO
2753	!
2754	!-- Gather start- and end indices
2755	start_index_h(l) = 1
2756	DO i = nxl, nxr
2757	DO j = nys, nyn
2758
2759	surf_h(l)%start_index(j,i) = start_index_h(l)
2760	surf_h(l)%end_index(j,i) = surf_h(l)%start_index(j,i) -1
2761
2762	DO m = surf_def_h(l)%start_index(j,i), &
2763	surf_def_h(l)%end_index(j,i)
2764	surf_h(l)%end_index(j,i) = surf_h(l)%end_index(j,i) + 1
2765	ENDDO
2766	IF ( l == 0 ) THEN
2767	DO m = surf_lsm_h%start_index(j,i), &
2768	surf_lsm_h%end_index(j,i)
2769	surf_h(l)%end_index(j,i) = surf_h(l)%end_index(j,i) + 1
2770	ENDDO
2771	DO m = surf_usm_h%start_index(j,i), &
2772	surf_usm_h%end_index(j,i)
2773	surf_h(l)%end_index(j,i) = surf_h(l)%end_index(j,i) + 1
2774	ENDDO
2775	ENDIF
2776
2777	start_index_h(l) = surf_h(l)%end_index(j,i) + 1
2778
2779	ENDDO
2780	ENDDO
2781	ENDDO
2782
2783
2784	mm(0:3) = 1
2785	DO l = 0, 3
2786	DO i = nxl, nxr
2787	DO j = nys, nyn
2788	DO m = surf_def_v(l)%start_index(j,i), &
2789	surf_def_v(l)%end_index(j,i)
2790	IF ( ALLOCATED( surf_def_v(l)%us ) ) &
2791	surf_v(l)%us(mm(l)) = surf_def_v(l)%us(m)
2792	IF ( ALLOCATED( surf_def_v(l)%ts ) ) &
2793	surf_v(l)%ts(mm(l)) = surf_def_v(l)%ts(m)
2794	IF ( ALLOCATED( surf_def_v(l)%qs ) ) &
2795	surf_v(l)%qs(mm(l)) = surf_def_v(l)%qs(m)
2796	IF ( ALLOCATED( surf_def_v(l)%ss ) ) &
2797	surf_v(l)%ss(mm(l)) = surf_def_v(l)%ss(m)
2798	IF ( ALLOCATED( surf_def_v(l)%qcs ) ) &
2799	surf_v(l)%qcs(mm(l)) = surf_def_v(l)%qcs(m)
2800	IF ( ALLOCATED( surf_def_v(l)%ncs ) ) &
2801	surf_v(l)%ncs(mm(l)) = surf_def_v(l)%ncs(m)
2802	IF ( ALLOCATED( surf_def_v(l)%qrs ) ) &
2803	surf_v(l)%qrs(mm(l)) = surf_def_v(l)%qrs(m)
2804	IF ( ALLOCATED( surf_def_v(l)%nrs ) ) &
2805	surf_v(l)%nrs(mm(l)) = surf_def_v(l)%nrs(m)
2806	IF ( ALLOCATED( surf_def_v(l)%ol ) ) &
2807	surf_v(l)%ol(mm(l)) = surf_def_v(l)%ol(m)
2808	IF ( ALLOCATED( surf_def_v(l)%rib ) ) &
2809	surf_v(l)%rib(mm(l)) = surf_def_v(l)%rib(m)
2810	IF ( ALLOCATED( surf_def_v(l)%pt_surface ) ) &
2811	surf_v(l)%pt_surface(mm(l)) = surf_def_v(l)%pt_surface(m)
2812	IF ( ALLOCATED( surf_def_v(l)%shf ) ) &
2813	surf_v(l)%shf(mm(l)) = surf_def_v(l)%shf(m)
2814	IF ( ALLOCATED( surf_def_v(l)%qsws ) ) &
2815	surf_v(l)%qsws(mm(l)) = surf_def_v(l)%qsws(m)
2816	IF ( ALLOCATED( surf_def_v(l)%ssws ) ) &
2817	surf_v(l)%ssws(mm(l)) = surf_def_v(l)%ssws(m)
2818	IF ( ALLOCATED( surf_def_v(l)%css ) ) THEN
2819	DO lsp = 1, nvar
2820	surf_v(l)%css(lsp,mm(l)) = surf_def_v(l)%css(lsp,m)
2821	ENDDO
2822	ENDIF
2823	IF ( ALLOCATED( surf_def_v(l)%cssws ) ) THEN
2824	DO lsp = 1, nvar
2825	surf_v(l)%cssws(lsp,mm(l)) = surf_def_v(l)%cssws(lsp,m)
2826	ENDDO
2827	ENDIF
2828	IF ( ALLOCATED( surf_def_v(l)%ncsws ) ) &
2829	surf_v(l)%ncsws(mm(l)) = surf_def_v(l)%ncsws(m)
2830	IF ( ALLOCATED( surf_def_v(l)%nrsws ) ) &
2831	surf_v(l)%nrsws(mm(l)) = surf_def_v(l)%nrsws(m)
2832	IF ( ALLOCATED( surf_def_v(l)%sasws ) ) &
2833	surf_v(l)%sasws(mm(l)) = surf_def_v(l)%sasws(m)
2834	IF ( ALLOCATED( surf_def_v(l)%mom_flux_uv) ) &
2835	surf_v(l)%mom_flux_uv(mm(l)) = surf_def_v(l)%mom_flux_uv(m)
2836	IF ( ALLOCATED( surf_def_v(l)%mom_flux_w) ) &
2837	surf_v(l)%mom_flux_w(mm(l)) = surf_def_v(l)%mom_flux_w(m)
2838	IF ( ALLOCATED( surf_def_v(l)%mom_flux_tke) ) &
2839	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_def_v(l)%mom_flux_tke(0:1,m)
2840
2841	mm(l) = mm(l) + 1
2842	ENDDO
2843
2844	DO m = surf_lsm_v(l)%start_index(j,i), &
2845	surf_lsm_v(l)%end_index(j,i)
2846	IF ( ALLOCATED( surf_lsm_v(l)%us ) ) &
2847	surf_v(l)%us(mm(l)) = surf_lsm_v(l)%us(m)
2848	IF ( ALLOCATED( surf_lsm_v(l)%ts ) ) &
2849	surf_v(l)%ts(mm(l)) = surf_lsm_v(l)%ts(m)
2850	IF ( ALLOCATED( surf_lsm_v(l)%qs ) ) &
2851	surf_v(l)%qs(mm(l)) = surf_lsm_v(l)%qs(m)
2852	IF ( ALLOCATED( surf_lsm_v(l)%ss ) ) &
2853	surf_v(l)%ss(mm(l)) = surf_lsm_v(l)%ss(m)
2854	IF ( ALLOCATED( surf_lsm_v(l)%qcs ) ) &
2855	surf_v(l)%qcs(mm(l)) = surf_lsm_v(l)%qcs(m)
2856	IF ( ALLOCATED( surf_lsm_v(l)%ncs ) ) &
2857	surf_v(l)%ncs(mm(l)) = surf_lsm_v(l)%ncs(m)
2858	IF ( ALLOCATED( surf_lsm_v(l)%qrs ) ) &
2859	surf_v(l)%qrs(mm(l)) = surf_lsm_v(l)%qrs(m)
2860	IF ( ALLOCATED( surf_lsm_v(l)%nrs ) ) &
2861	surf_v(l)%nrs(mm(l)) = surf_lsm_v(l)%nrs(m)
2862	IF ( ALLOCATED( surf_lsm_v(l)%ol ) ) &
2863	surf_v(l)%ol(mm(l)) = surf_lsm_v(l)%ol(m)
2864	IF ( ALLOCATED( surf_lsm_v(l)%rib ) ) &
2865	surf_v(l)%rib(mm(l)) = surf_lsm_v(l)%rib(m)
2866	IF ( ALLOCATED( surf_lsm_v(l)%pt_surface ) ) &
2867	surf_v(l)%pt_surface(mm(l)) = surf_lsm_v(l)%pt_surface(m)
2868	IF ( ALLOCATED( surf_lsm_v(l)%usws ) ) &
2869	surf_v(l)%usws(mm(l)) = surf_lsm_v(l)%usws(m)
2870	IF ( ALLOCATED( surf_lsm_v(l)%vsws ) ) &
2871	surf_v(l)%vsws(mm(l)) = surf_lsm_v(l)%vsws(m)
2872	IF ( ALLOCATED( surf_lsm_v(l)%shf ) ) &
2873	surf_v(l)%shf(mm(l)) = surf_lsm_v(l)%shf(m)
2874	IF ( ALLOCATED( surf_lsm_v(l)%qsws ) ) &
2875	surf_v(l)%qsws(mm(l)) = surf_lsm_v(l)%qsws(m)
2876	IF ( ALLOCATED( surf_lsm_v(l)%ssws ) ) &
2877	surf_v(l)%ssws(mm(l)) = surf_lsm_v(l)%ssws(m)
2878	IF ( ALLOCATED( surf_lsm_v(l)%css ) ) THEN
2879	DO lsp = 1, nvar
2880	surf_v(l)%css(lsp,mm(l)) = surf_lsm_v(l)%css(lsp,m)
2881	ENDDO
2882	ENDIF
2883	IF ( ALLOCATED( surf_lsm_v(l)%cssws ) ) THEN
2884	DO lsp = 1, nvar
2885	surf_v(l)%cssws(lsp,mm(l)) = surf_lsm_v(l)%cssws(lsp,m)
2886	ENDDO
2887	ENDIF
2888	IF ( ALLOCATED( surf_lsm_v(l)%ncsws ) ) &
2889	surf_v(l)%ncsws(mm(l)) = surf_lsm_v(l)%ncsws(m)
2890	IF ( ALLOCATED( surf_lsm_v(l)%nrsws ) ) &
2891	surf_v(l)%nrsws(mm(l)) = surf_lsm_v(l)%nrsws(m)
2892	IF ( ALLOCATED( surf_lsm_v(l)%sasws ) ) &
2893	surf_v(l)%sasws(mm(l)) = surf_lsm_v(l)%sasws(m)
2894	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_uv) ) &
2895	surf_v(l)%mom_flux_uv(mm(l)) = surf_lsm_v(l)%mom_flux_uv(m)
2896	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_w) ) &
2897	surf_v(l)%mom_flux_w(mm(l)) = surf_lsm_v(l)%mom_flux_w(m)
2898	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_tke) ) &
2899	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_lsm_v(l)%mom_flux_tke(0:1,m)
2900
2901	mm(l) = mm(l) + 1
2902	ENDDO
2903
2904	DO m = surf_usm_v(l)%start_index(j,i), &
2905	surf_usm_v(l)%end_index(j,i)
2906	IF ( ALLOCATED( surf_usm_v(l)%us ) ) &
2907	surf_v(l)%us(mm(l)) = surf_usm_v(l)%us(m)
2908	IF ( ALLOCATED( surf_usm_v(l)%ts ) ) &
2909	surf_v(l)%ts(mm(l)) = surf_usm_v(l)%ts(m)
2910	IF ( ALLOCATED( surf_usm_v(l)%qs ) ) &
2911	surf_v(l)%qs(mm(l)) = surf_usm_v(l)%qs(m)
2912	IF ( ALLOCATED( surf_usm_v(l)%ss ) ) &
2913	surf_v(l)%ss(mm(l)) = surf_usm_v(l)%ss(m)
2914	IF ( ALLOCATED( surf_usm_v(l)%qcs ) ) &
2915	surf_v(l)%qcs(mm(l)) = surf_usm_v(l)%qcs(m)
2916	IF ( ALLOCATED( surf_usm_v(l)%ncs ) ) &
2917	surf_v(l)%ncs(mm(l)) = surf_usm_v(l)%ncs(m)
2918	IF ( ALLOCATED( surf_usm_v(l)%qrs ) ) &
2919	surf_v(l)%qrs(mm(l)) = surf_usm_v(l)%qrs(m)
2920	IF ( ALLOCATED( surf_usm_v(l)%nrs ) ) &
2921	surf_v(l)%nrs(mm(l)) = surf_usm_v(l)%nrs(m)
2922	IF ( ALLOCATED( surf_usm_v(l)%ol ) ) &
2923	surf_v(l)%ol(mm(l)) = surf_usm_v(l)%ol(m)
2924	IF ( ALLOCATED( surf_usm_v(l)%rib ) ) &
2925	surf_v(l)%rib(mm(l)) = surf_usm_v(l)%rib(m)
2926	IF ( ALLOCATED( surf_usm_v(l)%pt_surface ) ) &
2927	surf_v(l)%pt_surface(mm(l)) = surf_usm_v(l)%pt_surface(m)
2928	IF ( ALLOCATED( surf_usm_v(l)%usws ) ) &
2929	surf_v(l)%usws(mm(l)) = surf_usm_v(l)%usws(m)
2930	IF ( ALLOCATED( surf_usm_v(l)%vsws ) ) &
2931	surf_v(l)%vsws(mm(l)) = surf_usm_v(l)%vsws(m)
2932	IF ( ALLOCATED( surf_usm_v(l)%shf ) ) &
2933	surf_v(l)%shf(mm(l)) = surf_usm_v(l)%shf(m)
2934	IF ( ALLOCATED( surf_usm_v(l)%qsws ) ) &
2935	surf_v(l)%qsws(mm(l)) = surf_usm_v(l)%qsws(m)
2936	IF ( ALLOCATED( surf_usm_v(l)%ssws ) ) &
2937	surf_v(l)%ssws(mm(l)) = surf_usm_v(l)%ssws(m)
2938	IF ( ALLOCATED( surf_usm_v(l)%css ) ) THEN
2939	DO lsp = 1, nvar
2940	surf_v(l)%css(lsp,mm(l)) = surf_usm_v(l)%css(lsp,m)
2941	ENDDO
2942	ENDIF
2943	IF ( ALLOCATED( surf_usm_v(l)%cssws ) ) THEN
2944	DO lsp = 1, nvar
2945	surf_v(l)%cssws(lsp,mm(l)) = surf_usm_v(l)%cssws(lsp,m)
2946	ENDDO
2947	ENDIF
2948	IF ( ALLOCATED( surf_usm_v(l)%ncsws ) ) &
2949	surf_v(l)%ncsws(mm(l)) = surf_usm_v(l)%ncsws(m)
2950	IF ( ALLOCATED( surf_usm_v(l)%nrsws ) ) &
2951	surf_v(l)%nrsws(mm(l)) = surf_usm_v(l)%nrsws(m)
2952	IF ( ALLOCATED( surf_usm_v(l)%sasws ) ) &
2953	surf_v(l)%sasws(mm(l)) = surf_usm_v(l)%sasws(m)
2954	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_uv) ) &
2955	surf_v(l)%mom_flux_uv(mm(l)) = surf_usm_v(l)%mom_flux_uv(m)
2956	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_w) ) &
2957	surf_v(l)%mom_flux_w(mm(l)) = surf_usm_v(l)%mom_flux_w(m)
2958	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_tke) ) &
2959	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_usm_v(l)%mom_flux_tke(0:1,m)
2960
2961	mm(l) = mm(l) + 1
2962	ENDDO
2963
2964	ENDDO
2965	ENDDO
2966	!
2967	!-- Gather start- and end indices
2968	start_index_v(l) = 1
2969	DO i = nxl, nxr
2970	DO j = nys, nyn
2971
2972	surf_v(l)%start_index(j,i) = start_index_v(l)
2973	surf_v(l)%end_index(j,i) = surf_v(l)%start_index(j,i) -1
2974
2975	DO m = surf_def_v(l)%start_index(j,i), &
2976	surf_def_v(l)%end_index(j,i)
2977	surf_v(l)%end_index(j,i) = surf_v(l)%end_index(j,i) + 1
2978	ENDDO
2979	DO m = surf_lsm_v(l)%start_index(j,i), &
2980	surf_lsm_v(l)%end_index(j,i)
2981	surf_v(l)%end_index(j,i) = surf_v(l)%end_index(j,i) + 1
2982	ENDDO
2983	DO m = surf_usm_v(l)%start_index(j,i), &
2984	surf_usm_v(l)%end_index(j,i)
2985	surf_v(l)%end_index(j,i) = surf_v(l)%end_index(j,i) + 1
2986	ENDDO
2987
2988	start_index_v(l) = surf_v(l)%end_index(j,i) + 1
2989	ENDDO
2990	ENDDO
2991
2992	ENDDO
2993
2994
2995	CALL wrd_write_string( 'ns_h_on_file' )
2996	WRITE ( 14 ) ns_h_on_file
2997
2998	CALL wrd_write_string( 'ns_v_on_file' )
2999	WRITE ( 14 ) ns_v_on_file
3000
3001	!
3002	!-- Write required restart data.
3003	!-- Start with horizontal surfaces (upward-, downward-facing, and model top)
3004	DO l = 0, 2
3005	WRITE( dum, '(I1)') l
3006
3007	CALL wrd_write_string( 'surf_h(' // dum // ')%start_index' )
3008	WRITE ( 14 ) surf_h(l)%start_index
3009
3010	CALL wrd_write_string( 'surf_h(' // dum // ')%end_index' )
3011	WRITE ( 14 ) surf_h(l)%end_index
3012
3013	IF ( ALLOCATED ( surf_h(l)%us ) ) THEN
3014	CALL wrd_write_string( 'surf_h(' // dum // ')%us' )
3015	WRITE ( 14 ) surf_h(l)%us
3016	ENDIF
3017
3018	IF ( ALLOCATED ( surf_h(l)%ts ) ) THEN
3019	CALL wrd_write_string( 'surf_h(' // dum // ')%ts' )
3020	WRITE ( 14 ) surf_h(l)%ts
3021	ENDIF
3022
3023	IF ( ALLOCATED ( surf_h(l)%qs ) ) THEN
3024	CALL wrd_write_string( 'surf_h(' // dum // ')%qs' )
3025	WRITE ( 14 ) surf_h(l)%qs
3026	ENDIF
3027
3028	IF ( ALLOCATED ( surf_h(l)%ss ) ) THEN
3029	CALL wrd_write_string( 'surf_h(' // dum // ')%ss' )
3030	WRITE ( 14 ) surf_h(l)%ss
3031	ENDIF
3032
3033	IF ( ALLOCATED ( surf_h(l)%qcs ) ) THEN
3034	CALL wrd_write_string( 'surf_h(' // dum // ')%qcs' )
3035	WRITE ( 14 ) surf_h(l)%qcs
3036	ENDIF
3037
3038	IF ( ALLOCATED ( surf_h(l)%ncs ) ) THEN
3039	CALL wrd_write_string( 'surf_h(' // dum // ')%ncs' )
3040	WRITE ( 14 ) surf_h(l)%ncs
3041	ENDIF
3042
3043	IF ( ALLOCATED ( surf_h(l)%qrs ) ) THEN
3044	CALL wrd_write_string( 'surf_h(' // dum // ')%qrs' )
3045	WRITE ( 14 ) surf_h(l)%qrs
3046	ENDIF
3047
3048	IF ( ALLOCATED ( surf_h(l)%nrs ) ) THEN
3049	CALL wrd_write_string( 'surf_h(' // dum // ')%nrs' )
3050	WRITE ( 14 ) surf_h(l)%nrs
3051	ENDIF
3052
3053	IF ( ALLOCATED ( surf_h(l)%ol ) ) THEN
3054	CALL wrd_write_string( 'surf_h(' // dum // ')%ol' )
3055	WRITE ( 14 ) surf_h(l)%ol
3056	ENDIF
3057
3058	IF ( ALLOCATED ( surf_h(l)%rib ) ) THEN
3059	CALL wrd_write_string( 'surf_h(' // dum // ')%rib' )
3060	WRITE ( 14 ) surf_h(l)%rib
3061	ENDIF
3062
3063	IF ( ALLOCATED ( surf_h(l)%pt_surface ) ) THEN
3064	CALL wrd_write_string( 'surf_h(' // dum // ')%pt_surface' )
3065	WRITE ( 14 ) surf_h(l)%pt_surface
3066	ENDIF
3067
3068	IF ( ALLOCATED ( surf_h(l)%usws ) ) THEN
3069	CALL wrd_write_string( 'surf_h(' // dum // ')%usws' )
3070	WRITE ( 14 ) surf_h(l)%usws
3071	ENDIF
3072
3073	IF ( ALLOCATED ( surf_h(l)%vsws ) ) THEN
3074	CALL wrd_write_string( 'surf_h(' // dum // ')%vsws' )
3075	WRITE ( 14 ) surf_h(l)%vsws
3076	ENDIF
3077
3078	IF ( ALLOCATED ( surf_h(l)%shf ) ) THEN
3079	CALL wrd_write_string( 'surf_h(' // dum // ')%shf' )
3080	WRITE ( 14 ) surf_h(l)%shf
3081	ENDIF
3082
3083	IF ( ALLOCATED ( surf_h(l)%qsws ) ) THEN
3084	CALL wrd_write_string( 'surf_h(' // dum // ')%qsws' )
3085	WRITE ( 14 ) surf_h(l)%qsws
3086	ENDIF
3087
3088	IF ( ALLOCATED ( surf_h(l)%ssws ) ) THEN
3089	CALL wrd_write_string( 'surf_h(' // dum // ')%ssws' )
3090	WRITE ( 14 ) surf_h(l)%ssws
3091	ENDIF
3092
3093	IF ( ALLOCATED ( surf_h(l)%css ) ) THEN
3094	CALL wrd_write_string( 'surf_h(' // dum // ')%css' )
3095	WRITE ( 14 ) surf_h(l)%css
3096	ENDIF
3097
3098	IF ( ALLOCATED ( surf_h(l)%cssws ) ) THEN
3099	CALL wrd_write_string( 'surf_h(' // dum // ')%cssws' )
3100	WRITE ( 14 ) surf_h(l)%cssws
3101	ENDIF
3102
3103	IF ( ALLOCATED ( surf_h(l)%qcsws ) ) THEN
3104	CALL wrd_write_string( 'surf_h(' // dum // ')%qcsws' )
3105	WRITE ( 14 ) surf_h(l)%qcsws
3106	ENDIF
3107
3108	IF ( ALLOCATED ( surf_h(l)%ncsws ) ) THEN
3109	CALL wrd_write_string( 'surf_h(' // dum // ')%ncsws' )
3110	WRITE ( 14 ) surf_h(l)%ncsws
3111	ENDIF
3112
3113	IF ( ALLOCATED ( surf_h(l)%qrsws ) ) THEN
3114	CALL wrd_write_string( 'surf_h(' // dum // ')%qrsws' )
3115	WRITE ( 14 ) surf_h(l)%qrsws
3116	ENDIF
3117
3118	IF ( ALLOCATED ( surf_h(l)%nrsws ) ) THEN
3119	CALL wrd_write_string( 'surf_h(' // dum // ')%nrsws' )
3120	WRITE ( 14 ) surf_h(l)%nrsws
3121	ENDIF
3122
3123	IF ( ALLOCATED ( surf_h(l)%sasws ) ) THEN
3124	CALL wrd_write_string( 'surf_h(' // dum // ')%sasws' )
3125	WRITE ( 14 ) surf_h(l)%sasws
3126	ENDIF
3127
3128	ENDDO
3129	!
3130	!-- Write vertical surfaces
3131	DO l = 0, 3
3132	WRITE( dum, '(I1)') l
3133
3134	CALL wrd_write_string( 'surf_v(' // dum // ')%start_index' )
3135	WRITE ( 14 ) surf_v(l)%start_index
3136
3137	CALL wrd_write_string( 'surf_v(' // dum // ')%end_index' )
3138	WRITE ( 14 ) surf_v(l)%end_index
3139
3140	IF ( ALLOCATED ( surf_v(l)%us ) ) THEN
3141	CALL wrd_write_string( 'surf_v(' // dum // ')%us' )
3142	WRITE ( 14 ) surf_v(l)%us
3143	ENDIF
3144
3145	IF ( ALLOCATED ( surf_v(l)%ts ) ) THEN
3146	CALL wrd_write_string( 'surf_v(' // dum // ')%ts' )
3147	WRITE ( 14 ) surf_v(l)%ts
3148	ENDIF
3149
3150	IF ( ALLOCATED ( surf_v(l)%qs ) ) THEN
3151	CALL wrd_write_string( 'surf_v(' // dum // ')%qs' )
3152	WRITE ( 14 ) surf_v(l)%qs
3153	ENDIF
3154
3155	IF ( ALLOCATED ( surf_v(l)%ss ) ) THEN
3156	CALL wrd_write_string( 'surf_v(' // dum // ')%ss' )
3157	WRITE ( 14 ) surf_v(l)%ss
3158	ENDIF
3159
3160	IF ( ALLOCATED ( surf_v(l)%qcs ) ) THEN
3161	CALL wrd_write_string( 'surf_v(' // dum // ')%qcs' )
3162	WRITE ( 14 ) surf_v(l)%qcs
3163	ENDIF
3164
3165	IF ( ALLOCATED ( surf_v(l)%ncs ) ) THEN
3166	CALL wrd_write_string( 'surf_v(' // dum // ')%ncs' )
3167	WRITE ( 14 ) surf_v(l)%ncs
3168	ENDIF
3169
3170	IF ( ALLOCATED ( surf_v(l)%qrs ) ) THEN
3171	CALL wrd_write_string( 'surf_v(' // dum // ')%qrs' )
3172	WRITE ( 14 ) surf_v(l)%qrs
3173	ENDIF
3174
3175	IF ( ALLOCATED ( surf_v(l)%nrs ) ) THEN
3176	CALL wrd_write_string( 'surf_v(' // dum // ')%nrs' )
3177	WRITE ( 14 ) surf_v(l)%nrs
3178	ENDIF
3179
3180	IF ( ALLOCATED ( surf_v(l)%ol ) ) THEN
3181	CALL wrd_write_string( 'surf_v(' // dum // ')%ol' )
3182	WRITE ( 14 ) surf_v(l)%ol
3183	ENDIF
3184
3185	IF ( ALLOCATED ( surf_v(l)%rib ) ) THEN
3186	CALL wrd_write_string( 'surf_v(' // dum // ')%rib' )
3187	WRITE ( 14 ) surf_v(l)%rib
3188	ENDIF
3189
3190	IF ( ALLOCATED ( surf_v(l)%pt_surface ) ) THEN
3191	CALL wrd_write_string( 'surf_v(' // dum // ')%pt_surface' )
3192	WRITE ( 14 ) surf_v(l)%pt_surface
3193	ENDIF
3194
3195	IF ( ALLOCATED ( surf_v(l)%shf ) ) THEN
3196	CALL wrd_write_string( 'surf_v(' // dum // ')%shf' )
3197	WRITE ( 14 ) surf_v(l)%shf
3198	ENDIF
3199
3200	IF ( ALLOCATED ( surf_v(l)%qsws ) ) THEN
3201	CALL wrd_write_string( 'surf_v(' // dum // ')%qsws' )
3202	WRITE ( 14 ) surf_v(l)%qsws
3203	ENDIF
3204
3205	IF ( ALLOCATED ( surf_v(l)%ssws ) ) THEN
3206	CALL wrd_write_string( 'surf_v(' // dum // ')%ssws' )
3207	WRITE ( 14 ) surf_v(l)%ssws
3208	ENDIF
3209
3210	IF ( ALLOCATED ( surf_v(l)%css ) ) THEN
3211	CALL wrd_write_string( 'surf_v(' // dum // ')%css' )
3212	WRITE ( 14 ) surf_v(l)%css
3213	ENDIF
3214
3215	IF ( ALLOCATED ( surf_v(l)%cssws ) ) THEN
3216	CALL wrd_write_string( 'surf_v(' // dum // ')%cssws' )
3217	WRITE ( 14 ) surf_v(l)%cssws
3218	ENDIF
3219
3220	IF ( ALLOCATED ( surf_v(l)%qcsws ) ) THEN
3221	CALL wrd_write_string( 'surf_v(' // dum // ')%qcsws' )
3222	WRITE ( 14 ) surf_v(l)%qcsws
3223	ENDIF
3224
3225	IF ( ALLOCATED ( surf_v(l)%ncsws ) ) THEN
3226	CALL wrd_write_string( 'surf_v(' // dum // ')%ncsws' )
3227	WRITE ( 14 ) surf_v(l)%ncsws
3228	ENDIF
3229
3230	IF ( ALLOCATED ( surf_v(l)%qrsws ) ) THEN
3231	CALL wrd_write_string( 'surf_v(' // dum // ')%qrsws' )
3232	WRITE ( 14 ) surf_v(l)%qrsws
3233	ENDIF
3234
3235	IF ( ALLOCATED ( surf_v(l)%nrsws ) ) THEN
3236	CALL wrd_write_string( 'surf_v(' // dum // ')%nrsws' )
3237	WRITE ( 14 ) surf_v(l)%nrsws
3238	ENDIF
3239
3240	IF ( ALLOCATED ( surf_v(l)%sasws ) ) THEN
3241	CALL wrd_write_string( 'surf_v(' // dum // ')%sasws' )
3242	WRITE ( 14 ) surf_v(l)%sasws
3243	ENDIF
3244
3245	IF ( ALLOCATED ( surf_v(l)%mom_flux_uv ) ) THEN
3246	CALL wrd_write_string( 'surf_v(' // dum // ')%mom_uv' )
3247	WRITE ( 14 ) surf_v(l)%mom_flux_uv
3248	ENDIF
3249
3250	IF ( ALLOCATED ( surf_v(l)%mom_flux_w ) ) THEN
3251	CALL wrd_write_string( 'surf_v(' // dum // ')%mom_w' )
3252	WRITE ( 14 ) surf_v(l)%mom_flux_w
3253	ENDIF
3254
3255	IF ( ALLOCATED ( surf_v(l)%mom_flux_tke ) ) THEN
3256	CALL wrd_write_string( 'surf_v(' // dum // ')%mom_tke' )
3257	WRITE ( 14 ) surf_v(l)%mom_flux_tke
3258	ENDIF
3259
3260	ENDDO
3261
3262
3263	END SUBROUTINE surface_wrd_local
3264
3265
3266	!------------------------------------------------------------------------------!
3267	! Description:
3268	! ------------
3269	!> Reads surface-related restart data. Please note, restart data for a certain
3270	!> surface orientation (e.g. horizontal upward-facing) is stored in one
3271	!> array, even if surface elements may belong to different surface types
3272	!> natural or urban for example). Surface elements are redistributed into its
3273	!> respective surface types within this routine. This allows e.g. changing the
3274	!> surface type after reading the restart data, which might be required in case
3275	!> of cyclic_fill mode.
3276	!------------------------------------------------------------------------------!
3277	SUBROUTINE surface_rrd_local( ii, kk, nxlf, nxlc, nxl_on_file, nxrf, nxrc, &
3278	nxr_on_file, nynf, nync, nyn_on_file, nysf, &
3279	nysc, nys_on_file, found )
3280
3281	USE pegrid, &
3282	ONLY: myid, numprocs_previous_run
3283
3284	IMPLICIT NONE
3285
3286	CHARACTER (LEN=1) :: dum !< dummy to create correct string for reading input variable
3287
3288	INTEGER(iwp) :: i !< running index along x-direction, refers to former domain size
3289	INTEGER(iwp) :: ic !< running index along x-direction, refers to current domain size
3290	INTEGER(iwp) :: j !< running index along y-direction, refers to former domain size
3291	INTEGER(iwp) :: jc !< running index along y-direction, refers to former domain size
3292	INTEGER(iwp) :: k !< running index along z-direction
3293	INTEGER(iwp) :: m !< running index for surface elements, refers to gathered array encompassing all surface types
3294	INTEGER(iwp) :: mm !< running index for surface elements, refers to individual surface types
3295	INTEGER(iwp) :: ii !< running index over input files
3296	INTEGER(iwp) :: kk !< running index over previous input files covering current local domain
3297	INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain
3298	INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain
3299	INTEGER(iwp) :: nxl_on_file !< index of left boundary on former local domain
3300	INTEGER(iwp) :: nxrc !< index of right boundary on current subdomain
3301	INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain
3302	INTEGER(iwp) :: nxr_on_file !< index of right boundary on former local domain
3303	INTEGER(iwp) :: nync !< index of north boundary on current subdomain
3304	INTEGER(iwp) :: nynf !< index of north boundary on former subdomain
3305	INTEGER(iwp) :: nyn_on_file !< index of norht boundary on former local domain
3306	INTEGER(iwp) :: nysc !< index of south boundary on current subdomain
3307	INTEGER(iwp) :: nysf !< index of south boundary on former subdomain
3308	INTEGER(iwp) :: nys_on_file !< index of south boundary on former local domain
3309
3310	INTEGER(iwp), SAVE :: l !< index variable for surface type
3311
3312	LOGICAL :: surf_match_def !< flag indicating that surface element is of default type
3313	LOGICAL :: surf_match_lsm !< flag indicating that surface element is of natural type
3314	LOGICAL :: surf_match_usm !< flag indicating that surface element is of urban type
3315
3316	LOGICAL, INTENT(OUT) :: found
3317
3318	LOGICAL, SAVE :: horizontal_surface !< flag indicating horizontal surfaces
3319	LOGICAL, SAVE :: vertical_surface !< flag indicating vertical surfaces
3320
3321	TYPE(surf_type), DIMENSION(0:2), SAVE :: surf_h !< horizontal surface type on file
3322	TYPE(surf_type), DIMENSION(0:3), SAVE :: surf_v !< vertical surface type on file
3323
3324
3325	found = .TRUE.
3326
3327	SELECT CASE ( restart_string(1:length) )
3328
3329	CASE ( 'ns_h_on_file' )
3330	IF ( kk == 1 ) THEN
3331	READ ( 13 ) ns_h_on_file
3332
3333	IF ( ALLOCATED( surf_h(0)%start_index ) ) &
3334	CALL deallocate_surface_attributes_h( surf_h(0) )
3335	IF ( ALLOCATED( surf_h(1)%start_index ) ) &
3336	CALL deallocate_surface_attributes_h( surf_h(1) )
3337	IF ( ALLOCATED( surf_h(2)%start_index ) ) &
3338	CALL deallocate_surface_attributes_h_top( surf_h(2) )
3339
3340	!-- Allocate memory for number of surface elements on file.
3341	!-- Please note, these number is not necessarily the same as
3342	!-- the final number of surface elements on local domain,
3343	!-- which is the case if processor topology changes during
3344	!-- restart runs.
3345	!-- Horizontal upward facing
3346	surf_h(0)%ns = ns_h_on_file(0)
3347	CALL allocate_surface_attributes_h( surf_h(0), &
3348	nys_on_file, nyn_on_file, &
3349	nxl_on_file, nxr_on_file )
3350
3351	!-- Horizontal downward facing
3352	surf_h(1)%ns = ns_h_on_file(1)
3353	CALL allocate_surface_attributes_h( surf_h(1), &
3354	nys_on_file, nyn_on_file, &
3355	nxl_on_file, nxr_on_file )
3356
3357	!-- Model top
3358	surf_h(2)%ns = ns_h_on_file(2)
3359	CALL allocate_surface_attributes_h_top( surf_h(2), &
3360	nys_on_file, nyn_on_file, &
3361	nxl_on_file, nxr_on_file )
3362
3363	!
3364	!-- Initial setting of flags for horizontal and vertical surfaces,
3365	!-- will be set after start- and end-indices are read.
3366	horizontal_surface = .FALSE.
3367	vertical_surface = .FALSE.
3368
3369	ENDIF
3370
3371	CASE ( 'ns_v_on_file' )
3372	IF ( kk == 1 ) THEN
3373	READ ( 13 ) ns_v_on_file
3374
3375	DO l = 0, 3
3376	IF ( ALLOCATED( surf_v(l)%start_index ) ) &
3377	CALL deallocate_surface_attributes_v( surf_v(l) )
3378	ENDDO
3379
3380	!-- Vertical surfaces
3381	DO l = 0, 3
3382	surf_v(l)%ns = ns_v_on_file(l)
3383	CALL allocate_surface_attributes_v( surf_v(l), &
3384	nys_on_file, nyn_on_file, &
3385	nxl_on_file, nxr_on_file )
3386	ENDDO
3387
3388	ENDIF
3389
3390	CASE ( 'surf_h(0)%start_index' )
3391	IF ( kk == 1 ) &
3392	READ ( 13 ) surf_h(0)%start_index
3393	l = 0
3394	CASE ( 'surf_h(0)%end_index' )
3395	IF ( kk == 1 ) &
3396	READ ( 13 ) surf_h(0)%end_index
3397	horizontal_surface = .TRUE.
3398	vertical_surface = .FALSE.
3399	CASE ( 'surf_h(0)%us' )
3400	IF ( ALLOCATED( surf_h(0)%us ) .AND. kk == 1 ) &
3401	READ ( 13 ) surf_h(0)%us
3402	CASE ( 'surf_h(0)%ts' )
3403	IF ( ALLOCATED( surf_h(0)%ts ) .AND. kk == 1 ) &
3404	READ ( 13 ) surf_h(0)%ts
3405	CASE ( 'surf_h(0)%qs' )
3406	IF ( ALLOCATED( surf_h(0)%qs ) .AND. kk == 1 ) &
3407	READ ( 13 ) surf_h(0)%qs
3408	CASE ( 'surf_h(0)%ss' )
3409	IF ( ALLOCATED( surf_h(0)%ss ) .AND. kk == 1 ) &
3410	READ ( 13 ) surf_h(0)%ss
3411	CASE ( 'surf_h(0)%qcs' )
3412	IF ( ALLOCATED( surf_h(0)%qcs ) .AND. kk == 1 ) &
3413	READ ( 13 ) surf_h(0)%qcs
3414	CASE ( 'surf_h(0)%ncs' )
3415	IF ( ALLOCATED( surf_h(0)%ncs ) .AND. kk == 1 ) &
3416	READ ( 13 ) surf_h(0)%ncs
3417	CASE ( 'surf_h(0)%qrs' )
3418	IF ( ALLOCATED( surf_h(0)%qrs ) .AND. kk == 1 ) &
3419	READ ( 13 ) surf_h(0)%qrs
3420	CASE ( 'surf_h(0)%nrs' )
3421	IF ( ALLOCATED( surf_h(0)%nrs ) .AND. kk == 1 ) &
3422	READ ( 13 ) surf_h(0)%nrs
3423	CASE ( 'surf_h(0)%ol' )
3424	IF ( ALLOCATED( surf_h(0)%ol ) .AND. kk == 1 ) &
3425	READ ( 13 ) surf_h(0)%ol
3426	CASE ( 'surf_h(0)%rib' )
3427	IF ( ALLOCATED( surf_h(0)%rib ) .AND. kk == 1 ) &
3428	READ ( 13 ) surf_h(0)%rib
3429	CASE ( 'surf_h(0)%pt_surface' )
3430	IF ( ALLOCATED( surf_h(0)%pt_surface ) .AND. kk == 1 ) &
3431	READ ( 13 ) surf_h(0)%pt_surface
3432	CASE ( 'surf_h(0)%usws' )
3433	IF ( ALLOCATED( surf_h(0)%usws ) .AND. kk == 1 ) &
3434	READ ( 13 ) surf_h(0)%usws
3435	CASE ( 'surf_h(0)%vsws' )
3436	IF ( ALLOCATED( surf_h(0)%vsws ) .AND. kk == 1 ) &
3437	READ ( 13 ) surf_h(0)%vsws
3438	CASE ( 'surf_h(0)%shf' )
3439	IF ( ALLOCATED( surf_h(0)%shf ) .AND. kk == 1 ) &
3440	READ ( 13 ) surf_h(0)%shf
3441	CASE ( 'surf_h(0)%qsws' )
3442	IF ( ALLOCATED( surf_h(0)%qsws ) .AND. kk == 1 ) &
3443	READ ( 13 ) surf_h(0)%qsws
3444	CASE ( 'surf_h(0)%ssws' )
3445	IF ( ALLOCATED( surf_h(0)%ssws ) .AND. kk == 1 ) &
3446	READ ( 13 ) surf_h(0)%ssws
3447	CASE ( 'surf_h(0)%css' )
3448	IF ( ALLOCATED( surf_h(0)%css ) .AND. kk == 1 ) &
3449	READ ( 13 ) surf_h(0)%css
3450	CASE ( 'surf_h(0)%cssws' )
3451	IF ( ALLOCATED( surf_h(0)%cssws ) .AND. kk == 1 ) &
3452	READ ( 13 ) surf_h(0)%cssws
3453	CASE ( 'surf_h(0)%qcsws' )
3454	IF ( ALLOCATED( surf_h(0)%qcsws ) .AND. kk == 1 ) &
3455	READ ( 13 ) surf_h(0)%qcsws
3456	CASE ( 'surf_h(0)%ncsws' )
3457	IF ( ALLOCATED( surf_h(0)%ncsws ) .AND. kk == 1 ) &
3458	READ ( 13 ) surf_h(0)%ncsws
3459	CASE ( 'surf_h(0)%qrsws' )
3460	IF ( ALLOCATED( surf_h(0)%qrsws ) .AND. kk == 1 ) &
3461	READ ( 13 ) surf_h(0)%qrsws
3462	CASE ( 'surf_h(0)%nrsws' )
3463	IF ( ALLOCATED( surf_h(0)%nrsws ) .AND. kk == 1 ) &
3464	READ ( 13 ) surf_h(0)%nrsws
3465	CASE ( 'surf_h(0)%sasws' )
3466	IF ( ALLOCATED( surf_h(0)%sasws ) .AND. kk == 1 ) &
3467	READ ( 13 ) surf_h(0)%sasws
3468
3469	CASE ( 'surf_h(1)%start_index' )
3470	IF ( kk == 1 ) &
3471	READ ( 13 ) surf_h(1)%start_index
3472	l = 1
3473	CASE ( 'surf_h(1)%end_index' )
3474	IF ( kk == 1 ) &
3475	READ ( 13 ) surf_h(1)%end_index
3476	CASE ( 'surf_h(1)%us' )
3477	IF ( ALLOCATED( surf_h(1)%us ) .AND. kk == 1 ) &
3478	READ ( 13 ) surf_h(1)%us
3479	CASE ( 'surf_h(1)%ts' )
3480	IF ( ALLOCATED( surf_h(1)%ts ) .AND. kk == 1 ) &
3481	READ ( 13 ) surf_h(1)%ts
3482	CASE ( 'surf_h(1)%qs' )
3483	IF ( ALLOCATED( surf_h(1)%qs ) .AND. kk == 1 ) &
3484	READ ( 13 ) surf_h(1)%qs
3485	CASE ( 'surf_h(1)%ss' )
3486	IF ( ALLOCATED( surf_h(1)%ss ) .AND. kk == 1 ) &
3487	READ ( 13 ) surf_h(1)%ss
3488	CASE ( 'surf_h(1)%qcs' )
3489	IF ( ALLOCATED( surf_h(1)%qcs ) .AND. kk == 1 ) &
3490	READ ( 13 ) surf_h(1)%qcs
3491	CASE ( 'surf_h(1)%ncs' )
3492	IF ( ALLOCATED( surf_h(1)%ncs ) .AND. kk == 1 ) &
3493	READ ( 13 ) surf_h(1)%ncs
3494	CASE ( 'surf_h(1)%qrs' )
3495	IF ( ALLOCATED( surf_h(1)%qrs ) .AND. kk == 1 ) &
3496	READ ( 13 ) surf_h(1)%qrs
3497	CASE ( 'surf_h(1)%nrs' )
3498	IF ( ALLOCATED( surf_h(1)%nrs ) .AND. kk == 1 ) &
3499	READ ( 13 ) surf_h(1)%nrs
3500	CASE ( 'surf_h(1)%ol' )
3501	IF ( ALLOCATED( surf_h(1)%ol ) .AND. kk == 1 ) &
3502	READ ( 13 ) surf_h(1)%ol
3503	CASE ( 'surf_h(1)%rib' )
3504	IF ( ALLOCATED( surf_h(1)%rib ) .AND. kk == 1 ) &
3505	READ ( 13 ) surf_h(1)%rib
3506	CASE ( 'surf_h(1)%pt_surface' )
3507	IF ( ALLOCATED( surf_h(1)%pt_surface ) .AND. kk == 1 ) &
3508	READ ( 13 ) surf_h(1)%pt_surface
3509	CASE ( 'surf_h(1)%usws' )
3510	IF ( ALLOCATED( surf_h(1)%usws ) .AND. kk == 1 ) &
3511	READ ( 13 ) surf_h(1)%usws
3512	CASE ( 'surf_h(1)%vsws' )
3513	IF ( ALLOCATED( surf_h(1)%vsws ) .AND. kk == 1 ) &
3514	READ ( 13 ) surf_h(1)%vsws
3515	CASE ( 'surf_h(1)%shf' )
3516	IF ( ALLOCATED( surf_h(1)%shf ) .AND. kk == 1 ) &
3517	READ ( 13 ) surf_h(1)%shf
3518	CASE ( 'surf_h(1)%qsws' )
3519	IF ( ALLOCATED( surf_h(1)%qsws ) .AND. kk == 1 ) &
3520	READ ( 13 ) surf_h(1)%qsws
3521	CASE ( 'surf_h(1)%ssws' )
3522	IF ( ALLOCATED( surf_h(1)%ssws ) .AND. kk == 1 ) &
3523	READ ( 13 ) surf_h(1)%ssws
3524	CASE ( 'surf_h(1)%css' )
3525	IF ( ALLOCATED( surf_h(1)%css ) .AND. kk == 1 ) &
3526	READ ( 13 ) surf_h(1)%css
3527	CASE ( 'surf_h(1)%cssws' )
3528	IF ( ALLOCATED( surf_h(1)%cssws ) .AND. kk == 1 ) &
3529	READ ( 13 ) surf_h(1)%cssws
3530	CASE ( 'surf_h(1)%qcsws' )
3531	IF ( ALLOCATED( surf_h(1)%qcsws ) .AND. kk == 1 ) &
3532	READ ( 13 ) surf_h(1)%qcsws
3533	CASE ( 'surf_h(1)%ncsws' )
3534	IF ( ALLOCATED( surf_h(1)%ncsws ) .AND. kk == 1 ) &
3535	READ ( 13 ) surf_h(1)%ncsws
3536	CASE ( 'surf_h(1)%qrsws' )
3537	IF ( ALLOCATED( surf_h(1)%qrsws ) .AND. kk == 1 ) &
3538	READ ( 13 ) surf_h(1)%qrsws
3539	CASE ( 'surf_h(1)%nrsws' )
3540	IF ( ALLOCATED( surf_h(1)%nrsws ) .AND. kk == 1 ) &
3541	READ ( 13 ) surf_h(1)%nrsws
3542	CASE ( 'surf_h(1)%sasws' )
3543	IF ( ALLOCATED( surf_h(1)%sasws ) .AND. kk == 1 ) &
3544	READ ( 13 ) surf_h(1)%sasws
3545
3546	CASE ( 'surf_h(2)%start_index' )
3547	IF ( kk == 1 ) &
3548	READ ( 13 ) surf_h(2)%start_index
3549	l = 2
3550	CASE ( 'surf_h(2)%end_index' )
3551	IF ( kk == 1 ) &
3552	READ ( 13 ) surf_h(2)%end_index
3553	CASE ( 'surf_h(2)%us' )
3554	IF ( ALLOCATED( surf_h(2)%us ) .AND. kk == 1 ) &
3555	READ ( 13 ) surf_h(2)%us
3556	CASE ( 'surf_h(2)%ts' )
3557	IF ( ALLOCATED( surf_h(2)%ts ) .AND. kk == 1 ) &
3558	READ ( 13 ) surf_h(2)%ts
3559	CASE ( 'surf_h(2)%qs' )
3560	IF ( ALLOCATED( surf_h(2)%qs ) .AND. kk == 1 ) &
3561	READ ( 13 ) surf_h(2)%qs
3562	CASE ( 'surf_h(2)%ss' )
3563	IF ( ALLOCATED( surf_h(2)%ss ) .AND. kk == 1 ) &
3564	READ ( 13 ) surf_h(2)%ss
3565	CASE ( 'surf_h(2)%qcs' )
3566	IF ( ALLOCATED( surf_h(2)%qcs ) .AND. kk == 1 ) &
3567	READ ( 13 ) surf_h(2)%qcs
3568	CASE ( 'surf_h(2)%ncs' )
3569	IF ( ALLOCATED( surf_h(2)%ncs ) .AND. kk == 1 ) &
3570	READ ( 13 ) surf_h(2)%ncs
3571	CASE ( 'surf_h(2)%qrs' )
3572	IF ( ALLOCATED( surf_h(2)%qrs ) .AND. kk == 1 ) &
3573	READ ( 13 ) surf_h(2)%qrs
3574	CASE ( 'surf_h(2)%nrs' )
3575	IF ( ALLOCATED( surf_h(2)%nrs ) .AND. kk == 1 ) &
3576	READ ( 13 ) surf_h(2)%nrs
3577	CASE ( 'surf_h(2)%ol' )
3578	IF ( ALLOCATED( surf_h(2)%ol ) .AND. kk == 1 ) &
3579	READ ( 13 ) surf_h(2)%ol
3580	CASE ( 'surf_h(2)%rib' )
3581	IF ( ALLOCATED( surf_h(2)%rib ) .AND. kk == 1 ) &
3582	READ ( 13 ) surf_h(2)%rib
3583	CASE ( 'surf_h(2)%pt_surface' )
3584	IF ( ALLOCATED( surf_h(2)%pt_surface ) .AND. kk == 1 ) &
3585	READ ( 13 ) surf_h(2)%pt_surface
3586	CASE ( 'surf_h(2)%usws' )
3587	IF ( ALLOCATED( surf_h(2)%usws ) .AND. kk == 1 ) &
3588	READ ( 13 ) surf_h(2)%usws
3589	CASE ( 'surf_h(2)%vsws' )
3590	IF ( ALLOCATED( surf_h(2)%vsws ) .AND. kk == 1 ) &
3591	READ ( 13 ) surf_h(2)%vsws
3592	CASE ( 'surf_h(2)%shf' )
3593	IF ( ALLOCATED( surf_h(2)%shf ) .AND. kk == 1 ) &
3594	READ ( 13 ) surf_h(2)%shf
3595	CASE ( 'surf_h(2)%qsws' )
3596	IF ( ALLOCATED( surf_h(2)%qsws ) .AND. kk == 1 ) &
3597	READ ( 13 ) surf_h(2)%qsws
3598	CASE ( 'surf_h(2)%ssws' )
3599	IF ( ALLOCATED( surf_h(2)%ssws ) .AND. kk == 1 ) &
3600	READ ( 13 ) surf_h(2)%ssws
3601	CASE ( 'surf_h(2)%css' )
3602	IF ( ALLOCATED( surf_h(2)%css ) .AND. kk == 1 ) &
3603	READ ( 13 ) surf_h(2)%css
3604	CASE ( 'surf_h(2)%cssws' )
3605	IF ( ALLOCATED( surf_h(2)%cssws ) .AND. kk == 1 ) &
3606	READ ( 13 ) surf_h(2)%cssws
3607	CASE ( 'surf_h(2)%qcsws' )
3608	IF ( ALLOCATED( surf_h(2)%qcsws ) .AND. kk == 1 ) &
3609	READ ( 13 ) surf_h(2)%qcsws
3610	CASE ( 'surf_h(2)%ncsws' )
3611	IF ( ALLOCATED( surf_h(2)%ncsws ) .AND. kk == 1 ) &
3612	READ ( 13 ) surf_h(2)%ncsws
3613	CASE ( 'surf_h(2)%qrsws' )
3614	IF ( ALLOCATED( surf_h(2)%qrsws ) .AND. kk == 1 ) &
3615	READ ( 13 ) surf_h(2)%qrsws
3616	CASE ( 'surf_h(2)%nrsws' )
3617	IF ( ALLOCATED( surf_h(2)%nrsws ) .AND. kk == 1 ) &
3618	READ ( 13 ) surf_h(2)%nrsws
3619	CASE ( 'surf_h(2)%sasws' )
3620	IF ( ALLOCATED( surf_h(2)%sasws ) .AND. kk == 1 ) &
3621	READ ( 13 ) surf_h(2)%sasws
3622
3623	CASE ( 'surf_v(0)%start_index' )
3624	IF ( kk == 1 ) &
3625	READ ( 13 ) surf_v(0)%start_index
3626	l = 0
3627	horizontal_surface = .FALSE.
3628	vertical_surface = .TRUE.
3629	CASE ( 'surf_v(0)%end_index' )
3630	IF ( kk == 1 ) &
3631	READ ( 13 ) surf_v(0)%end_index
3632	CASE ( 'surf_v(0)%us' )
3633	IF ( ALLOCATED( surf_v(0)%us ) .AND. kk == 1 ) &
3634	READ ( 13 ) surf_v(0)%us
3635	CASE ( 'surf_v(0)%ts' )
3636	IF ( ALLOCATED( surf_v(0)%ts ) .AND. kk == 1 ) &
3637	READ ( 13 ) surf_v(0)%ts
3638	CASE ( 'surf_v(0)%qs' )
3639	IF ( ALLOCATED( surf_v(0)%qs ) .AND. kk == 1 ) &
3640	READ ( 13 ) surf_v(0)%qs
3641	CASE ( 'surf_v(0)%ss' )
3642	IF ( ALLOCATED( surf_v(0)%ss ) .AND. kk == 1 ) &
3643	READ ( 13 ) surf_v(0)%ss
3644	CASE ( 'surf_v(0)%qcs' )
3645	IF ( ALLOCATED( surf_v(0)%qcs ) .AND. kk == 1 ) &
3646	READ ( 13 ) surf_v(0)%qcs
3647	CASE ( 'surf_v(0)%ncs' )
3648	IF ( ALLOCATED( surf_v(0)%ncs ) .AND. kk == 1 ) &
3649	READ ( 13 ) surf_v(0)%ncs
3650	CASE ( 'surf_v(0)%qrs' )
3651	IF ( ALLOCATED( surf_v(0)%qrs ) .AND. kk == 1 ) &
3652	READ ( 13 ) surf_v(0)%qrs
3653	CASE ( 'surf_v(0)%nrs' )
3654	IF ( ALLOCATED( surf_v(0)%nrs ) .AND. kk == 1 ) &
3655	READ ( 13 ) surf_v(0)%nrs
3656	CASE ( 'surf_v(0)%ol' )
3657	IF ( ALLOCATED( surf_v(0)%ol ) .AND. kk == 1 ) &
3658	READ ( 13 ) surf_v(0)%ol
3659	CASE ( 'surf_v(0)%rib' )
3660	IF ( ALLOCATED( surf_v(0)%rib ) .AND. kk == 1 ) &
3661	READ ( 13 ) surf_v(0)%rib
3662	CASE ( 'surf_v(0)%pt_surface' )
3663	IF ( ALLOCATED( surf_v(0)%pt_surface ) .AND. kk == 1 ) &
3664	READ ( 13 ) surf_v(0)%pt_surface
3665	CASE ( 'surf_v(0)%shf' )
3666	IF ( ALLOCATED( surf_v(0)%shf ) .AND. kk == 1 ) &
3667	READ ( 13 ) surf_v(0)%shf
3668	CASE ( 'surf_v(0)%qsws' )
3669	IF ( ALLOCATED( surf_v(0)%qsws ) .AND. kk == 1 ) &
3670	READ ( 13 ) surf_v(0)%qsws
3671	CASE ( 'surf_v(0)%ssws' )
3672	IF ( ALLOCATED( surf_v(0)%ssws ) .AND. kk == 1 ) &
3673	READ ( 13 ) surf_v(0)%ssws
3674	CASE ( 'surf_v(0)%css' )
3675	IF ( ALLOCATED( surf_v(0)%css ) .AND. kk == 1 ) &
3676	READ ( 13 ) surf_v(0)%css
3677	CASE ( 'surf_v(0)%cssws' )
3678	IF ( ALLOCATED( surf_v(0)%cssws ) .AND. kk == 1 ) &
3679	READ ( 13 ) surf_v(0)%cssws
3680	CASE ( 'surf_v(0)%qcsws' )
3681	IF ( ALLOCATED( surf_v(0)%qcsws ) .AND. kk == 1 ) &
3682	READ ( 13 ) surf_v(0)%qcsws
3683	CASE ( 'surf_v(0)%ncsws' )
3684	IF ( ALLOCATED( surf_v(0)%ncsws ) .AND. kk == 1 ) &
3685	READ ( 13 ) surf_v(0)%ncsws
3686	CASE ( 'surf_v(0)%qrsws' )
3687	IF ( ALLOCATED( surf_v(0)%qrsws ) .AND. kk == 1 ) &
3688	READ ( 13 ) surf_v(0)%qrsws
3689	CASE ( 'surf_v(0)%nrsws' )
3690	IF ( ALLOCATED( surf_v(0)%nrsws ) .AND. kk == 1 ) &
3691	READ ( 13 ) surf_v(0)%nrsws
3692	CASE ( 'surf_v(0)%sasws' )
3693	IF ( ALLOCATED( surf_v(0)%sasws ) .AND. kk == 1 ) &
3694	READ ( 13 ) surf_v(0)%sasws
3695	CASE ( 'surf_v(0)%mom_uv' )
3696	IF ( ALLOCATED( surf_v(0)%mom_flux_uv ) .AND. kk == 1 ) &
3697	READ ( 13 ) surf_v(0)%mom_flux_uv
3698	CASE ( 'surf_v(0)%mom_w' )
3699	IF ( ALLOCATED( surf_v(0)%mom_flux_w ) .AND. kk == 1 ) &
3700	READ ( 13 ) surf_v(0)%mom_flux_w
3701	CASE ( 'surf_v(0)%mom_tke' )
3702	IF ( ALLOCATED( surf_v(0)%mom_flux_tke ) .AND. kk == 1 ) &
3703	READ ( 13 ) surf_v(0)%mom_flux_tke
3704
3705	CASE ( 'surf_v(1)%start_index' )
3706	IF ( kk == 1 ) &
3707	READ ( 13 ) surf_v(1)%start_index
3708	l = 1
3709	CASE ( 'surf_v(1)%end_index' )
3710	IF ( kk == 1 ) &
3711	READ ( 13 ) surf_v(1)%end_index
3712	CASE ( 'surf_v(1)%us' )
3713	IF ( ALLOCATED( surf_v(1)%us ) .AND. kk == 1 ) &
3714	READ ( 13 ) surf_v(1)%us
3715	CASE ( 'surf_v(1)%ts' )
3716	IF ( ALLOCATED( surf_v(1)%ts ) .AND. kk == 1 ) &
3717	READ ( 13 ) surf_v(1)%ts
3718	CASE ( 'surf_v(1)%qs' )
3719	IF ( ALLOCATED( surf_v(1)%qs ) .AND. kk == 1 ) &
3720	READ ( 13 ) surf_v(1)%qs
3721	CASE ( 'surf_v(1)%ss' )
3722	IF ( ALLOCATED( surf_v(1)%ss ) .AND. kk == 1 ) &
3723	READ ( 13 ) surf_v(1)%ss
3724	CASE ( 'surf_v(1)%qcs' )
3725	IF ( ALLOCATED( surf_v(1)%qcs ) .AND. kk == 1 ) &
3726	READ ( 13 ) surf_v(1)%qcs
3727	CASE ( 'surf_v(1)%ncs' )
3728	IF ( ALLOCATED( surf_v(1)%ncs ) .AND. kk == 1 ) &
3729	READ ( 13 ) surf_v(1)%ncs
3730	CASE ( 'surf_v(1)%qrs' )
3731	IF ( ALLOCATED( surf_v(1)%qrs ) .AND. kk == 1 ) &
3732	READ ( 13 ) surf_v(1)%qrs
3733	CASE ( 'surf_v(1)%nrs' )
3734	IF ( ALLOCATED( surf_v(1)%nrs ) .AND. kk == 1 ) &
3735	READ ( 13 ) surf_v(1)%nrs
3736	CASE ( 'surf_v(1)%ol' )
3737	IF ( ALLOCATED( surf_v(1)%ol ) .AND. kk == 1 ) &
3738	READ ( 13 ) surf_v(1)%ol
3739	CASE ( 'surf_v(1)%rib' )
3740	IF ( ALLOCATED( surf_v(1)%rib ) .AND. kk == 1 ) &
3741	READ ( 13 ) surf_v(1)%rib
3742	CASE ( 'surf_v(1)%pt_surface' )
3743	IF ( ALLOCATED( surf_v(1)%pt_surface ) .AND. kk == 1 ) &
3744	READ ( 13 ) surf_v(1)%pt_surface
3745	CASE ( 'surf_v(1)%shf' )
3746	IF ( ALLOCATED( surf_v(1)%shf ) .AND. kk == 1 ) &
3747	READ ( 13 ) surf_v(1)%shf
3748	CASE ( 'surf_v(1)%qsws' )
3749	IF ( ALLOCATED( surf_v(1)%qsws ) .AND. kk == 1 ) &
3750	READ ( 13 ) surf_v(1)%qsws
3751	CASE ( 'surf_v(1)%ssws' )
3752	IF ( ALLOCATED( surf_v(1)%ssws ) .AND. kk == 1 ) &
3753	READ ( 13 ) surf_v(1)%ssws
3754	CASE ( 'surf_v(1)%css' )
3755	IF ( ALLOCATED( surf_v(1)%css ) .AND. kk == 1 ) &
3756	READ ( 13 ) surf_v(1)%css
3757	CASE ( 'surf_v(1)%cssws' )
3758	IF ( ALLOCATED( surf_v(1)%cssws ) .AND. kk == 1 ) &
3759	READ ( 13 ) surf_v(1)%cssws
3760	CASE ( 'surf_v(1)%qcsws' )
3761	IF ( ALLOCATED( surf_v(1)%qcsws ) .AND. kk == 1 ) &
3762	READ ( 13 ) surf_v(1)%qcsws
3763	CASE ( 'surf_v(1)%ncsws' )
3764	IF ( ALLOCATED( surf_v(1)%ncsws ) .AND. kk == 1 ) &
3765	READ ( 13 ) surf_v(1)%ncsws
3766	CASE ( 'surf_v(1)%qrsws' )
3767	IF ( ALLOCATED( surf_v(1)%qrsws ) .AND. kk == 1 ) &
3768	READ ( 13 ) surf_v(1)%qrsws
3769	CASE ( 'surf_v(1)%nrsws' )
3770	IF ( ALLOCATED( surf_v(1)%nrsws ) .AND. kk == 1 ) &
3771	READ ( 13 ) surf_v(1)%nrsws
3772	CASE ( 'surf_v(1)%sasws' )
3773	IF ( ALLOCATED( surf_v(1)%sasws ) .AND. kk == 1 ) &
3774	READ ( 13 ) surf_v(1)%sasws
3775	CASE ( 'surf_v(1)%mom_uv' )
3776	IF ( ALLOCATED( surf_v(1)%mom_flux_uv ) .AND. kk == 1 ) &
3777	READ ( 13 ) surf_v(1)%mom_flux_uv
3778	CASE ( 'surf_v(1)%mom_w' )
3779	IF ( ALLOCATED( surf_v(1)%mom_flux_w ) .AND. kk == 1 ) &
3780	READ ( 13 ) surf_v(1)%mom_flux_w
3781	CASE ( 'surf_v(1)%mom_tke' )
3782	IF ( ALLOCATED( surf_v(1)%mom_flux_tke ) .AND. kk == 1 ) &
3783	READ ( 13 ) surf_v(1)%mom_flux_tke
3784
3785	CASE ( 'surf_v(2)%start_index' )
3786	IF ( kk == 1 ) &
3787	READ ( 13 ) surf_v(2)%start_index
3788	l = 2
3789	CASE ( 'surf_v(2)%end_index' )
3790	IF ( kk == 1 ) &
3791	READ ( 13 ) surf_v(2)%end_index
3792	CASE ( 'surf_v(2)%us' )
3793	IF ( ALLOCATED( surf_v(2)%us ) .AND. kk == 1 ) &
3794	READ ( 13 ) surf_v(2)%us
3795	CASE ( 'surf_v(2)%ts' )
3796	IF ( ALLOCATED( surf_v(2)%ts ) .AND. kk == 1 ) &
3797	READ ( 13 ) surf_v(2)%ts
3798	CASE ( 'surf_v(2)%qs' )
3799	IF ( ALLOCATED( surf_v(2)%qs ) .AND. kk == 1 ) &
3800	READ ( 13 ) surf_v(2)%qs
3801	CASE ( 'surf_v(2)%ss' )
3802	IF ( ALLOCATED( surf_v(2)%ss ) .AND. kk == 1 ) &
3803	READ ( 13 ) surf_v(2)%ss
3804	CASE ( 'surf_v(2)%qcs' )
3805	IF ( ALLOCATED( surf_v(2)%qcs ) .AND. kk == 1 ) &
3806	READ ( 13 ) surf_v(2)%qcs
3807	CASE ( 'surf_v(2)%ncs' )
3808	IF ( ALLOCATED( surf_v(2)%ncs ) .AND. kk == 1 ) &
3809	READ ( 13 ) surf_v(2)%ncs
3810	CASE ( 'surf_v(2)%qrs' )
3811	IF ( ALLOCATED( surf_v(2)%qrs ) .AND. kk == 1 ) &
3812	READ ( 13 ) surf_v(2)%qrs
3813	CASE ( 'surf_v(2)%nrs' )
3814	IF ( ALLOCATED( surf_v(2)%nrs ) .AND. kk == 1 ) &
3815	READ ( 13 ) surf_v(2)%nrs
3816	CASE ( 'surf_v(2)%ol' )
3817	IF ( ALLOCATED( surf_v(2)%ol ) .AND. kk == 1 ) &
3818	READ ( 13 ) surf_v(2)%ol
3819	CASE ( 'surf_v(2)%rib' )
3820	IF ( ALLOCATED( surf_v(2)%rib ) .AND. kk == 1 ) &
3821	READ ( 13 ) surf_v(2)%rib
3822	CASE ( 'surf_v(2)%pt_surface' )
3823	IF ( ALLOCATED( surf_v(2)%pt_surface ) .AND. kk == 1 ) &
3824	READ ( 13 ) surf_v(2)%pt_surface
3825	CASE ( 'surf_v(2)%shf' )
3826	IF ( ALLOCATED( surf_v(2)%shf ) .AND. kk == 1 ) &
3827	READ ( 13 ) surf_v(2)%shf
3828	CASE ( 'surf_v(2)%qsws' )
3829	IF ( ALLOCATED( surf_v(2)%qsws ) .AND. kk == 1 ) &
3830	READ ( 13 ) surf_v(2)%qsws
3831	CASE ( 'surf_v(2)%ssws' )
3832	IF ( ALLOCATED( surf_v(2)%ssws ) .AND. kk == 1 ) &
3833	READ ( 13 ) surf_v(2)%ssws
3834	CASE ( 'surf_v(2)%css' )
3835	IF ( ALLOCATED( surf_v(2)%css ) .AND. kk == 1 ) &
3836	READ ( 13 ) surf_v(2)%css
3837	CASE ( 'surf_v(2)%cssws' )
3838	IF ( ALLOCATED( surf_v(2)%cssws ) .AND. kk == 1 ) &
3839	READ ( 13 ) surf_v(2)%cssws
3840	CASE ( 'surf_v(2)%qcsws' )
3841	IF ( ALLOCATED( surf_v(2)%qcsws ) .AND. kk == 1 ) &
3842	READ ( 13 ) surf_v(2)%qcsws
3843	CASE ( 'surf_v(2)%ncsws' )
3844	IF ( ALLOCATED( surf_v(2)%ncsws ) .AND. kk == 1 ) &
3845	READ ( 13 ) surf_v(2)%ncsws
3846	CASE ( 'surf_v(2)%qrsws' )
3847	IF ( ALLOCATED( surf_v(2)%qrsws ) .AND. kk == 1 ) &
3848	READ ( 13 ) surf_v(2)%qrsws
3849	CASE ( 'surf_v(2)%nrsws' )
3850	IF ( ALLOCATED( surf_v(2)%nrsws ) .AND. kk == 1 ) &
3851	READ ( 13 ) surf_v(2)%nrsws
3852	CASE ( 'surf_v(2)%sasws' )
3853	IF ( ALLOCATED( surf_v(2)%sasws ) .AND. kk == 1 ) &
3854	READ ( 13 ) surf_v(2)%sasws
3855	CASE ( 'surf_v(2)%mom_uv' )
3856	IF ( ALLOCATED( surf_v(2)%mom_flux_uv ) .AND. kk == 1 ) &
3857	READ ( 13 ) surf_v(2)%mom_flux_uv
3858	CASE ( 'surf_v(2)%mom_w' )
3859	IF ( ALLOCATED( surf_v(2)%mom_flux_w ) .AND. kk == 1 ) &
3860	READ ( 13 ) surf_v(2)%mom_flux_w
3861	CASE ( 'surf_v(2)%mom_tke' )
3862	IF ( ALLOCATED( surf_v(2)%mom_flux_tke ) .AND. kk == 1 ) &
3863	READ ( 13 ) surf_v(2)%mom_flux_tke
3864
3865	CASE ( 'surf_v(3)%start_index' )
3866	IF ( kk == 1 ) &
3867	READ ( 13 ) surf_v(3)%start_index
3868	l = 3
3869	CASE ( 'surf_v(3)%end_index' )
3870	IF ( kk == 1 ) &
3871	READ ( 13 ) surf_v(3)%end_index
3872	CASE ( 'surf_v(3)%us' )
3873	IF ( ALLOCATED( surf_v(3)%us ) .AND. kk == 1 ) &
3874	READ ( 13 ) surf_v(3)%us
3875	CASE ( 'surf_v(3)%ts' )
3876	IF ( ALLOCATED( surf_v(3)%ts ) .AND. kk == 1 ) &
3877	READ ( 13 ) surf_v(3)%ts
3878	CASE ( 'surf_v(3)%qs' )
3879	IF ( ALLOCATED( surf_v(3)%qs ) .AND. kk == 1 ) &
3880	READ ( 13 ) surf_v(3)%qs
3881	CASE ( 'surf_v(3)%ss' )
3882	IF ( ALLOCATED( surf_v(3)%ss ) .AND. kk == 1 ) &
3883	READ ( 13 ) surf_v(3)%ss
3884	CASE ( 'surf_v(3)%qcs' )
3885	IF ( ALLOCATED( surf_v(3)%qcs ) .AND. kk == 1 ) &
3886	READ ( 13 ) surf_v(3)%qcs
3887	CASE ( 'surf_v(3)%ncs' )
3888	IF ( ALLOCATED( surf_v(3)%ncs ) .AND. kk == 1 ) &
3889	READ ( 13 ) surf_v(3)%ncs
3890	CASE ( 'surf_v(3)%qrs' )
3891	IF ( ALLOCATED( surf_v(3)%qrs ) .AND. kk == 1 ) &
3892	READ ( 13 ) surf_v(3)%qrs
3893	CASE ( 'surf_v(3)%nrs' )
3894	IF ( ALLOCATED( surf_v(3)%nrs ) .AND. kk == 1 ) &
3895	READ ( 13 ) surf_v(3)%nrs
3896	CASE ( 'surf_v(3)%ol' )
3897	IF ( ALLOCATED( surf_v(3)%ol ) .AND. kk == 1 ) &
3898	READ ( 13 ) surf_v(3)%ol
3899	CASE ( 'surf_v(3)%rib' )
3900	IF ( ALLOCATED( surf_v(3)%rib ) .AND. kk == 1 ) &
3901	READ ( 13 ) surf_v(3)%rib
3902	CASE ( 'surf_v(3)%pt_surface' )
3903	IF ( ALLOCATED( surf_v(3)%pt_surface ) .AND. kk == 1 ) &
3904	READ ( 13 ) surf_v(3)%pt_surface
3905	CASE ( 'surf_v(3)%shf' )
3906	IF ( ALLOCATED( surf_v(3)%shf ) .AND. kk == 1 ) &
3907	READ ( 13 ) surf_v(3)%shf
3908	CASE ( 'surf_v(3)%qsws' )
3909	IF ( ALLOCATED( surf_v(3)%qsws ) .AND. kk == 1 ) &
3910	READ ( 13 ) surf_v(3)%qsws
3911	CASE ( 'surf_v(3)%ssws' )
3912	IF ( ALLOCATED( surf_v(3)%ssws ) .AND. kk == 1 ) &
3913	READ ( 13 ) surf_v(3)%ssws
3914	CASE ( 'surf_v(3)%css' )
3915	IF ( ALLOCATED( surf_v(3)%css ) .AND. kk == 1 ) &
3916	READ ( 13 ) surf_v(3)%css
3917	CASE ( 'surf_v(3)%cssws' )
3918	IF ( ALLOCATED( surf_v(3)%cssws ) .AND. kk == 1 ) &
3919	READ ( 13 ) surf_v(3)%cssws
3920	CASE ( 'surf_v(3)%qcsws' )
3921	IF ( ALLOCATED( surf_v(3)%qcsws ) .AND. kk == 1 ) &
3922	READ ( 13 ) surf_v(3)%qcsws
3923	CASE ( 'surf_v(3)%ncsws' )
3924	IF ( ALLOCATED( surf_v(3)%ncsws ) .AND. kk == 1 ) &
3925	READ ( 13 ) surf_v(3)%ncsws
3926	CASE ( 'surf_v(3)%qrsws' )
3927	IF ( ALLOCATED( surf_v(3)%qrsws ) .AND. kk == 1 ) &
3928	READ ( 13 ) surf_v(3)%qrsws
3929	CASE ( 'surf_v(3)%nrsws' )
3930	IF ( ALLOCATED( surf_v(3)%nrsws ) .AND. kk == 1 ) &
3931	READ ( 13 ) surf_v(3)%nrsws
3932	CASE ( 'surf_v(3)%sasws' )
3933	IF ( ALLOCATED( surf_v(3)%sasws ) .AND. kk == 1 ) &
3934	READ ( 13 ) surf_v(3)%sasws
3935	CASE ( 'surf_v(3)%mom_uv' )
3936	IF ( ALLOCATED( surf_v(3)%mom_flux_uv ) .AND. kk == 1 ) &
3937	READ ( 13 ) surf_v(3)%mom_flux_uv
3938	CASE ( 'surf_v(3)%mom_w' )
3939	IF ( ALLOCATED( surf_v(3)%mom_flux_w ) .AND. kk == 1 ) &
3940	READ ( 13 ) surf_v(3)%mom_flux_w
3941	CASE ( 'surf_v(3)%mom_tke' )
3942	IF ( ALLOCATED( surf_v(3)%mom_flux_tke ) .AND. kk == 1 ) &
3943	READ ( 13 ) surf_v(3)%mom_flux_tke
3944
3945	CASE DEFAULT
3946
3947	found = .FALSE.
3948
3949	END SELECT
3950	!
3951	!-- Redistribute surface elements on its respective type.
3952	IF ( horizontal_surface .AND. &
3953	.NOT. INDEX( restart_string(1:length), '%start_index' ) /= 0 ) &
3954	THEN
3955
3956	ic = nxlc
3957	DO i = nxlf, nxrf
3958	jc = nysc
3959	DO j = nysf, nynf
3960
3961	surf_match_def = surf_def_h(l)%end_index(jc,ic) >= &
3962	surf_def_h(l)%start_index(jc,ic)
3963	surf_match_lsm = ( surf_lsm_h%end_index(jc,ic) >= &
3964	surf_lsm_h%start_index(jc,ic) ) &
3965	.AND. l == 0
3966	surf_match_usm = ( surf_usm_h%end_index(jc,ic) >= &
3967	surf_usm_h%start_index(jc,ic) ) &
3968	.AND. l == 0
3969
3970	IF ( surf_match_def ) THEN
3971	mm = surf_def_h(l)%start_index(jc,ic)
3972	DO m = surf_h(l)%start_index(j,i), &
3973	surf_h(l)%end_index(j,i)
3974	IF ( surf_def_h(l)%end_index(jc,ic) >= mm ) &
3975	CALL restore_surface_elements( surf_def_h(l), &
3976	mm, surf_h(l), m )
3977	mm = mm + 1
3978	ENDDO
3979	ENDIF
3980
3981	IF ( surf_match_lsm ) THEN
3982	mm = surf_lsm_h%start_index(jc,ic)
3983	DO m = surf_h(l)%start_index(j,i), &
3984	surf_h(l)%end_index(j,i)
3985	IF ( surf_lsm_h%end_index(jc,ic) >= mm ) &
3986	CALL restore_surface_elements( surf_lsm_h, &
3987	mm, surf_h(l), m )
3988	mm = mm + 1
3989	ENDDO
3990	ENDIF
3991
3992	IF ( surf_match_usm ) THEN
3993	mm = surf_usm_h%start_index(jc,ic)
3994	DO m = surf_h(l)%start_index(j,i), &
3995	surf_h(l)%end_index(j,i)
3996	IF ( surf_usm_h%end_index(jc,ic) >= mm ) &
3997	CALL restore_surface_elements( surf_usm_h, &
3998	mm, surf_h(l), m )
3999	mm = mm + 1
4000	ENDDO
4001	ENDIF
4002
4003	jc = jc + 1
4004	ENDDO
4005	ic = ic + 1
4006	ENDDO
4007	ELSEIF ( vertical_surface .AND. &
4008	.NOT. INDEX( restart_string(1:length), '%start_index' ) /= 0 ) &
4009	THEN
4010	ic = nxlc
4011	DO i = nxlf, nxrf
4012	jc = nysc
4013	DO j = nysf, nynf
4014
4015	surf_match_def = surf_def_v(l)%end_index(jc,ic) >= &
4016	surf_def_v(l)%start_index(jc,ic)
4017	surf_match_lsm = surf_lsm_v(l)%end_index(jc,ic) >= &
4018	surf_lsm_v(l)%start_index(jc,ic)
4019	surf_match_usm = surf_usm_v(l)%end_index(jc,ic) >= &
4020	surf_usm_v(l)%start_index(jc,ic)
4021
4022	IF ( surf_match_def ) THEN
4023	mm = surf_def_v(l)%start_index(jc,ic)
4024	DO m = surf_v(l)%start_index(j,i), &
4025	surf_v(l)%end_index(j,i)
4026	IF ( surf_def_v(l)%end_index(jc,ic) >= mm ) &
4027	CALL restore_surface_elements( surf_def_v(l), &
4028	mm, surf_v(l), m )
4029	mm = mm + 1
4030	ENDDO
4031	ENDIF
4032
4033	IF ( surf_match_lsm ) THEN
4034	mm = surf_lsm_v(l)%start_index(jc,ic)
4035	DO m = surf_v(l)%start_index(j,i), &
4036	surf_v(l)%end_index(j,i)
4037	IF ( surf_lsm_v(l)%end_index(jc,ic) >= mm ) &
4038	CALL restore_surface_elements( surf_lsm_v(l), &
4039	mm, surf_v(l), m )
4040	mm = mm + 1
4041	ENDDO
4042	ENDIF
4043
4044	IF ( surf_match_usm ) THEN
4045	mm = surf_usm_v(l)%start_index(jc,ic)
4046	DO m = surf_v(l)%start_index(j,i), &
4047	surf_v(l)%end_index(j,i)
4048	IF ( surf_usm_v(l)%end_index(jc,ic) >= mm ) &
4049	CALL restore_surface_elements( surf_usm_v(l), &
4050	mm, surf_v(l), m )
4051	mm = mm + 1
4052	ENDDO
4053	ENDIF
4054
4055	jc = jc + 1
4056	ENDDO
4057	ic = ic + 1
4058	ENDDO
4059	ENDIF
4060
4061
4062	CONTAINS
4063	!------------------------------------------------------------------------------!
4064	! Description:
4065	! ------------
4066	!> Restores surfacle elements back on its respective type.
4067	!------------------------------------------------------------------------------!
4068	SUBROUTINE restore_surface_elements( surf_target, m_target, &
4069	surf_file, m_file )
4070
4071	IMPLICIT NONE
4072
4073	INTEGER(iwp) :: m_file !< respective surface-element index of current surface array
4074	INTEGER(iwp) :: m_target !< respecitve surface-element index of surface array on file
4075	INTEGER(iwp) :: lsp !< running index chemical species
4076
4077	TYPE( surf_type ) :: surf_target !< target surface type
4078	TYPE( surf_type ) :: surf_file !< surface type on file
4079
4080
4081	IF ( INDEX( restart_string(1:length), '%us' ) /= 0 ) THEN
4082	IF ( ALLOCATED( surf_target%us ) .AND. &
4083	ALLOCATED( surf_file%us ) ) &
4084	surf_target%us(m_target) = surf_file%us(m_file)
4085	ENDIF
4086
4087	IF ( INDEX( restart_string(1:length), '%ol' ) /= 0 ) THEN
4088	IF ( ALLOCATED( surf_target%ol ) .AND. &
4089	ALLOCATED( surf_file%ol ) ) &
4090	surf_target%ol(m_target) = surf_file%ol(m_file)
4091	ENDIF
4092
4093	IF ( INDEX( restart_string(1:length), '%pt_surface' ) /= 0 ) THEN
4094	IF ( ALLOCATED( surf_target%pt_surface ) .AND. &
4095	ALLOCATED( surf_file%pt_surface ) ) &
4096	surf_target%pt_surface(m_target) = surf_file%pt_surface(m_file)
4097	ENDIF
4098
4099	IF ( INDEX( restart_string(1:length), '%usws' ) /= 0 ) THEN
4100	IF ( ALLOCATED( surf_target%usws ) .AND. &
4101	ALLOCATED( surf_file%usws ) ) &
4102	surf_target%usws(m_target) = surf_file%usws(m_file)
4103	ENDIF
4104
4105	IF ( INDEX( restart_string(1:length), '%vsws' ) /= 0 ) THEN
4106	IF ( ALLOCATED( surf_target%vsws ) .AND. &
4107	ALLOCATED( surf_file%vsws ) ) &
4108	surf_target%vsws(m_target) = surf_file%vsws(m_file)
4109	ENDIF
4110
4111	IF ( INDEX( restart_string(1:length), '%ts' ) /= 0 ) THEN
4112	IF ( ALLOCATED( surf_target%ts ) .AND. &
4113	ALLOCATED( surf_file%ts ) ) &
4114	surf_target%ts(m_target) = surf_file%ts(m_file)
4115	ENDIF
4116
4117	IF ( INDEX( restart_string(1:length), '%shf' ) /= 0 ) THEN
4118	IF ( ALLOCATED( surf_target%shf ) .AND. &
4119	ALLOCATED( surf_file%shf ) ) &
4120	surf_target%shf(m_target) = surf_file%shf(m_file)
4121	ENDIF
4122
4123	IF ( INDEX( restart_string(1:length), '%qs' ) /= 0 ) THEN
4124	IF ( ALLOCATED( surf_target%qs ) .AND. &
4125	ALLOCATED( surf_file%qs ) ) &
4126	surf_target%qs(m_target) = surf_file%qs(m_file)
4127	ENDIF
4128
4129	IF ( INDEX( restart_string(1:length), '%qsws' ) /= 0 ) THEN
4130	IF ( ALLOCATED( surf_target%qsws ) .AND. &
4131	ALLOCATED( surf_file%qsws ) ) &
4132	surf_target%qsws(m_target) = surf_file%qsws(m_file)
4133	ENDIF
4134
4135	IF ( INDEX( restart_string(1:length), '%ss' ) /= 0 ) THEN
4136	IF ( ALLOCATED( surf_target%ss ) .AND. &
4137	ALLOCATED( surf_file%ss ) ) &
4138	surf_target%ss(m_target) = surf_file%ss(m_file)
4139	ENDIF
4140
4141	IF ( INDEX( restart_string(1:length), '%ssws' ) /= 0 ) THEN
4142	IF ( ALLOCATED( surf_target%ssws ) .AND. &
4143	ALLOCATED( surf_file%ssws ) ) &
4144	surf_target%ssws(m_target) = surf_file%ssws(m_file)
4145	ENDIF
4146
4147	IF ( INDEX( restart_string(1:length), '%css' ) /= 0 ) THEN
4148	IF ( ALLOCATED( surf_target%css ) .AND. &
4149	ALLOCATED( surf_file%css ) ) THEN
4150	DO lsp = 1, nvar
4151	surf_target%css(lsp,m_target) = surf_file%css(lsp,m_file)
4152	ENDDO
4153	ENDIF
4154	ENDIF
4155	IF ( INDEX( restart_string(1:length), '%cssws' ) /= 0 ) THEN
4156	IF ( ALLOCATED( surf_target%cssws ) .AND. &
4157	ALLOCATED( surf_file%cssws ) ) THEN
4158	DO lsp = 1, nvar
4159	surf_target%cssws(lsp,m_target) = surf_file%cssws(lsp,m_file)
4160	ENDDO
4161	ENDIF
4162	ENDIF
4163
4164	IF ( INDEX( restart_string(1:length), '%qcs' ) /= 0 ) THEN
4165	IF ( ALLOCATED( surf_target%qcs ) .AND. &
4166	ALLOCATED( surf_file%qcs ) ) &
4167	surf_target%qcs(m_target) = surf_file%qcs(m_file)
4168	ENDIF
4169
4170	IF ( INDEX( restart_string(1:length), '%qcsws' ) /= 0 ) THEN
4171	IF ( ALLOCATED( surf_target%qcsws ) .AND. &
4172	ALLOCATED( surf_file%qcsws ) ) &
4173	surf_target%qcsws(m_target) = surf_file%qcsws(m_file)
4174	ENDIF
4175
4176	IF ( INDEX( restart_string(1:length), '%ncs' ) /= 0 ) THEN
4177	IF ( ALLOCATED( surf_target%ncs ) .AND. &
4178	ALLOCATED( surf_file%ncs ) ) &
4179	surf_target%ncs(m_target) = surf_file%ncs(m_file)
4180	ENDIF
4181
4182	IF ( INDEX( restart_string(1:length), '%ncsws' ) /= 0 ) THEN
4183	IF ( ALLOCATED( surf_target%ncsws ) .AND. &
4184	ALLOCATED( surf_file%ncsws ) ) &
4185	surf_target%ncsws(m_target) = surf_file%ncsws(m_file)
4186	ENDIF
4187
4188	IF ( INDEX( restart_string(1:length), '%qrs' ) /= 0 ) THEN
4189	IF ( ALLOCATED( surf_target%qrs ) .AND. &
4190	ALLOCATED( surf_file%qrs ) ) &
4191	surf_target%qrs(m_target) = surf_file%qrs(m_file)
4192	ENDIF
4193
4194	IF ( INDEX( restart_string(1:length), '%qrsws' ) /= 0 ) THEN
4195	IF ( ALLOCATED( surf_target%qrsws ) .AND. &
4196	ALLOCATED( surf_file%qrsws ) ) &
4197	surf_target%qrsws(m_target) = surf_file%qrsws(m_file)
4198	ENDIF
4199
4200	IF ( INDEX( restart_string(1:length), '%nrs' ) /= 0 ) THEN
4201	IF ( ALLOCATED( surf_target%nrs ) .AND. &
4202	ALLOCATED( surf_file%nrs ) ) &
4203	surf_target%nrs(m_target) = surf_file%nrs(m_file)
4204	ENDIF
4205
4206	IF ( INDEX( restart_string(1:length), '%nrsws' ) /= 0 ) THEN
4207	IF ( ALLOCATED( surf_target%nrsws ) .AND. &
4208	ALLOCATED( surf_file%nrsws ) ) &
4209	surf_target%nrsws(m_target) = surf_file%nrsws(m_file)
4210	ENDIF
4211
4212	IF ( INDEX( restart_string(1:length), '%sasws' ) /= 0 ) THEN
4213	IF ( ALLOCATED( surf_target%sasws ) .AND. &
4214	ALLOCATED( surf_file%sasws ) ) &
4215	surf_target%sasws(m_target) = surf_file%sasws(m_file)
4216	ENDIF
4217
4218	IF ( INDEX( restart_string(1:length), '%mom_uv' ) /= 0 ) THEN
4219	IF ( ALLOCATED( surf_target%mom_flux_uv ) .AND. &
4220	ALLOCATED( surf_file%mom_flux_uv ) ) &
4221	surf_target%mom_flux_uv(m_target) = &
4222	surf_file%mom_flux_uv(m_file)
4223	ENDIF
4224
4225	IF ( INDEX( restart_string(1:length), '%mom_w' ) /= 0 ) THEN
4226	IF ( ALLOCATED( surf_target%mom_flux_w ) .AND. &
4227	ALLOCATED( surf_file%mom_flux_w ) ) &
4228	surf_target%mom_flux_w(m_target) = &
4229	surf_file%mom_flux_w(m_file)
4230	ENDIF
4231
4232	IF ( INDEX( restart_string(1:length), '%mom_tke' ) /= 0 ) THEN
4233	IF ( ALLOCATED( surf_target%mom_flux_tke ) .AND. &
4234	ALLOCATED( surf_file%mom_flux_tke ) ) &
4235	surf_target%mom_flux_tke(0:1,m_target) = &
4236	surf_file%mom_flux_tke(0:1,m_file)
4237	ENDIF
4238
4239
4240	END SUBROUTINE restore_surface_elements
4241
4242
4243	END SUBROUTINE surface_rrd_local
4244
4245
4246	!------------------------------------------------------------------------------!
4247	! Description:
4248	! ------------
4249	!> Counts the number of surface elements with the same facing, required for
4250	!> reading and writing restart data.
4251	!------------------------------------------------------------------------------!
4252	SUBROUTINE surface_last_actions
4253
4254	IMPLICIT NONE
4255	!
4256	!-- Horizontal surfaces
4257	ns_h_on_file(0) = surf_def_h(0)%ns + surf_lsm_h%ns + surf_usm_h%ns
4258	ns_h_on_file(1) = surf_def_h(1)%ns
4259	ns_h_on_file(2) = surf_def_h(2)%ns
4260	!
4261	!-- Vertical surfaces
4262	ns_v_on_file(0) = surf_def_v(0)%ns + surf_lsm_v(0)%ns + surf_usm_v(0)%ns
4263	ns_v_on_file(1) = surf_def_v(1)%ns + surf_lsm_v(1)%ns + surf_usm_v(1)%ns
4264	ns_v_on_file(2) = surf_def_v(2)%ns + surf_lsm_v(2)%ns + surf_usm_v(2)%ns
4265	ns_v_on_file(3) = surf_def_v(3)%ns + surf_lsm_v(3)%ns + surf_usm_v(3)%ns
4266
4267	END SUBROUTINE surface_last_actions
4268
4269	!------------------------------------------------------------------------------!
4270	! Description:
4271	! ------------
4272	!> Routine maps surface data read from file after restart - 1D arrays
4273	!------------------------------------------------------------------------------!
4274	SUBROUTINE surface_restore_elements_1d( surf_target, surf_file, &
4275	start_index_c, &
4276	start_index_on_file, &
4277	end_index_on_file, &
4278	nxlc, nysc, nxlf, nxrf, nysf, nynf,&
4279	nys_on_file, nyn_on_file, &
4280	nxl_on_file,nxr_on_file )
4281
4282	IMPLICIT NONE
4283
4284	INTEGER(iwp) :: i !< running index along x-direction, refers to former domain size
4285	INTEGER(iwp) :: ic !< running index along x-direction, refers to current domain size
4286	INTEGER(iwp) :: j !< running index along y-direction, refers to former domain size
4287	INTEGER(iwp) :: jc !< running index along y-direction, refers to former domain size
4288	INTEGER(iwp) :: m !< surface-element index on file
4289	INTEGER(iwp) :: mm !< surface-element index on current subdomain
4290	INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain
4291	INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain
4292	INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain
4293	INTEGER(iwp) :: nysc !< index of north boundary on current subdomain
4294	INTEGER(iwp) :: nynf !< index of north boundary on former subdomain
4295	INTEGER(iwp) :: nysf !< index of south boundary on former subdomain
4296
4297	INTEGER(iwp) :: nxl_on_file !< leftmost index on file
4298	INTEGER(iwp) :: nxr_on_file !< rightmost index on file
4299	INTEGER(iwp) :: nyn_on_file !< northmost index on file
4300	INTEGER(iwp) :: nys_on_file !< southmost index on file
4301
4302	INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) :: start_index_c
4303	INTEGER(iwp), DIMENSION(nys_on_file:nyn_on_file,nxl_on_file:nxr_on_file) :: &
4304	start_index_on_file !< start index of surface elements on file
4305	INTEGER(iwp), DIMENSION(nys_on_file:nyn_on_file,nxl_on_file:nxr_on_file) :: &
4306	end_index_on_file !< end index of surface elements on file
4307
4308	REAL(wp), DIMENSION(:) :: surf_target !< target surface type
4309	REAL(wp), DIMENSION(:) :: surf_file !< surface type on file
4310
4311	ic = nxlc
4312	DO i = nxlf, nxrf
4313	jc = nysc
4314	DO j = nysf, nynf
4315
4316	mm = start_index_c(jc,ic)
4317	DO m = start_index_on_file(j,i), end_index_on_file(j,i)
4318	surf_target(mm) = surf_file(m)
4319	mm = mm + 1
4320	ENDDO
4321
4322	jc = jc + 1
4323	ENDDO
4324	ic = ic + 1
4325	ENDDO
4326
4327
4328	END SUBROUTINE surface_restore_elements_1d
4329
4330	!------------------------------------------------------------------------------!
4331	! Description:
4332	! ------------
4333	!> Routine maps surface data read from file after restart - 2D arrays
4334	!------------------------------------------------------------------------------!
4335	SUBROUTINE surface_restore_elements_2d( surf_target, surf_file, &
4336	start_index_c, &
4337	start_index_on_file, &
4338	end_index_on_file, &
4339	nxlc, nysc, nxlf, nxrf, nysf, nynf,&
4340	nys_on_file, nyn_on_file, &
4341	nxl_on_file,nxr_on_file )
4342
4343	IMPLICIT NONE
4344
4345	INTEGER(iwp) :: i !< running index along x-direction, refers to former domain size
4346	INTEGER(iwp) :: ic !< running index along x-direction, refers to current domain size
4347	INTEGER(iwp) :: j !< running index along y-direction, refers to former domain size
4348	INTEGER(iwp) :: jc !< running index along y-direction, refers to former domain size
4349	INTEGER(iwp) :: m !< surface-element index on file
4350	INTEGER(iwp) :: mm !< surface-element index on current subdomain
4351	INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain
4352	INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain
4353	INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain
4354	INTEGER(iwp) :: nysc !< index of north boundary on current subdomain
4355	INTEGER(iwp) :: nynf !< index of north boundary on former subdomain
4356	INTEGER(iwp) :: nysf !< index of south boundary on former subdomain
4357
4358	INTEGER(iwp) :: nxl_on_file !< leftmost index on file
4359	INTEGER(iwp) :: nxr_on_file !< rightmost index on file
4360	INTEGER(iwp) :: nyn_on_file !< northmost index on file
4361	INTEGER(iwp) :: nys_on_file !< southmost index on file
4362
4363	INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) :: start_index_c !< start index of surface type
4364	INTEGER(iwp), DIMENSION(nys_on_file:nyn_on_file,nxl_on_file:nxr_on_file) :: &
4365	start_index_on_file !< start index of surface elements on file
4366	INTEGER(iwp), DIMENSION(nys_on_file:nyn_on_file,nxl_on_file:nxr_on_file) :: &
4367	end_index_on_file !< end index of surface elements on file
4368
4369	REAL(wp), DIMENSION(:,:) :: surf_target !< target surface type
4370	REAL(wp), DIMENSION(:,:) :: surf_file !< surface type on file
4371
4372	ic = nxlc
4373	DO i = nxlf, nxrf
4374	jc = nysc
4375	DO j = nysf, nynf
4376	mm = start_index_c(jc,ic)
4377	DO m = start_index_on_file(j,i), end_index_on_file(j,i)
4378	surf_target(:,mm) = surf_file(:,m)
4379	mm = mm + 1
4380	ENDDO
4381
4382	jc = jc + 1
4383	ENDDO
4384	ic = ic + 1
4385	ENDDO
4386
4387	END SUBROUTINE surface_restore_elements_2d
4388
4389
4390	END MODULE surface_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |