Home

Context Navigation

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

Last change on this file since 2501 was 2478, checked in by suehring, 7 years ago
Bugfixes concerning top fluxes and TKE production
Property svn:keywords set to `Id`
File size: 161.8 KB

Rev	Line
[2232]	1	!> @file surface_mod.f90
	2	!------------------------------------------------------------------------------!
	3	! This file is part of PALM.
	4	!
	5	! PALM is free software: you can redistribute it and/or modify it under the
	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	!
[2233]	17	! Copyright 1997-2017 Leibniz Universitaet Hannover
[2232]	18	!
	19	!------------------------------------------------------------------------------!
	20	! Current revisions:
	21	! ------------------
	22	!
[2318]	23	!
[2232]	24	! Former revisions:
	25	! -----------------
[2256]	26	! $Id: surface_mod.f90 2478 2017-09-18 13:37:24Z raasch $
[2478]	27	! Bugfixes in initializing model top
	28	!
	29	! 2378 2017-08-31 13:57:27Z suehring
[2378]	30	! Bugfix in write restart data
	31	!
	32	! 2339 2017-08-07 13:55:26Z gronemeier
[2339]	33	! corrected timestamp in header
	34	!
	35	! 2338 2017-08-07 12:15:38Z gronemeier
[2338]	36	! Modularize 1D model
	37	!
[2339]	38	! 2318 2017-07-20 17:27:44Z suehring
[2318]	39	! New function to obtain topography top index.
	40	!
	41	! 2317 2017-07-20 17:27:19Z suehring
[2292]	42	! Implementation of new microphysic scheme: cloud_scheme = 'morrison'
	43	! includes two more prognostic equations for cloud drop concentration (nc)
	44	! and cloud water content (qc).
	45	!
	46	! 2270 2017-06-09 12:18:47Z maronga
[2270]	47	! Parameters removed/added due to changes in the LSM
	48	!
	49	! 2269 2017-06-09 11:57:32Z suehring
[2269]	50	! Formatting and description adjustments
	51	!
	52	! 2256 2017-06-07 13:58:08Z suehring
[2256]	53	! Enable heating at downward-facing surfaces
	54	!
	55	! 2233 2017-05-30 18:08:54Z suehring
[2232]	56	! Initial revision
	57	!
	58	!
	59	! Description:
	60	! ------------
[2269]	61	!> Surface module defines derived data structures to treat surface-
[2232]	62	!> bounded grid cells. Three different types of surfaces are defined:
[2269]	63	!> default surfaces, natural surfaces, and urban surfaces. The module
	64	!> encompasses the allocation and initialization of surface arrays, and handles
	65	!> reading and writing restart data.
[2232]	66	!> In addition, a further derived data structure is defined, in order to set
[2269]	67	!> boundary conditions at surfaces.
[2232]	68	!------------------------------------------------------------------------------!
	69	MODULE surface_mod
	70
	71	USE arrays_3d, &
	72	ONLY: zu, zw, heatflux_input_conversion, waterflux_input_conversion, &
	73	momentumflux_input_conversion
	74
	75	USE control_parameters
	76
	77	USE indices, &
	78	ONLY: nxl, nxlg, nxr, nxrg, nys, nysg, nyn, nyng, nzb, nzt, wall_flags_0
	79
	80	USE grid_variables, &
	81	ONLY: dx, dy
	82
	83	USE kinds
	84
[2338]	85	USE model_1d_mod, &
[2232]	86	ONLY: rif1d, us1d, usws1d, vsws1d
	87
[2338]	88
[2232]	89	IMPLICIT NONE
	90
	91	!
	92	!-- Data type used to identify grid-points where horizontal boundary conditions
	93	!-- are applied
	94	TYPE bc_type
	95
	96	INTEGER(iwp) :: ns !< number of surface elements on the PE
	97
	98	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: i !< x-index linking to the PALM 3D-grid
	99	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: j !< y-index linking to the PALM 3D-grid
	100	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: k !< z-index linking to the PALM 3D-grid
	101
	102	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: start_index !< start index within surface data type for given (j,i)
	103	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: end_index !< end index within surface data type for given (j,i)
	104
	105	END TYPE bc_type
	106	!
	107	!-- Data type used to identify and treat surface-bounded grid points
	108	TYPE surf_type
	109
	110	INTEGER(iwp) :: ns !< number of surface elements on the PE
	111	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: i !< x-index linking to the PALM 3D-grid
	112	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: j !< y-index linking to the PALM 3D-grid
	113	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: k !< z-index linking to the PALM 3D-grid
	114
	115	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: facing !< Bit indicating surface orientation
	116
	117	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: start_index !< Start index within surface data type for given (j,i)
	118	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: end_index !< End index within surface data type for given (j,i)
	119
	120	REAL(wp), DIMENSION(:), ALLOCATABLE :: z_mo !< surface-layer height
	121	REAL(wp), DIMENSION(:), ALLOCATABLE :: uvw_abs !< absolute surface-parallel velocity
	122	REAL(wp), DIMENSION(:), ALLOCATABLE :: us !< friction velocity
	123	REAL(wp), DIMENSION(:), ALLOCATABLE :: ts !< scaling parameter temerature
	124	REAL(wp), DIMENSION(:), ALLOCATABLE :: qs !< scaling parameter humidity
	125	REAL(wp), DIMENSION(:), ALLOCATABLE :: ss !< scaling parameter passive scalar
[2292]	126	REAL(wp), DIMENSION(:), ALLOCATABLE :: qcs !< scaling parameter qc
	127	REAL(wp), DIMENSION(:), ALLOCATABLE :: ncs !< scaling parameter nc
[2232]	128	REAL(wp), DIMENSION(:), ALLOCATABLE :: qrs !< scaling parameter qr
	129	REAL(wp), DIMENSION(:), ALLOCATABLE :: nrs !< scaling parameter nr
	130
	131	REAL(wp), DIMENSION(:), ALLOCATABLE :: ol !< Obukhov length
	132	REAL(wp), DIMENSION(:), ALLOCATABLE :: rib !< Richardson bulk number
	133
	134	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0 !< roughness length for momentum
	135	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0h !< roughness length for heat
	136	REAL(wp), DIMENSION(:), ALLOCATABLE :: z0q !< roughness length for humidity
	137
	138	REAL(wp), DIMENSION(:), ALLOCATABLE :: pt1 !< Specific humidity at first grid level (required for cloud_physics = .T.)
	139	REAL(wp), DIMENSION(:), ALLOCATABLE :: qv1 !< Potential temperature at first grid level (required for cloud_physics = .T.)
	140	!
	141	!-- Define arrays for surface fluxes
	142	REAL(wp), DIMENSION(:), ALLOCATABLE :: usws !< vertical momentum flux for u-component at horizontal surfaces
	143	REAL(wp), DIMENSION(:), ALLOCATABLE :: vsws !< vertical momentum flux for v-component at horizontal surfaces
	144
	145	REAL(wp), DIMENSION(:), ALLOCATABLE :: shf !< surface flux sensible heat
	146	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws !< surface flux latent heat
	147	REAL(wp), DIMENSION(:), ALLOCATABLE :: ssws !< surface flux passive scalar
[2292]	148	REAL(wp), DIMENSION(:), ALLOCATABLE :: qcsws !< surface flux qc
	149	REAL(wp), DIMENSION(:), ALLOCATABLE :: ncsws !< surface flux nc
[2232]	150	REAL(wp), DIMENSION(:), ALLOCATABLE :: qrsws !< surface flux qr
	151	REAL(wp), DIMENSION(:), ALLOCATABLE :: nrsws !< surface flux nr
	152	REAL(wp), DIMENSION(:), ALLOCATABLE :: sasws !< surface flux salinity
	153	!
	154	!-- Required for horizontal walls in production_e
	155	REAL(wp), DIMENSION(:), ALLOCATABLE :: u_0 !< virtual velocity component (see production_e_init for further explanation)
	156	REAL(wp), DIMENSION(:), ALLOCATABLE :: v_0 !< virtual velocity component (see production_e_init for further explanation)
	157
	158	REAL(wp), DIMENSION(:), ALLOCATABLE :: mom_flux_uv !< momentum flux usvs and vsus at vertical surfaces (used in diffusion_u and diffusion_v)
	159	REAL(wp), DIMENSION(:), ALLOCATABLE :: mom_flux_w !< momentum flux wsus and wsvs at vertical surfaces (used in diffusion_w)
	160	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: mom_flux_tke !< momentum flux usvs, vsus, wsus, wsvs at vertical surfaces at grid center (used in production_e)
	161	!
	162	!-- Variables required for LSM as well as for USM
	163	REAL(wp), DIMENSION(:), ALLOCATABLE :: pt_surface !< skin-surface temperature
	164	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net_l !< net radiation
	165	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_h !< heat conductivity of soil (W/m/K)
[2270]	166	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_h_def !< default heat conductivity of soil (W/m/K)
	167
[2232]	168	LOGICAL, DIMENSION(:), ALLOCATABLE :: building_surface !< flag parameter indicating that the surface element is covered by buildings (no LSM actions, not implemented yet)
[2270]	169	LOGICAL, DIMENSION(:), ALLOCATABLE :: pavement_surface !< flag parameter for pavements
	170	LOGICAL, DIMENSION(:), ALLOCATABLE :: water_surface !< flag parameter for water surfaces
	171	LOGICAL, DIMENSION(:), ALLOCATABLE :: vegetation_surface !< flag parameter for natural land surfaces
[2232]	172
	173	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_liq !< liquid water coverage (of vegetated area)
	174	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_veg !< vegetation coverage
	175	REAL(wp), DIMENSION(:), ALLOCATABLE :: f_sw_in !< fraction of absorbed shortwave radiation by the surface layer (not implemented yet)
[2270]	176	REAL(wp), DIMENSION(:), ALLOCATABLE :: ghf !< ground heat flux
[2232]	177	REAL(wp), DIMENSION(:), ALLOCATABLE :: g_d !< coefficient for dependence of r_canopy on water vapour pressure deficit
	178	REAL(wp), DIMENSION(:), ALLOCATABLE :: lai !< leaf area index
[2270]	179	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surface_u !< coupling between surface and soil (depends on vegetation type) (W/m2/K)
	180	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surface_s !< coupling between surface and soil (depends on vegetation type) (W/m2/K)
	181	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_liq !< surface flux of latent heat (liquid water portion)
	182	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_soil !< surface flux of latent heat (soil portion)
	183	REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_veg !< surface flux of latent heat (vegetation portion)
[2232]	184	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_a !< aerodynamic resistance
	185	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_canopy !< canopy resistance
	186	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_soil !< soil resistance
	187	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_soil_min !< minimum soil resistance
	188	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_s !< total surface resistance (combination of r_soil and r_canopy)
	189	REAL(wp), DIMENSION(:), ALLOCATABLE :: r_canopy_min !< minimum canopy (stomatal) resistance
	190
[2270]	191	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: alpha_vg !< coef. of Van Genuchten
[2232]	192	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: lambda_w !< hydraulic diffusivity of soil (?)
	193	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: gamma_w !< hydraulic conductivity of soil (W/m/K)
[2270]	194	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: gamma_w_sat !< hydraulic conductivity at saturation
	195	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: l_vg !< coef. of Van Genuchten
	196	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_fc !< soil moisture at field capacity (m3/m3)
	197	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_res !< residual soil moisture
	198	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_sat !< saturation soil moisture (m3/m3)
	199	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: m_wilt !< soil moisture at permanent wilting point (m3/m3)
	200	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: n_vg !< coef. Van Genuchten
	201	REAL(wp), DIMENSION(:), ALLOCATABLE :: rho_c_def !< default volumetric heat capacity of the (soil) layer (J/m3/K)
	202	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_total !< volumetric heat capacity of the actual soil matrix (J/m3/K)
[2232]	203	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: root_fr !< root fraction within the soil layers
	204	!
	205	!-- Urban surface variables
	206	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: surface_types !< array of types of wall parameters
	207
	208	LOGICAL, DIMENSION(:), ALLOCATABLE :: isroof_surf !< flag indication roof surfaces
	209
	210	REAL(wp), DIMENSION(:), ALLOCATABLE :: albedo_surf !< albedo of the surface
[2270]	211	REAL(wp), DIMENSION(:), ALLOCATABLE :: c_surface !< heat capacity of the wall surface skin (J/m2/K)
[2232]	212	REAL(wp), DIMENSION(:), ALLOCATABLE :: emiss_surf !< emissivity of the wall surface
[2270]	213	REAL(wp), DIMENSION(:), ALLOCATABLE :: lambda_surf !< heat conductivity between air and surface (W/m2/K)
[2232]	214	REAL(wp), DIMENSION(:), ALLOCATABLE :: roughness_wall !< roughness relative to concrete
	215	REAL(wp), DIMENSION(:), ALLOCATABLE :: thickness_wall !< thickness of the wall, roof and soil layers
	216
	217	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsl !< reflected shortwave radiation for local surface in i-th reflection
	218	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutll !< reflected + emitted longwave radiation for local surface in i-th reflection
	219	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf !< total radiation flux incoming to minus outgoing from local surface
	220
	221	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_m !< surface temperature tendency (K)
	222	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf !< kinematic wall heat flux of sensible heat (actually no longer needed)
	223	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf_eb !< wall heat flux of sensible heat in wall normal direction
	224
	225
	226	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb !< wall ground heat flux
	227
	228	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_in_sw !< incoming shortwave radiation
	229	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_out_sw !< emitted shortwave radiation
	230	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_in_lw !< incoming longwave radiation
	231	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_out_lw !< emitted longwave radiation
	232
	233	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinsw !< shortwave radiation falling to local surface including radiation from reflections
	234	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsw !< total shortwave radiation outgoing from nonvirtual surfaces surfaces after all reflection
	235	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlw !< longwave radiation falling to local surface including radiation from reflections
	236	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw !< total longwave radiation outgoing from nonvirtual surfaces surfaces after all reflection
	237
	238
	239	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: rho_c_wall !< volumetric heat capacity of the material ( J m-3 K-1 ) (= 2.19E6)
	240	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_wall !< wall grid spacing (center-center)
	241	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_wall !< 1/dz_wall
	242	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: dz_wall_stag !< wall grid spacing (edge-edge)
	243	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ddz_wall_stag !< 1/dz_wall_stag
	244	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_wall_m !< t_wall prognostic array
	245	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zw !< wall layer depths (m)
	246
	247
	248	!-- arrays for time averages
	249	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net_av !< average of rad_net_l
	250	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinsw_av !< average of sw radiation falling to local surface including radiation from reflections
	251	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlw_av !< average of lw radiation falling to local surface including radiation from reflections
	252	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdir_av !< average of direct sw radiation falling to local surface
	253	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdif_av !< average of diffuse sw radiation from sky and model boundary falling to local surface
	254	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwdif_av !< average of diffuse lw radiation from sky and model boundary falling to local surface
	255	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswref_av !< average of sw radiation falling to surface from reflections
	256	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwref_av !< average of lw radiation falling to surface from reflections
	257	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsw_av !< average of total sw radiation outgoing from nonvirtual surfaces surfaces after all reflection
	258	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw_av !< average of total lw radiation outgoing from nonvirtual surfaces surfaces after all reflection
	259	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfins_av !< average of array of residua of sw radiation absorbed in surface after last reflection
	260	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinl_av !< average of array of residua of lw radiation absorbed in surface after last reflection
	261	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf_av !< average of total radiation flux incoming to minus outgoing from local surface
	262	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_av !< average of wghf_eb
	263	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf_eb_av !< average of wshf_eb
	264	REAL(wp), DIMENSION(:), ALLOCATABLE :: t_surf_av !< average of wall surface temperature (K)
	265
	266	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_wall_av !< Average of t_wall
	267
	268	END TYPE surf_type
	269
	270	TYPE (bc_type), DIMENSION(0:1) :: bc_h !< boundary condition data type, horizontal upward- and downward facing surfaces
	271
	272	TYPE (surf_type), DIMENSION(0:2), TARGET :: surf_def_h !< horizontal default surfaces (Up, Down, and Top)
	273	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_def_v !< vertical default surfaces (North, South, West, East)
	274	TYPE (surf_type) , TARGET :: surf_lsm_h !< horizontal natural land surfaces, so far only upward-facing
	275	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_lsm_v !< vertical land surfaces (North, South, West, East)
	276	TYPE (surf_type) , TARGET :: surf_usm_h !< horizontal urban surfaces, so far only upward-facing
	277	TYPE (surf_type), DIMENSION(0:3), TARGET :: surf_usm_v !< vertical urban surfaces (North, South, West, East)
	278
	279	INTEGER(iwp) :: ns_h_on_file(0:2) !< total number of horizontal surfaces with the same facing, required for writing restart data
	280	INTEGER(iwp) :: ns_v_on_file(0:3) !< total number of vertical surfaces with the same facing, required for writing restart data
	281
	282
	283	SAVE
	284
	285	PRIVATE
	286
[2317]	287	INTERFACE get_topography_top_index
	288	MODULE PROCEDURE get_topography_top_index
	289	END INTERFACE get_topography_top_index
	290
[2232]	291	INTERFACE init_bc
	292	MODULE PROCEDURE init_bc
	293	END INTERFACE init_bc
	294
	295	INTERFACE init_surfaces
	296	MODULE PROCEDURE init_surfaces
	297	END INTERFACE init_surfaces
	298
	299	INTERFACE init_surface_arrays
	300	MODULE PROCEDURE init_surface_arrays
	301	END INTERFACE init_surface_arrays
	302
	303	INTERFACE surface_read_restart_data
	304	MODULE PROCEDURE surface_read_restart_data
	305	END INTERFACE surface_read_restart_data
	306
	307	INTERFACE surface_write_restart_data
	308	MODULE PROCEDURE surface_write_restart_data
	309	END INTERFACE surface_write_restart_data
	310
	311	INTERFACE surface_last_actions
	312	MODULE PROCEDURE surface_last_actions
	313	END INTERFACE surface_last_actions
	314
[2317]	315	!
	316	!-- Public variables
	317	PUBLIC bc_h, ns_h_on_file, ns_v_on_file, surf_def_h, surf_def_v, &
	318	surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v, surf_type
	319	!
	320	!-- Public subroutines and functions
	321	PUBLIC get_topography_top_index, init_bc, init_surfaces, &
	322	init_surface_arrays, surface_read_restart_data, &
	323	surface_write_restart_data, surface_last_actions
[2232]	324
[2317]	325
[2232]	326	CONTAINS
	327
	328	!------------------------------------------------------------------------------!
	329	! Description:
	330	! ------------
	331	!> Initialize data type for setting boundary conditions at horizontal surfaces.
	332	!------------------------------------------------------------------------------!
	333	SUBROUTINE init_bc
	334
	335	IMPLICIT NONE
	336
	337	INTEGER(iwp) :: i !<
	338	INTEGER(iwp) :: j !<
	339	INTEGER(iwp) :: k !<
	340
	341	INTEGER(iwp), DIMENSION(0:1) :: num_h !<
	342	INTEGER(iwp), DIMENSION(0:1) :: num_h_kji !<
	343	INTEGER(iwp), DIMENSION(0:1) :: start_index_h !<
	344
	345	!
	346	!-- First of all, count the number of upward- and downward-facing surfaces
	347	num_h = 0
	348	DO i = nxlg, nxrg
	349	DO j = nysg, nyng
	350	DO k = nzb+1, nzt
	351	!
	352	!-- Check if current gridpoint belongs to the atmosphere
	353	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
	354	!
	355	!-- Upward-facing
	356	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) &
	357	num_h(0) = num_h(0) + 1
	358	!
	359	!-- Downward-facing
	360	IF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) &
	361	num_h(1) = num_h(1) + 1
	362	ENDIF
	363	ENDDO
	364	ENDDO
	365	ENDDO
	366	!
	367	!-- Save the number of surface elements
	368	bc_h(0)%ns = num_h(0)
	369	bc_h(1)%ns = num_h(1)
	370	!
	371	!-- ALLOCATE data type variables
	372	!-- Upward facing
	373	ALLOCATE( bc_h(0)%i(1:bc_h(0)%ns) )
	374	ALLOCATE( bc_h(0)%j(1:bc_h(0)%ns) )
	375	ALLOCATE( bc_h(0)%k(1:bc_h(0)%ns) )
	376	ALLOCATE( bc_h(0)%start_index(nysg:nyng,nxlg:nxrg) )
	377	ALLOCATE( bc_h(0)%end_index(nysg:nyng,nxlg:nxrg) )
	378	bc_h(0)%start_index = 1
	379	bc_h(0)%end_index = 0
	380	!
	381	!-- Downward facing
	382	ALLOCATE( bc_h(1)%i(1:bc_h(1)%ns) )
	383	ALLOCATE( bc_h(1)%j(1:bc_h(1)%ns) )
	384	ALLOCATE( bc_h(1)%k(1:bc_h(1)%ns) )
	385	ALLOCATE( bc_h(1)%start_index(nysg:nyng,nxlg:nxrg) )
	386	ALLOCATE( bc_h(1)%end_index(nysg:nyng,nxlg:nxrg) )
	387	bc_h(1)%start_index = 1
	388	bc_h(1)%end_index = 0
	389	!
	390	!-- Store the respective indices on data type
	391	num_h(0:1) = 1
	392	start_index_h(0:1) = 1
	393	DO i = nxlg, nxrg
	394	DO j = nysg, nyng
	395
	396	num_h_kji(0:1) = 0
	397	DO k = nzb+1, nzt
	398	!
	399	!-- Check if current gridpoint belongs to the atmosphere
	400	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
	401	!
	402	!-- Upward-facing
	403	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
	404	bc_h(0)%i(num_h(0)) = i
	405	bc_h(0)%j(num_h(0)) = j
	406	bc_h(0)%k(num_h(0)) = k
	407	num_h_kji(0) = num_h_kji(0) + 1
	408	num_h(0) = num_h(0) + 1
	409	ENDIF
	410	!
	411	!-- Downward-facing
	412	IF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
	413	bc_h(1)%i(num_h(1)) = i
	414	bc_h(1)%j(num_h(1)) = j
	415	bc_h(1)%k(num_h(1)) = k
	416	num_h_kji(1) = num_h_kji(1) + 1
	417	num_h(1) = num_h(1) + 1
	418	ENDIF
	419	ENDIF
	420	ENDDO
	421	bc_h(0)%start_index(j,i) = start_index_h(0)
	422	bc_h(0)%end_index(j,i) = bc_h(0)%start_index(j,i) + num_h_kji(0) - 1
	423	start_index_h(0) = bc_h(0)%end_index(j,i) + 1
	424
	425	bc_h(1)%start_index(j,i) = start_index_h(1)
	426	bc_h(1)%end_index(j,i) = bc_h(1)%start_index(j,i) + num_h_kji(1) - 1
	427	start_index_h(1) = bc_h(1)%end_index(j,i) + 1
	428	ENDDO
	429	ENDDO
	430
	431
	432	END SUBROUTINE init_bc
	433
	434
	435	!------------------------------------------------------------------------------!
	436	! Description:
	437	! ------------
	438	!> Initialize horizontal and vertical surfaces. Counts the number of default-,
	439	!> natural and urban surfaces and allocates memory, respectively.
	440	!------------------------------------------------------------------------------!
	441	SUBROUTINE init_surface_arrays
	442
	443	IMPLICIT NONE
	444
	445	INTEGER(iwp) :: i !< running index x-direction
	446	INTEGER(iwp) :: j !< running index y-direction
	447	INTEGER(iwp) :: k !< running index z-direction
	448	INTEGER(iwp) :: l !< index variable for surface facing
	449	INTEGER(iwp) :: num_lsm_h !< number of horizontally-aligned natural surfaces
	450	INTEGER(iwp) :: num_usm_h !< number of horizontally-aligned urban surfaces
	451
	452	INTEGER(iwp), DIMENSION(0:2) :: num_def_h !< number of horizontally-aligned default surfaces
	453	INTEGER(iwp), DIMENSION(0:3) :: num_def_v !< number of vertically-aligned default surfaces
	454	INTEGER(iwp), DIMENSION(0:3) :: num_lsm_v !< number of vertically-aligned natural surfaces
	455	INTEGER(iwp), DIMENSION(0:3) :: num_usm_v !< number of vertically-aligned urban surfaces
	456
	457
	458	num_def_h = 0
	459	num_def_v = 0
	460	num_lsm_h = 0
	461	num_lsm_v = 0
	462	num_usm_h = 0
	463	num_usm_v = 0
	464	!
	465	!-- Count number of horizontal surfaces on local domain
	466	DO i = nxl, nxr
	467	DO j = nys, nyn
	468	DO k = nzb+1, nzt
	469	!
	470	!-- Check if current gridpoint belongs to the atmosphere
	471	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
	472	!
	473	!-- Check if grid point adjoins to any upward-facing horizontal
	474	!-- surface, e.g. the Earth surface, plane roofs, or ceilings.
	475	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
	476	!
	477	!-- Land-surface type
	478	IF ( land_surface ) THEN
	479	num_lsm_h = num_lsm_h + 1
	480	!
	481	!-- Urban surface tpye
	482	ELSEIF ( urban_surface ) THEN
	483	num_usm_h = num_usm_h + 1
	484	!
	485	!-- Default-surface type
	486	ELSE
	487	num_def_h(0) = num_def_h(0) + 1
	488	ENDIF
	489
	490	ENDIF
	491	!
	492	!-- Check for top-fluxes
	493	IF ( k == nzt .AND. use_top_fluxes ) THEN
	494	num_def_h(2) = num_def_h(2) + 1
	495	!
	496	!-- Check for any other downward-facing surface. So far only for
	497	!-- default surface type.
	498	ELSEIF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
	499	num_def_h(1) = num_def_h(1) + 1
	500	ENDIF
	501
	502	ENDIF
	503	ENDDO
	504	ENDDO
	505	ENDDO
	506	!
	507	!-- Count number of vertical surfaces on local domain
	508	DO i = nxl, nxr
	509	DO j = nys, nyn
	510	DO k = nzb+1, nzt
	511	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
	512	!
	513	!-- Northward-facing
	514	IF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 0 ) ) THEN
	515	IF ( urban_surface ) THEN
	516	num_usm_v(0) = num_usm_v(0) + 1
	517	ELSEIF ( land_surface ) THEN
	518	num_lsm_v(0) = num_lsm_v(0) + 1
	519	ELSE
	520	num_def_v(0) = num_def_v(0) + 1
	521	ENDIF
	522	ENDIF
	523	!
	524	!-- Southward-facing
	525	IF ( .NOT. BTEST( wall_flags_0(k,j+1,i), 0 ) ) THEN
	526	IF ( urban_surface ) THEN
	527	num_usm_v(1) = num_usm_v(1) + 1
	528	ELSEIF ( land_surface ) THEN
	529	num_lsm_v(1) = num_lsm_v(1) + 1
	530	ELSE
	531	num_def_v(1) = num_def_v(1) + 1
	532	ENDIF
	533	ENDIF
	534	!
	535	!-- Eastward-facing
	536	IF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 0 ) ) THEN
	537	IF ( urban_surface ) THEN
	538	num_usm_v(2) = num_usm_v(2) + 1
	539	ELSEIF ( land_surface ) THEN
	540	num_lsm_v(2) = num_lsm_v(2) + 1
	541	ELSE
	542	num_def_v(2) = num_def_v(2) + 1
	543	ENDIF
	544	ENDIF
	545	!
	546	!-- Westward-facing
	547	IF ( .NOT. BTEST( wall_flags_0(k,j,i+1), 0 ) ) THEN
	548	IF ( urban_surface ) THEN
	549	num_usm_v(3) = num_usm_v(3) + 1
	550	ELSEIF ( land_surface ) THEN
	551	num_lsm_v(3) = num_lsm_v(3) + 1
	552	ELSE
	553	num_def_v(3) = num_def_v(3) + 1
	554	ENDIF
	555	ENDIF
	556	ENDIF
	557	ENDDO
	558	ENDDO
	559	ENDDO
	560
	561	!
	562	!-- Store number of surfaces per core.
	563	!-- Horizontal surface, default type, upward facing
	564	surf_def_h(0)%ns = num_def_h(0)
	565	!
	566	!-- Horizontal surface, default type, downward facing
	567	surf_def_h(1)%ns = num_def_h(1)
	568	!
	569	!-- Horizontal surface, default type, top downward facing
	570	surf_def_h(2)%ns = num_def_h(2)
	571	!
	572	!-- Horizontal surface, natural type, so far only upward-facing
	573	surf_lsm_h%ns = num_lsm_h
	574	!
	575	!-- Horizontal surface, urban type, so far only upward-facing
	576	surf_usm_h%ns = num_usm_h
	577	!
	578	!-- Vertical surface, default type, northward facing
	579	surf_def_v(0)%ns = num_def_v(0)
	580	!
	581	!-- Vertical surface, default type, southward facing
	582	surf_def_v(1)%ns = num_def_v(1)
	583	!
	584	!-- Vertical surface, default type, eastward facing
	585	surf_def_v(2)%ns = num_def_v(2)
	586	!
	587	!-- Vertical surface, default type, westward facing
	588	surf_def_v(3)%ns = num_def_v(3)
	589	!
	590	!-- Vertical surface, natural type, northward facing
	591	surf_lsm_v(0)%ns = num_lsm_v(0)
	592	!
	593	!-- Vertical surface, natural type, southward facing
	594	surf_lsm_v(1)%ns = num_lsm_v(1)
	595	!
	596	!-- Vertical surface, natural type, eastward facing
	597	surf_lsm_v(2)%ns = num_lsm_v(2)
	598	!
	599	!-- Vertical surface, natural type, westward facing
	600	surf_lsm_v(3)%ns = num_lsm_v(3)
	601	!
	602	!-- Vertical surface, urban type, northward facing
	603	surf_usm_v(0)%ns = num_usm_v(0)
	604	!
	605	!-- Vertical surface, urban type, southward facing
	606	surf_usm_v(1)%ns = num_usm_v(1)
	607	!
	608	!-- Vertical surface, urban type, eastward facing
	609	surf_usm_v(2)%ns = num_usm_v(2)
	610	!
	611	!-- Vertical surface, urban type, westward facing
	612	surf_usm_v(3)%ns = num_usm_v(3)
	613	!
	614	!-- Allocate required attributes for horizontal surfaces - default type.
	615	!-- Upward-facing (l=0) and downward-facing (l=1).
	616	DO l = 0, 1
	617	CALL allocate_surface_attributes_h ( surf_def_h(l), nys, nyn, nxl, nxr )
	618	ENDDO
	619	!
	620	!-- Allocate required attributes for model top
	621	CALL allocate_surface_attributes_h_top ( surf_def_h(2), nys, nyn, nxl, nxr )
	622	!
	623	!-- Allocate required attributes for horizontal surfaces - natural type.
	624	CALL allocate_surface_attributes_h ( surf_lsm_h, nys, nyn, nxl, nxr )
	625	!
	626	!-- Allocate required attributes for horizontal surfaces - urban type.
	627	CALL allocate_surface_attributes_h ( surf_usm_h, nys, nyn, nxl, nxr )
	628
	629	!
	630	!-- Allocate required attributes for vertical surfaces.
	631	!-- Northward-facing (l=0), southward-facing (l=1), eastward-facing (l=2)
	632	!-- and westward-facing (l=3).
	633	!-- Default type.
	634	DO l = 0, 3
	635	CALL allocate_surface_attributes_v ( surf_def_v(l), .FALSE., &
	636	nys, nyn, nxl, nxr )
	637	ENDDO
	638	!
	639	!-- Natural type
	640	DO l = 0, 3
	641	CALL allocate_surface_attributes_v ( surf_lsm_v(l), .TRUE., &
	642	nys, nyn, nxl, nxr )
	643	ENDDO
	644	!
	645	!-- Urban type
	646	DO l = 0, 3
	647	CALL allocate_surface_attributes_v ( surf_usm_v(l), .FALSE., &
	648	nys, nyn, nxl, nxr )
	649	ENDDO
	650
	651	END SUBROUTINE init_surface_arrays
	652
	653	!------------------------------------------------------------------------------!
	654	! Description:
	655	! ------------
	656	!> Allocating memory for upward and downward-facing horizontal surface types,
	657	!> except for top fluxes.
	658	!------------------------------------------------------------------------------!
	659	SUBROUTINE allocate_surface_attributes_h( surfaces, &
	660	nys_l, nyn_l, nxl_l, nxr_l )
	661
	662	IMPLICIT NONE
	663
	664	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	665	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	666	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	667	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	668
	669	TYPE(surf_type) :: surfaces !< respective surface type
	670
	671	!
	672	!-- Allocate arrays for start and end index of horizontal surface type
	673	!-- for each (j,i)-grid point. This is required e.g. in diffion_x, which is
	674	!-- called for each (j,i). In order to find the location where the
	675	!-- respective flux is store within the surface-type, start- and end-
	676	!-- index are stored for each (j,i). For example, each (j,i) can have
	677	!-- several entries where fluxes for horizontal surfaces might be stored,
	678	!-- e.g. for overhanging structures where several upward-facing surfaces
	679	!-- might exist for given (j,i).
	680	!-- If no surface of respective type exist at current (j,i), set indicies
	681	!-- such that loop in diffusion routines will not be entered.
	682	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
	683	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
	684	surfaces%start_index = 0
	685	surfaces%end_index = -1
	686	!
	687	!-- Indices to locate surface element
	688	ALLOCATE ( surfaces%i(1:surfaces%ns) )
	689	ALLOCATE ( surfaces%j(1:surfaces%ns) )
	690	ALLOCATE ( surfaces%k(1:surfaces%ns) )
	691	!
	692	!-- Surface-layer height
	693	ALLOCATE ( surfaces%z_mo(1:surfaces%ns) )
	694	!
	695	!-- Surface orientation
	696	ALLOCATE ( surfaces%facing(1:surfaces%ns) )
	697	!
	698	!-- Surface-parallel wind velocity
	699	ALLOCATE ( surfaces%uvw_abs(1:surfaces%ns) )
	700	! ALLOCATE ( surfaces%pt_surface(1:surfaces%ns) )
	701	!
	702	!-- Roughness
	703	ALLOCATE ( surfaces%z0(1:surfaces%ns) )
	704	ALLOCATE ( surfaces%z0h(1:surfaces%ns) )
	705	ALLOCATE ( surfaces%z0q(1:surfaces%ns) )
	706	!
	707	!-- Friction velocity
	708	ALLOCATE ( surfaces%us(1:surfaces%ns) )
	709	!
	710	!-- Stability parameter
	711	ALLOCATE ( surfaces%ol(1:surfaces%ns) )
	712	!
	713	!-- Bulk Richardson number
	714	ALLOCATE ( surfaces%rib(1:surfaces%ns) )
	715	!
	716	!-- Vertical momentum fluxes of u and v
	717	ALLOCATE ( surfaces%usws(1:surfaces%ns) )
	718	ALLOCATE ( surfaces%vsws(1:surfaces%ns) )
	719	!
	720	!-- Required in production_e
	721	IF ( .NOT. constant_diffusion ) THEN
	722	ALLOCATE ( surfaces%u_0(1:surfaces%ns) )
	723	ALLOCATE ( surfaces%v_0(1:surfaces%ns) )
	724	ENDIF
	725	!
	726	!-- Characteristic temperature and surface flux of sensible heat
	727	ALLOCATE ( surfaces%ts(1:surfaces%ns) )
	728	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
	729	!
	730	!-- Characteristic humidity and surface flux of latent heat
	731	IF ( humidity ) THEN
	732	ALLOCATE ( surfaces%qs(1:surfaces%ns) )
	733	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
	734	ENDIF
	735	!
	736	!-- Characteristic scalar and surface flux of scalar
	737	IF ( passive_scalar ) THEN
	738	ALLOCATE ( surfaces%ss(1:surfaces%ns) )
	739	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
	740	ENDIF
	741	!
	742	!-- When cloud physics is used, arrays for storing potential temperature and
	743	!-- specific humidity at first grid level are required.
	744	IF ( cloud_physics ) THEN
	745	ALLOCATE ( surfaces%pt1(1:surfaces%ns) )
	746	ALLOCATE ( surfaces%qv1(1:surfaces%ns) )
	747	ENDIF
	748	!
	749	!--
[2292]	750	IF ( cloud_physics .AND. microphysics_morrison) THEN
	751	ALLOCATE ( surfaces%qcs(1:surfaces%ns) )
	752	ALLOCATE ( surfaces%ncs(1:surfaces%ns) )
	753	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
	754	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
	755	ENDIF
	756	!
	757	!--
[2232]	758	IF ( cloud_physics .AND. microphysics_seifert) THEN
	759	ALLOCATE ( surfaces%qrs(1:surfaces%ns) )
	760	ALLOCATE ( surfaces%nrs(1:surfaces%ns) )
	761	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
	762	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
	763	ENDIF
	764	!
	765	!-- Salinity surface flux
	766	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
	767
	768	END SUBROUTINE allocate_surface_attributes_h
	769
	770
	771	!------------------------------------------------------------------------------!
	772	! Description:
	773	! ------------
	774	!> Allocating memory for model-top fluxes
	775	!------------------------------------------------------------------------------!
	776	SUBROUTINE allocate_surface_attributes_h_top( surfaces, &
	777	nys_l, nyn_l, nxl_l, nxr_l )
	778
	779	IMPLICIT NONE
	780
	781	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	782	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	783	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	784	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	785
	786	TYPE(surf_type) :: surfaces !< respective surface type
	787
	788	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
	789	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
	790	surfaces%start_index = 0
	791	surfaces%end_index = -1
	792	!
	793	!-- Indices to locate surface (model-top) element
	794	ALLOCATE ( surfaces%i(1:surfaces%ns) )
	795	ALLOCATE ( surfaces%j(1:surfaces%ns) )
	796	ALLOCATE ( surfaces%k(1:surfaces%ns) )
	797	!
	798	!-- Vertical momentum fluxes of u and v
	799	ALLOCATE ( surfaces%usws(1:surfaces%ns) )
	800	ALLOCATE ( surfaces%vsws(1:surfaces%ns) )
	801	!
	802	!-- Sensible heat flux
	803	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
	804	!
	805	!-- Latent heat flux
	806	IF ( humidity ) THEN
	807	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
	808	ENDIF
	809	!
	810	!-- Scalar flux
	811	IF ( passive_scalar ) THEN
	812	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
	813	ENDIF
	814	!
	815	!--
[2292]	816	IF ( cloud_physics .AND. microphysics_morrison) THEN
	817	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
	818	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
	819	ENDIF
	820	!
	821	!--
[2232]	822	IF ( cloud_physics .AND. microphysics_seifert) THEN
	823	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
	824	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
	825	ENDIF
	826	!
	827	!-- Salinity flux
	828	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
	829
	830	END SUBROUTINE allocate_surface_attributes_h_top
	831
	832	!------------------------------------------------------------------------------!
	833	! Description:
	834	! ------------
	835	!> Allocating memory for vertical surface types.
	836	!------------------------------------------------------------------------------!
	837	SUBROUTINE allocate_surface_attributes_v( surfaces, lsm, &
	838	nys_l, nyn_l, nxl_l, nxr_l )
	839
	840	IMPLICIT NONE
	841
	842	INTEGER(iwp) :: nyn_l !< north bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	843	INTEGER(iwp) :: nys_l !< south bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	844	INTEGER(iwp) :: nxl_l !< west bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	845	INTEGER(iwp) :: nxr_l !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
	846
	847	LOGICAL :: lsm !< flag indicating data type of natural land surface
	848
	849	TYPE(surf_type) :: surfaces !< respective surface type
	850
	851	!
	852	!-- Allocate arrays for start and end index of vertical surface type
	853	!-- for each (j,i)-grid point. This is required in diffion_x, which is
	854	!-- called for each (j,i). In order to find the location where the
	855	!-- respective flux is store within the surface-type, start- and end-
	856	!-- index are stored for each (j,i). For example, each (j,i) can have
	857	!-- several entries where fluxes for vertical surfaces might be stored.
	858	!-- In the flat case, where no vertical walls exit, set indicies such
	859	!-- that loop in diffusion routines will not be entered.
	860	ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
	861	ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
	862	surfaces%start_index = 0
	863	surfaces%end_index = -1
	864	!
	865	!-- Indices to locate surface element.
	866	ALLOCATE ( surfaces%i(1:surfaces%ns) )
	867	ALLOCATE ( surfaces%j(1:surfaces%ns) )
	868	ALLOCATE ( surfaces%k(1:surfaces%ns) )
	869	!
	870	!-- Surface-layer height
	871	ALLOCATE ( surfaces%z_mo(1:surfaces%ns) )
	872	!
	873	!-- Surface orientation
	874	ALLOCATE ( surfaces%facing(1:surfaces%ns) )
	875	!
	876	!-- Surface parallel wind velocity
	877	ALLOCATE ( surfaces%uvw_abs(1:surfaces%ns) )
	878	!
	879	!-- Roughness
	880	ALLOCATE ( surfaces%z0(1:surfaces%ns) )
	881	ALLOCATE ( surfaces%z0h(1:surfaces%ns) )
	882	ALLOCATE ( surfaces%z0q(1:surfaces%ns) )
	883
	884	!
	885	!-- Friction velocity
	886	ALLOCATE ( surfaces%us(1:surfaces%ns) )
	887	!
	888	!-- Allocate Obukhov length and bulk Richardson number. Only required
	889	!-- for natural land surfaces
	890	IF ( lsm ) THEN
	891	ALLOCATE( surfaces%ol(1:surfaces%ns) )
	892	ALLOCATE( surfaces%rib(1:surfaces%ns) )
	893	ENDIF
	894	!
	895	!-- Allocate arrays for surface momentum fluxes for u and v. For u at north-
	896	!-- and south-facing surfaces, for v at east- and west-facing surfaces.
	897	ALLOCATE ( surfaces%mom_flux_uv(1:surfaces%ns) )
	898	!
	899	!-- Allocate array for surface momentum flux for w - wsus and wsvs
	900	ALLOCATE ( surfaces%mom_flux_w(1:surfaces%ns) )
	901	!
	902	!-- Allocate array for surface momentum flux for subgrid-scale tke wsus and
	903	!-- wsvs; first index usvs or vsws, second index for wsus or wsvs, depending
	904	!-- on surface.
	905	ALLOCATE ( surfaces%mom_flux_tke(0:1,1:surfaces%ns) )
	906	!
	907	!-- Characteristic temperature and surface flux of sensible heat
	908	ALLOCATE ( surfaces%ts(1:surfaces%ns) )
	909	ALLOCATE ( surfaces%shf(1:surfaces%ns) )
	910	!
	911	!-- Characteristic humidity and surface flux of latent heat
	912	IF ( humidity ) THEN
	913	ALLOCATE ( surfaces%qs(1:surfaces%ns) )
	914	ALLOCATE ( surfaces%qsws(1:surfaces%ns) )
	915	ENDIF
	916	!
	917	!-- Characteristic scalar and surface flux of scalar
	918	IF ( passive_scalar ) THEN
	919	ALLOCATE ( surfaces%ss(1:surfaces%ns) )
	920	ALLOCATE ( surfaces%ssws(1:surfaces%ns) )
	921	ENDIF
	922
	923	IF ( cloud_physics .AND. microphysics_seifert) THEN
[2292]	924	ALLOCATE ( surfaces%qcs(1:surfaces%ns) )
	925	ALLOCATE ( surfaces%ncs(1:surfaces%ns) )
	926	ALLOCATE ( surfaces%qcsws(1:surfaces%ns) )
	927	ALLOCATE ( surfaces%ncsws(1:surfaces%ns) )
	928	ENDIF
	929
	930	IF ( cloud_physics .AND. microphysics_seifert) THEN
[2232]	931	ALLOCATE ( surfaces%qrs(1:surfaces%ns) )
	932	ALLOCATE ( surfaces%nrs(1:surfaces%ns) )
	933	ALLOCATE ( surfaces%qrsws(1:surfaces%ns) )
	934	ALLOCATE ( surfaces%nrsws(1:surfaces%ns) )
	935	ENDIF
	936	!
	937	!-- Salinity surface flux
	938	IF ( ocean ) ALLOCATE ( surfaces%sasws(1:surfaces%ns) )
	939
	940	END SUBROUTINE allocate_surface_attributes_v
	941
	942	!------------------------------------------------------------------------------!
	943	! Description:
	944	! ------------
	945	!> Initialize surface elements.
	946	!------------------------------------------------------------------------------!
	947	SUBROUTINE init_surfaces
	948
	949	IMPLICIT NONE
	950
	951	INTEGER(iwp) :: i !< running index x-direction
	952	INTEGER(iwp) :: j !< running index y-direction
	953	INTEGER(iwp) :: k !< running index z-direction
	954	INTEGER(iwp) :: l !< index variable used to distinguish surface facing
	955	INTEGER(iwp) :: m !< running index surface elements
	956
	957	INTEGER(iwp) :: start_index_lsm_h !< dummy to determing local start index in surface type for given (j,i), for horizontal natural surfaces
	958	INTEGER(iwp) :: start_index_usm_h !< dummy to determing local start index in surface type for given (j,i), for horizontal urban surfaces
	959
	960	INTEGER(iwp) :: num_lsm_h !< current number of horizontal surface element, natural type
	961	INTEGER(iwp) :: num_lsm_h_kji !< dummy to determing local end index in surface type for given (j,i), for for horizonal natural surfaces
	962	INTEGER(iwp) :: num_usm_h !< current number of horizontal surface element, urban type
	963	INTEGER(iwp) :: num_usm_h_kji !< dummy to determing local end index in surface type for given (j,i), for for horizonal urban surfaces
	964
	965	INTEGER(iwp), DIMENSION(0:2) :: num_def_h !< current number of horizontal surface element, default type
	966	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
	967	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
	968
	969	INTEGER(iwp), DIMENSION(0:3) :: num_def_v !< current number of vertical surface element, default type
	970	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
	971	INTEGER(iwp), DIMENSION(0:3) :: num_lsm_v !< current number of vertical surface element, natural type
	972	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
	973	INTEGER(iwp), DIMENSION(0:3) :: num_usm_v !< current number of vertical surface element, urban type
	974	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
	975
	976	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
	977	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
	978	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
	979
	980
	981	!
	982	!-- Initialize surface attributes, store indicies, surfaces orientation, etc.,
	983	num_def_h(0:2) = 1
	984	num_def_v(0:3) = 1
	985
	986	num_lsm_h = 1
	987	num_lsm_v(0:3) = 1
	988
	989	num_usm_h = 1
	990	num_usm_v(0:3) = 1
	991
	992	start_index_def_h(0:2) = 1
	993	start_index_def_v(0:3) = 1
	994
	995	start_index_lsm_h = 1
	996	start_index_lsm_v(0:3) = 1
	997
	998	start_index_usm_h = 1
	999	start_index_usm_v(0:3) = 1
	1000
	1001	DO i = nxl, nxr
	1002	DO j = nys, nyn
	1003
	1004	num_def_h_kji = 0
	1005	num_def_v_kji = 0
	1006	num_lsm_h_kji = 0
	1007	num_lsm_v_kji = 0
	1008	num_usm_h_kji = 0
	1009	num_usm_v_kji = 0
	1010
	1011	DO k = nzb+1, nzt
	1012	!
	1013	!-- Check if current gridpoint belongs to the atmosphere
	1014	IF ( BTEST( wall_flags_0(k,j,i), 0 ) ) THEN
	1015	!
	1016	!-- Upward-facing surface. Distinguish between differet surface types.
	1017	!-- To do, think about method to flag natural and non-natural
	1018	!-- surfaces. Only to ask for land_surface or urban surface
	1019	!-- is just a work-around.
	1020	IF ( .NOT. BTEST( wall_flags_0(k-1,j,i), 0 ) ) THEN
	1021	!
	1022	!-- Natural surface type
	1023	IF ( land_surface ) THEN
	1024	CALL initialize_horizontal_surfaces( k, j, i, &
	1025	surf_lsm_h, &
	1026	num_lsm_h, &
	1027	num_lsm_h_kji, &
	1028	.TRUE., .FALSE. )
	1029	!
	1030	!-- Urban surface tpye
	1031	ELSEIF ( urban_surface ) THEN
	1032	CALL initialize_horizontal_surfaces( k, j, i, &
	1033	surf_usm_h, &
	1034	num_usm_h, &
	1035	num_usm_h_kji, &
	1036	.TRUE., .FALSE. )
	1037	!
	1038	!-- Default surface type
	1039	ELSE
	1040	CALL initialize_horizontal_surfaces( k, j, i, &
	1041	surf_def_h(0), &
	1042	num_def_h(0), &
	1043	num_def_h_kji(0),&
	1044	.TRUE., .FALSE. )
	1045	ENDIF
	1046	ENDIF
	1047	!
	1048	!-- downward-facing surface, first, model top
	1049	IF ( k == nzt .AND. use_top_fluxes ) THEN
	1050	CALL initialize_top( k, j, i, surf_def_h(2), &
	1051	num_def_h(2), num_def_h_kji(2) )
	1052	!
	1053	!-- Check for any other downward-facing surface. So far only for
	1054	!-- default surface type.
	1055	ELSEIF ( .NOT. BTEST( wall_flags_0(k+1,j,i), 0 ) ) THEN
	1056	CALL initialize_horizontal_surfaces( k, j, i, &
	1057	surf_def_h(1), &
	1058	num_def_h(1), &
	1059	num_def_h_kji(1), &
	1060	.FALSE., .TRUE. )
	1061	ENDIF
	1062	!
	1063	!-- Check for vertical walls and, if required, initialize it.
	1064	! Start with northward-facing surface.
	1065	IF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 0 ) ) THEN
	1066	IF ( urban_surface ) THEN
	1067	CALL initialize_vertical_surfaces( 0, k, j, i, &
	1068	surf_usm_v(0), &
	1069	num_usm_v(0), &
	1070	num_usm_v_kji(0), &
	1071	.FALSE., .FALSE., &
	1072	.FALSE., .TRUE. )
	1073	ELSEIF ( land_surface ) THEN
	1074	CALL initialize_vertical_surfaces( 0, k, j, i, &
	1075	surf_lsm_v(0), &
	1076	num_lsm_v(0), &
	1077	num_lsm_v_kji(0), &
	1078	.FALSE., .FALSE., &
	1079	.FALSE., .TRUE. )
	1080	ELSE
	1081	CALL initialize_vertical_surfaces( 0, k, j, i, &
	1082	surf_def_v(0), &
	1083	num_def_v(0), &
	1084	num_def_v_kji(0), &
	1085	.FALSE., .FALSE., &
	1086	.FALSE., .TRUE. )
	1087	ENDIF
	1088	ENDIF
	1089	!
	1090	!-- southward-facing surface
	1091	IF ( .NOT. BTEST( wall_flags_0(k,j+1,i), 0 ) ) THEN
	1092	IF ( urban_surface ) THEN
	1093	CALL initialize_vertical_surfaces( 1, k, j, i, &
	1094	surf_usm_v(1), &
	1095	num_usm_v(1), &
	1096	num_usm_v_kji(1), &
	1097	.FALSE., .FALSE., &
	1098	.TRUE., .FALSE. )
	1099	ELSEIF ( land_surface ) THEN
	1100	CALL initialize_vertical_surfaces( 1, k, j, i, &
	1101	surf_lsm_v(1), &
	1102	num_lsm_v(1), &
	1103	num_lsm_v_kji(1), &
	1104	.FALSE., .FALSE., &
	1105	.TRUE., .FALSE. )
	1106	ELSE
	1107	CALL initialize_vertical_surfaces( 1, k, j, i, &
	1108	surf_def_v(1), &
	1109	num_def_v(1), &
	1110	num_def_v_kji(1), &
	1111	.FALSE., .FALSE., &
	1112	.TRUE., .FALSE. )
	1113	ENDIF
	1114	ENDIF
	1115	!
	1116	!-- eastward-facing surface
	1117	IF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 0 ) ) THEN
	1118	IF ( urban_surface ) THEN
	1119	CALL initialize_vertical_surfaces( 2, k, j, i, &
	1120	surf_usm_v(2), &
	1121	num_usm_v(2), &
	1122	num_usm_v_kji(2), &
	1123	.TRUE., .FALSE., &
	1124	.FALSE., .FALSE. )
	1125	ELSEIF ( land_surface ) THEN
	1126	CALL initialize_vertical_surfaces( 2, k, j, i, &
	1127	surf_lsm_v(2), &
	1128	num_lsm_v(2), &
	1129	num_lsm_v_kji(2), &
	1130	.TRUE., .FALSE., &
	1131	.FALSE., .FALSE. )
	1132	ELSE
	1133	CALL initialize_vertical_surfaces( 2, k, j, i, &
	1134	surf_def_v(2), &
	1135	num_def_v(2), &
	1136	num_def_v_kji(2), &
	1137	.TRUE., .FALSE., &
	1138	.FALSE., .FALSE. )
	1139	ENDIF
	1140	ENDIF
	1141	!
	1142	!-- westward-facing surface
	1143	IF ( .NOT. BTEST( wall_flags_0(k,j,i+1), 0 ) ) THEN
	1144	IF ( urban_surface ) THEN
	1145	CALL initialize_vertical_surfaces( 3, k, j, i, &
	1146	surf_usm_v(3), &
	1147	num_usm_v(3), &
	1148	num_usm_v_kji(3), &
	1149	.FALSE., .TRUE., &
	1150	.FALSE., .FALSE. )
	1151	ELSEIF ( land_surface ) THEN
	1152	CALL initialize_vertical_surfaces( 3, k, j, i, &
	1153	surf_lsm_v(3), &
	1154	num_lsm_v(3), &
	1155	num_lsm_v_kji(3), &
	1156	.FALSE., .TRUE., &
	1157	.FALSE., .FALSE. )
	1158	ELSE
	1159	CALL initialize_vertical_surfaces( 3, k, j, i, &
	1160	surf_def_v(3), &
	1161	num_def_v(3), &
	1162	num_def_v_kji(3), &
	1163	.FALSE., .TRUE., &
	1164	.FALSE., .FALSE. )
	1165	ENDIF
	1166	ENDIF
	1167	ENDIF
	1168
	1169
	1170	ENDDO
	1171	!
	1172	!-- Determine start- and end-index at grid point (j,i). Also, for
	1173	!-- horizontal surfaces more than 1 horizontal surface element can
	1174	!-- exist at grid point (j,i) if overhanging structures are present.
	1175	!-- Upward-facing surfaces
	1176	surf_def_h(0)%start_index(j,i) = start_index_def_h(0)
	1177	surf_def_h(0)%end_index(j,i) = surf_def_h(0)%start_index(j,i) + &
	1178	num_def_h_kji(0) - 1
	1179	start_index_def_h(0) = surf_def_h(0)%end_index(j,i) + 1
	1180	!
	1181	!-- Downward-facing surfaces, except model top
	1182	surf_def_h(1)%start_index(j,i) = start_index_def_h(1)
	1183	surf_def_h(1)%end_index(j,i) = surf_def_h(1)%start_index(j,i) + &
	1184	num_def_h_kji(1) - 1
	1185	start_index_def_h(1) = surf_def_h(1)%end_index(j,i) + 1
	1186	!
	1187	!-- Downward-facing surfaces -- model top fluxes
	1188	surf_def_h(2)%start_index(j,i) = start_index_def_h(2)
	1189	surf_def_h(2)%end_index(j,i) = surf_def_h(2)%start_index(j,i) + &
	1190	num_def_h_kji(2) - 1
	1191	start_index_def_h(2) = surf_def_h(2)%end_index(j,i) + 1
	1192	!
	1193	!-- Horizontal natural land surfaces
	1194	surf_lsm_h%start_index(j,i) = start_index_lsm_h
	1195	surf_lsm_h%end_index(j,i) = surf_lsm_h%start_index(j,i) + &
	1196	num_lsm_h_kji - 1
	1197	start_index_lsm_h = surf_lsm_h%end_index(j,i) + 1
	1198	!
	1199	!-- Horizontal urban surfaces
	1200	surf_usm_h%start_index(j,i) = start_index_usm_h
	1201	surf_usm_h%end_index(j,i) = surf_usm_h%start_index(j,i) + &
	1202	num_usm_h_kji - 1
	1203	start_index_usm_h = surf_usm_h%end_index(j,i) + 1
	1204
	1205	!
	1206	!-- Vertical surfaces - Default type
	1207	surf_def_v(0)%start_index(j,i) = start_index_def_v(0)
	1208	surf_def_v(1)%start_index(j,i) = start_index_def_v(1)
	1209	surf_def_v(2)%start_index(j,i) = start_index_def_v(2)
	1210	surf_def_v(3)%start_index(j,i) = start_index_def_v(3)
	1211	surf_def_v(0)%end_index(j,i) = start_index_def_v(0) + &
	1212	num_def_v_kji(0) - 1
	1213	surf_def_v(1)%end_index(j,i) = start_index_def_v(1) + &
	1214	num_def_v_kji(1) - 1
	1215	surf_def_v(2)%end_index(j,i) = start_index_def_v(2) + &
	1216	num_def_v_kji(2) - 1
	1217	surf_def_v(3)%end_index(j,i) = start_index_def_v(3) + &
	1218	num_def_v_kji(3) - 1
	1219	start_index_def_v(0) = surf_def_v(0)%end_index(j,i) + 1
	1220	start_index_def_v(1) = surf_def_v(1)%end_index(j,i) + 1
	1221	start_index_def_v(2) = surf_def_v(2)%end_index(j,i) + 1
	1222	start_index_def_v(3) = surf_def_v(3)%end_index(j,i) + 1
	1223	!
	1224	!-- Natural type
	1225	surf_lsm_v(0)%start_index(j,i) = start_index_lsm_v(0)
	1226	surf_lsm_v(1)%start_index(j,i) = start_index_lsm_v(1)
	1227	surf_lsm_v(2)%start_index(j,i) = start_index_lsm_v(2)
	1228	surf_lsm_v(3)%start_index(j,i) = start_index_lsm_v(3)
	1229	surf_lsm_v(0)%end_index(j,i) = start_index_lsm_v(0) + &
	1230	num_lsm_v_kji(0) - 1
	1231	surf_lsm_v(1)%end_index(j,i) = start_index_lsm_v(1) + &
	1232	num_lsm_v_kji(1) - 1
	1233	surf_lsm_v(2)%end_index(j,i) = start_index_lsm_v(2) + &
	1234	num_lsm_v_kji(2) - 1
	1235	surf_lsm_v(3)%end_index(j,i) = start_index_lsm_v(3) + &
	1236	num_lsm_v_kji(3) - 1
	1237	start_index_lsm_v(0) = surf_lsm_v(0)%end_index(j,i) + 1
	1238	start_index_lsm_v(1) = surf_lsm_v(1)%end_index(j,i) + 1
	1239	start_index_lsm_v(2) = surf_lsm_v(2)%end_index(j,i) + 1
	1240	start_index_lsm_v(3) = surf_lsm_v(3)%end_index(j,i) + 1
	1241	!
	1242	!-- Urban type
	1243	surf_usm_v(0)%start_index(j,i) = start_index_usm_v(0)
	1244	surf_usm_v(1)%start_index(j,i) = start_index_usm_v(1)
	1245	surf_usm_v(2)%start_index(j,i) = start_index_usm_v(2)
	1246	surf_usm_v(3)%start_index(j,i) = start_index_usm_v(3)
	1247	surf_usm_v(0)%end_index(j,i) = start_index_usm_v(0) + &
	1248	num_usm_v_kji(0) - 1
	1249	surf_usm_v(1)%end_index(j,i) = start_index_usm_v(1) + &
	1250	num_usm_v_kji(1) - 1
	1251	surf_usm_v(2)%end_index(j,i) = start_index_usm_v(2) + &
	1252	num_usm_v_kji(2) - 1
	1253	surf_usm_v(3)%end_index(j,i) = start_index_usm_v(3) + &
	1254	num_usm_v_kji(3) - 1
	1255	start_index_usm_v(0) = surf_usm_v(0)%end_index(j,i) + 1
	1256	start_index_usm_v(1) = surf_usm_v(1)%end_index(j,i) + 1
	1257	start_index_usm_v(2) = surf_usm_v(2)%end_index(j,i) + 1
	1258	start_index_usm_v(3) = surf_usm_v(3)%end_index(j,i) + 1
	1259
	1260
	1261	ENDDO
	1262	ENDDO
	1263
	1264	CONTAINS
	1265
	1266	!------------------------------------------------------------------------------!
	1267	! Description:
	1268	! ------------
	1269	!> Initialize horizontal surface elements, upward- and downward-facing.
	1270	!> Note, horizontal surface type alsw comprises model-top fluxes, which are,
	1271	!> initialized in a different routine.
	1272	!------------------------------------------------------------------------------!
	1273	SUBROUTINE initialize_horizontal_surfaces( k, j, i, surf, num_h, &
	1274	num_h_kji, upward_facing, &
	1275	downward_facing )
	1276
	1277	IMPLICIT NONE
	1278
	1279	INTEGER(iwp) :: i !< running index x-direction
	1280	INTEGER(iwp) :: j !< running index y-direction
	1281	INTEGER(iwp) :: k !< running index z-direction
	1282	INTEGER(iwp) :: num_h !< current number of surface element
	1283	INTEGER(iwp) :: num_h_kji !< dummy increment
	1284
	1285	LOGICAL :: upward_facing !< flag indicating upward-facing surface
	1286	LOGICAL :: downward_facing !< flag indicating downward-facing surface
	1287
	1288	TYPE( surf_type ) :: surf !< respective surface type
	1289	!
	1290	!-- Store indices of respective surface element
	1291	surf%i(num_h) = i
	1292	surf%j(num_h) = j
	1293	surf%k(num_h) = k
	1294	!
	1295	!-- Surface orientation, bit 0 is set to 1 for upward-facing surfaces,
	1296	!-- bit 1 is for downward-facing surfaces.
	1297	IF ( upward_facing ) surf%facing(num_h) = IBSET( surf%facing(num_h), 0 )
	1298	IF ( downward_facing ) surf%facing(num_h) = IBSET( surf%facing(num_h), 1 )
	1299	!
	1300	!-- Initialize surface-layer height
	1301	IF ( upward_facing ) THEN
	1302	surf%z_mo(num_h) = zu(k) - zw(k-1)
	1303	ELSE
	1304	surf%z_mo(num_h) = zw(k) - zu(k)
	1305	ENDIF
	1306
	1307	surf%z0(num_h) = roughness_length
	1308	surf%z0h(num_h) = z0h_factor * roughness_length
	1309	surf%z0q(num_h) = z0h_factor * roughness_length
	1310	!
	1311	!-- Initialization in case of 1D pre-cursor run
	1312	IF ( INDEX( initializing_actions, 'set_1d-model_profiles' ) /= 0 )&
	1313	THEN
	1314	IF ( .NOT. constant_diffusion ) THEN
	1315	IF ( constant_flux_layer ) THEN
	1316	surf%ol(num_h) = surf%z_mo(num_h) / &
	1317	( rif1d(nzb+1) + 1.0E-20_wp )
	1318	surf%us(num_h) = us1d
	1319	surf%usws(num_h) = usws1d
	1320	surf%vsws(num_h) = vsws1d
	1321	ELSE
	1322	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
	1323	surf%us(num_h) = 0.0_wp
	1324	surf%usws(num_h) = 0.0_wp
	1325	surf%vsws(num_h) = 0.0_wp
	1326	ENDIF
	1327	ELSE
	1328	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
	1329	surf%us(num_h) = 0.0_wp
	1330	surf%usws(num_h) = 0.0_wp
	1331	surf%vsws(num_h) = 0.0_wp
	1332	ENDIF
	1333	!
	1334	!-- Initialization in case of constant profiles
	1335	ELSEIF ( INDEX(initializing_actions, 'set_constant_profiles') /= 0 )&
	1336	THEN
	1337
	1338	surf%ol(num_h) = surf%z_mo(num_h) / zeta_min
	1339	!
	1340	!-- Very small number is required for calculation of Obukhov length
	1341	!-- at first timestep
	1342	surf%us(num_h) = 1E-30_wp
	1343	surf%usws(num_h) = 0.0_wp
	1344	surf%vsws(num_h) = 0.0_wp
	1345
	1346	ENDIF
	1347
	1348	surf%rib(num_h) = 0.0_wp
	1349	surf%uvw_abs(num_h) = 0.0_wp
	1350
	1351	IF ( .NOT. constant_diffusion ) THEN
	1352	surf%u_0(num_h) = 0.0_wp
	1353	surf%v_0(num_h) = 0.0_wp
	1354	ENDIF
	1355
	1356	surf%ts(num_h) = 0.0_wp
	1357
	1358	IF ( humidity ) THEN
	1359	surf%qs(num_h) = 0.0_wp
[2292]	1360	IF ( cloud_physics .AND. microphysics_morrison) THEN
	1361	surf%qcs(num_h) = 0.0_wp
	1362	surf%ncs(num_h) = 0.0_wp
	1363
	1364	surf%qcsws(num_h) = 0.0_wp
	1365	surf%ncsws(num_h) = 0.0_wp
	1366
	1367	ENDIF
[2232]	1368	IF ( cloud_physics .AND. microphysics_seifert) THEN
	1369	surf%qrs(num_h) = 0.0_wp
	1370	surf%nrs(num_h) = 0.0_wp
	1371
	1372	surf%qrsws(num_h) = 0.0_wp
	1373	surf%nrsws(num_h) = 0.0_wp
	1374
	1375	surf%pt1(num_h) = 0.0_wp
	1376	surf%qv1(num_h) = 0.0_wp
	1377
	1378	ENDIF
	1379	ENDIF
	1380
	1381	IF ( passive_scalar ) surf%ss(num_h) = 0.0_wp
	1382	!
	1383	!-- Inititalize surface fluxes of sensible and latent heat, as well as
	1384	!-- passive scalar
	1385	IF ( use_surface_fluxes ) THEN
	1386
	1387	IF ( upward_facing ) THEN
	1388	IF ( constant_heatflux ) THEN
	1389	!
	1390	!-- Initialize surface heatflux. However, skip this for now if
	1391	!-- if random_heatflux is set. This case, shf is initialized later.
	1392	IF ( .NOT. random_heatflux ) THEN
	1393	surf%shf(num_h) = surface_heatflux * &
	1394	heatflux_input_conversion(nzb)
	1395	!
	1396	!-- Check if surface heat flux might be replaced by
	1397	!-- prescribed wall heatflux
	1398	IF ( k-1 /= 0 ) THEN
	1399	surf%shf(num_h) = wall_heatflux(0) * &
	1400	heatflux_input_conversion(k-1)
	1401	ENDIF
	1402	!
	1403	!-- Initialize shf with data from external file LSF_DATA. Will
	1404	!-- be done directly in ls_foring_surf
	1405	!-- Attention: Just a workaround, need to be revised!!!
	1406	IF ( large_scale_forcing .AND. lsf_surf ) THEN
	1407	! CALL ls_forcing_surf ( simulated_time )
	1408	! surf%shf(num_h) = shf(j,i)
	1409	ENDIF
	1410	ENDIF
	1411	ELSE
	1412	surf%shf(num_h) = 0.0_wp
	1413	ENDIF
	1414	!
[2256]	1415	!-- Set heat-flux at downward-facing surfaces
[2232]	1416	ELSE
[2256]	1417	surf%shf(num_h) = wall_heatflux(5) * &
	1418	heatflux_input_conversion(k)
[2232]	1419	ENDIF
	1420
	1421	IF ( humidity ) THEN
	1422	IF ( upward_facing ) THEN
	1423	IF ( constant_waterflux ) THEN
	1424	surf%qsws(num_h) = surface_waterflux * &
	1425	waterflux_input_conversion(nzb)
	1426	IF ( k-1 /= 0 ) THEN
	1427	surf%qsws(num_h) = wall_humidityflux(0) * &
	1428	waterflux_input_conversion(k-1)
	1429	ENDIF
	1430	ELSE
	1431	surf%qsws(num_h) = 0.0_wp
	1432	ENDIF
	1433	ELSE
[2256]	1434	surf%qsws(num_h) = wall_humidityflux(5) * &
	1435	heatflux_input_conversion(k)
[2232]	1436	ENDIF
	1437	ENDIF
	1438
	1439	IF ( passive_scalar ) THEN
	1440	IF ( upward_facing ) THEN
	1441	IF ( constant_scalarflux ) THEN
	1442	surf%ssws(num_h) = surface_scalarflux
	1443
	1444	IF ( k-1 /= 0 ) &
[2256]	1445	surf%ssws(num_h) = wall_scalarflux(0)
[2232]	1446
	1447	ELSE
	1448	surf%ssws(num_h) = 0.0_wp
	1449	ENDIF
	1450	ELSE
[2256]	1451	surf%ssws(num_h) = wall_scalarflux(5)
[2232]	1452	ENDIF
	1453	ENDIF
	1454
	1455	IF ( ocean ) THEN
	1456	IF ( upward_facing ) THEN
	1457	surf%sasws(num_h) = bottom_salinityflux
	1458	ELSE
	1459	surf%sasws(num_h) = 0.0_wp
	1460	ENDIF
	1461	ENDIF
	1462	ENDIF
	1463	!
	1464	!-- Increment surface indices
	1465	num_h = num_h + 1
	1466	num_h_kji = num_h_kji + 1
	1467
	1468
	1469	END SUBROUTINE initialize_horizontal_surfaces
	1470
	1471
	1472	!------------------------------------------------------------------------------!
	1473	! Description:
	1474	! ------------
	1475	!> Initialize model-top fluxes. Currently, only the heatflux and salinity flux
	1476	!> can be prescribed, latent flux is zero in this case!
	1477	!------------------------------------------------------------------------------!
	1478	SUBROUTINE initialize_top( k, j, i, surf, num_h, num_h_kji )
	1479
	1480	IMPLICIT NONE
	1481
	1482	INTEGER(iwp) :: i !< running index x-direction
	1483	INTEGER(iwp) :: j !< running index y-direction
	1484	INTEGER(iwp) :: k !< running index z-direction
	1485	INTEGER(iwp) :: num_h !< current number of surface element
	1486	INTEGER(iwp) :: num_h_kji !< dummy increment
	1487
	1488	TYPE( surf_type ) :: surf !< respective surface type
	1489	!
	1490	!-- Store indices of respective surface element
	1491	surf%i(num_h) = i
	1492	surf%j(num_h) = j
	1493	surf%k(num_h) = k
	1494	!
	1495	!-- Initialize top heat flux
	1496	IF ( constant_top_heatflux ) &
[2478]	1497	surf%shf(num_h) = top_heatflux * heatflux_input_conversion(nzt+1)
[2232]	1498	!
	1499	!-- Initialization in case of a coupled model run
	1500	IF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
[2478]	1501	surf%shf(num_h) = 0.0_wp
[2232]	1502	ENDIF
	1503	!
	1504	!-- Prescribe latent heat flux at the top
	1505	IF ( humidity ) THEN
[2478]	1506	surf%qsws(num_h) = 0.0_wp
[2292]	1507	IF ( cloud_physics .AND. microphysics_morrison ) THEN
[2478]	1508	surf%ncsws(num_h) = 0.0_wp
	1509	surf%qcsws(num_h) = 0.0_wp
[2292]	1510	ENDIF
[2232]	1511	IF ( cloud_physics .AND. microphysics_seifert ) THEN
[2478]	1512	surf%nrsws(num_h) = 0.0_wp
	1513	surf%qrsws(num_h) = 0.0_wp
[2232]	1514	ENDIF
	1515	ENDIF
	1516	!
	1517	!-- Prescribe top scalar flux
	1518	IF ( passive_scalar .AND. constant_top_scalarflux ) &
[2478]	1519	surf%ssws(num_h) = top_scalarflux
[2232]	1520	!
	1521	!-- Prescribe top salinity flux
[2478]	1522	IF ( ocean .AND. constant_top_salinityflux) &
	1523	surf%sasws(num_h) = top_salinityflux
[2232]	1524	!
	1525	!-- Initialization in case of a coupled model run
	1526	IF ( coupling_mode == 'ocean_to_atmosphere' ) THEN
[2478]	1527	surf%shf(num_h) = 0.0_wp
[2232]	1528	ENDIF
	1529	!
	1530	!-- Top momentum fluxes
[2478]	1531	surf%usws(num_h) = top_momentumflux_u * &
	1532	momentumflux_input_conversion(nzt+1)
	1533	surf%vsws(num_h) = top_momentumflux_v * &
	1534	momentumflux_input_conversion(nzt+1)
[2232]	1535	!
	1536	!-- Increment surface indices
	1537	num_h = num_h + 1
	1538	num_h_kji = num_h_kji + 1
	1539
	1540
	1541	END SUBROUTINE initialize_top
	1542
	1543
	1544	!------------------------------------------------------------------------------!
	1545	! Description:
	1546	! ------------
	1547	!> Initialize vertical surface elements.
	1548	!------------------------------------------------------------------------------!
	1549	SUBROUTINE initialize_vertical_surfaces( l, k, j, i, surf, num_v, &
	1550	num_v_kji, east_facing, &
	1551	west_facing, south_facing, &
	1552	north_facing )
	1553
	1554	IMPLICIT NONE
	1555
	1556	INTEGER(iwp) :: component !<
	1557	INTEGER(iwp) :: i !< running index x-direction
	1558	INTEGER(iwp) :: j !< running index x-direction
	1559	INTEGER(iwp) :: k !< running index x-direction
	1560	INTEGER(iwp) :: l !< index variable for the surface type, indicating the facing
	1561	INTEGER(iwp) :: num_v !< current number of surface element
	1562	INTEGER(iwp) :: num_v_kji !< current number of surface element at (j,i)
	1563
	1564
	1565	LOGICAL :: east_facing !< flag indicating east-facing surfaces
	1566	LOGICAL :: north_facing !< flag indicating north-facing surfaces
	1567	LOGICAL :: south_facing !< flag indicating south-facing surfaces
	1568	LOGICAL :: west_facing !< flag indicating west-facing surfaces
	1569
	1570	TYPE( surf_type ) :: surf !< respective surface type
	1571
	1572	!
	1573	!-- Store indices of respective wall element
	1574	surf%i(num_v) = i
	1575	surf%j(num_v) = j
	1576	surf%k(num_v) = k
	1577	!
	1578	!-- Initialize surface-layer height, or more precisely, distance to surface
	1579	IF ( north_facing .OR. south_facing ) THEN
	1580	surf%z_mo(num_v) = 0.5_wp * dy
	1581	ELSE
	1582	surf%z_mo(num_v) = 0.5_wp * dx
	1583	ENDIF
	1584
	1585	surf%facing(num_v) = 0
	1586	!
	1587	!-- Surface orientation. Moreover, set component id to map wall_heatflux,
	1588	!-- etc., on surface type (further below)
	1589	IF ( north_facing ) THEN
	1590	surf%facing(num_v) = IBSET( surf%facing(num_v), 0 )
	1591	component = 4
	1592	ENDIF
	1593
	1594	IF ( south_facing ) THEN
	1595	surf%facing(num_v) = IBSET( surf%facing(num_v), 1 )
	1596	component = 3
	1597	ENDIF
	1598
	1599	IF ( east_facing ) THEN
	1600	surf%facing(num_v) = IBSET( surf%facing(num_v), 2 )
	1601	component = 2
	1602	ENDIF
	1603
	1604	IF ( west_facing ) THEN
	1605	surf%facing(num_v) = IBSET( surf%facing(num_v), 3 )
	1606	component = 1
	1607	ENDIF
	1608
	1609
	1610	surf%z0(num_v) = roughness_length
	1611	surf%z0h(num_v) = z0h_factor * roughness_length
	1612	surf%z0q(num_v) = z0h_factor * roughness_length
	1613
	1614	surf%us(num_v) = 0.0_wp
	1615	!
	1616	!-- If required, initialize Obukhov length
	1617	IF ( ALLOCATED( surf%ol ) ) &
	1618	surf%ol(num_v) = surf%z_mo(num_v) / zeta_min
	1619
	1620	surf%uvw_abs(num_v) = 0.0_wp
	1621
	1622	surf%mom_flux_uv(num_v) = 0.0_wp
	1623	surf%mom_flux_w(num_v) = 0.0_wp
	1624	surf%mom_flux_tke(0:1,num_v) = 0.0_wp
	1625
	1626	surf%ts(num_v) = 0.0_wp
	1627	surf%shf(num_v) = wall_heatflux(component)
	1628
	1629	IF ( humidity ) THEN
	1630	surf%qs(num_v) = 0.0_wp
	1631	surf%qsws(num_v) = wall_humidityflux(component)
	1632	!
	1633	!-- Following wall fluxes are assumed to be zero
[2292]	1634	IF ( cloud_physics .AND. microphysics_morrison) THEN
	1635	surf%qcs(num_v) = 0.0_wp
	1636	surf%ncs(num_v) = 0.0_wp
	1637
	1638	surf%qcsws(num_v) = 0.0_wp
	1639	surf%ncsws(num_v) = 0.0_wp
	1640	ENDIF
[2232]	1641	IF ( cloud_physics .AND. microphysics_seifert) THEN
	1642	surf%qrs(num_v) = 0.0_wp
	1643	surf%nrs(num_v) = 0.0_wp
	1644
	1645	surf%qrsws(num_v) = 0.0_wp
	1646	surf%nrsws(num_v) = 0.0_wp
	1647	ENDIF
	1648	ENDIF
	1649
	1650	IF ( passive_scalar ) THEN
	1651	surf%ss(num_v) = 0.0_wp
	1652	surf%ssws(num_v) = wall_scalarflux(component)
	1653	ENDIF
	1654	!
	1655	!-- So far, salinityflux at vertical surfaces is simply zero
	1656	!-- at the moment
	1657	IF ( ocean ) surf%sasws(num_v) = wall_salinityflux(component)
	1658	!
	1659	!-- Increment wall indices
	1660	num_v = num_v + 1
	1661	num_v_kji = num_v_kji + 1
	1662
	1663	END SUBROUTINE initialize_vertical_surfaces
	1664
	1665	END SUBROUTINE init_surfaces
	1666
[2317]	1667
[2232]	1668	!------------------------------------------------------------------------------!
	1669	! Description:
	1670	! ------------
[2317]	1671	!> Determines topography-top index at given (j,i)-position.
	1672	!------------------------------------------------------------------------------!
	1673	FUNCTION get_topography_top_index( j, i, grid )
	1674
	1675	USE kinds
	1676
	1677	IMPLICIT NONE
	1678
	1679	CHARACTER(LEN=*) :: grid !< flag to distinquish between staggered grids
	1680	INTEGER(iwp) :: i !< grid index in x-dimension
	1681	INTEGER(iwp) :: ibit !< bit position where topography information is stored on respective grid
	1682	INTEGER(iwp) :: j !< grid index in y-dimension
	1683	INTEGER(iwp) :: get_topography_top_index !< topography top index
	1684
	1685	SELECT CASE ( TRIM( grid ) )
	1686
	1687	CASE ( 's' )
	1688	ibit = 12
	1689	CASE ( 'u' )
	1690	ibit = 14
	1691	CASE ( 'v' )
	1692	ibit = 16
	1693	CASE ( 'w' )
	1694	ibit = 18
	1695	CASE ( 's_out' )
	1696	ibit = 24
	1697	CASE ( 'u_out' )
	1698	ibit = 26
	1699	CASE ( 'v_out' )
	1700	ibit = 27
	1701	CASE ( 'w_out' )
	1702	ibit = 28
	1703	CASE DEFAULT
	1704	!
	1705	!-- Set default to scalar grid
	1706	ibit = 12
	1707
	1708	END SELECT
	1709
	1710	get_topography_top_index = MAXLOC( &
	1711	MERGE( 1, 0, &
	1712	BTEST( wall_flags_0(:,j,i), ibit ) &
	1713	), DIM = 1 &
	1714	) - 1
	1715
	1716	RETURN
	1717
	1718	END FUNCTION get_topography_top_index
	1719
	1720	!------------------------------------------------------------------------------!
	1721	! Description:
	1722	! ------------
[2232]	1723	!> Gathers all surface elements with the same facing (but possibly different
	1724	!> type) onto a surface type, and writes binary data into restart files.
	1725	!------------------------------------------------------------------------------!
	1726	SUBROUTINE surface_write_restart_data
	1727
	1728	IMPLICIT NONE
	1729
[2269]	1730	CHARACTER(LEN=1) :: dum !< dummy string to create output-variable name
[2232]	1731
[2269]	1732	INTEGER(iwp) :: i !< running index x-direction
	1733	INTEGER(iwp) :: j !< running index y-direction
	1734	INTEGER(iwp) :: l !< index surface type orientation
	1735	INTEGER(iwp) :: m !< running index for surface elements on individual surface array
	1736	INTEGER(iwp), DIMENSION(0:3) :: mm !< running index for surface elements on gathered surface array
[2232]	1737
[2269]	1738	TYPE(surf_type), DIMENSION(0:2) :: surf_h !< gathered horizontal surfaces, contains all surface types
	1739	TYPE(surf_type), DIMENSION(0:3) :: surf_v !< gathered vertical surfaces, contains all surface types
[2232]	1740
	1741	!
	1742	!-- Determine total number of horizontal and vertical surface elements before
	1743	!-- writing var_list
	1744	CALL surface_last_actions
	1745	!
	1746	!-- Count number of grid points with same facing and allocate attributes respectively
	1747	!-- Horizontal upward facing
	1748	surf_h(0)%ns = ns_h_on_file(0)
	1749	CALL allocate_surface_attributes_h( surf_h(0), nys, nyn, nxl, nxr )
	1750	!
	1751	!-- Horizontal downward facing
	1752	surf_h(1)%ns = ns_h_on_file(1)
	1753	CALL allocate_surface_attributes_h( surf_h(1), nys, nyn, nxl, nxr )
	1754	!
	1755	!-- Model top
	1756	surf_h(2)%ns = ns_h_on_file(2)
	1757	CALL allocate_surface_attributes_h_top( surf_h(2), nys, nyn, nxl, nxr )
	1758	!
	1759	!-- Vertical surfaces
	1760	DO l = 0, 3
	1761	surf_v(l)%ns = ns_v_on_file(l)
	1762	CALL allocate_surface_attributes_v( surf_v(l), .FALSE., &
	1763	nys, nyn, nxl, nxr )
	1764	ENDDO
	1765	!
	1766	!-- In the following, gather data from surfaces elements with the same
	1767	!-- facing (but possibly differt type) on 1 data-type array.
	1768	mm(0:2) = 1
	1769	DO l = 0, 2
	1770	DO i = nxl, nxr
	1771	DO j = nys, nyn
	1772	DO m = surf_def_h(l)%start_index(j,i), &
	1773	surf_def_h(l)%end_index(j,i)
	1774	IF ( ALLOCATED( surf_def_h(l)%us ) ) &
	1775	surf_h(l)%us(mm(l)) = surf_def_h(l)%us(m)
	1776	IF ( ALLOCATED( surf_def_h(l)%ts ) ) &
	1777	surf_h(l)%ts(mm(l)) = surf_def_h(l)%ts(m)
	1778	IF ( ALLOCATED( surf_def_h(l)%qs ) ) &
	1779	surf_h(l)%qs(mm(l)) = surf_def_h(l)%qs(m)
	1780	IF ( ALLOCATED( surf_def_h(l)%ss ) ) &
	1781	surf_h(l)%ss(mm(l)) = surf_def_h(l)%ss(m)
[2292]	1782	IF ( ALLOCATED( surf_def_h(l)%qcs ) ) &
	1783	surf_h(l)%qcs(mm(l)) = surf_def_h(l)%qcs(m)
	1784	IF ( ALLOCATED( surf_def_h(l)%ncs ) ) &
	1785	surf_h(l)%ncs(mm(l)) = surf_def_h(l)%ncs(m)
[2232]	1786	IF ( ALLOCATED( surf_def_h(l)%qrs ) ) &
	1787	surf_h(l)%qrs(mm(l)) = surf_def_h(l)%qrs(m)
	1788	IF ( ALLOCATED( surf_def_h(l)%nrs ) ) &
	1789	surf_h(l)%nrs(mm(l)) = surf_def_h(l)%nrs(m)
	1790	IF ( ALLOCATED( surf_def_h(l)%ol ) ) &
	1791	surf_h(l)%ol(mm(l)) = surf_def_h(l)%ol(m)
	1792	IF ( ALLOCATED( surf_def_h(l)%rib ) ) &
	1793	surf_h(l)%rib(mm(l)) = surf_def_h(l)%rib(m)
	1794	IF ( ALLOCATED( surf_def_h(l)%usws ) ) &
	1795	surf_h(l)%usws(mm(l)) = surf_def_h(l)%usws(m)
	1796	IF ( ALLOCATED( surf_def_h(l)%vsws ) ) &
	1797	surf_h(l)%vsws(mm(l)) = surf_def_h(l)%vsws(m)
	1798	IF ( ALLOCATED( surf_def_h(l)%shf ) ) &
	1799	surf_h(l)%shf(mm(l)) = surf_def_h(l)%shf(m)
	1800	IF ( ALLOCATED( surf_def_h(l)%qsws ) ) &
	1801	surf_h(l)%qsws(mm(l)) = surf_def_h(l)%qsws(m)
	1802	IF ( ALLOCATED( surf_def_h(l)%ssws ) ) &
[2378]	1803	surf_h(l)%ssws(mm(l)) = surf_def_h(l)%ssws(m)
[2292]	1804	IF ( ALLOCATED( surf_def_h(l)%ncsws ) ) &
	1805	surf_h(l)%ncsws(mm(l)) = surf_def_h(l)%ncsws(m)
[2232]	1806	IF ( ALLOCATED( surf_def_h(l)%nrsws ) ) &
	1807	surf_h(l)%nrsws(mm(l)) = surf_def_h(l)%nrsws(m)
	1808	IF ( ALLOCATED( surf_def_h(l)%sasws ) ) &
	1809	surf_h(l)%sasws(mm(l)) = surf_def_h(l)%sasws(m)
	1810
	1811	mm(l) = mm(l) + 1
	1812	ENDDO
	1813
	1814	IF ( l == 0 ) THEN
	1815	DO m = surf_lsm_h%start_index(j,i), &
	1816	surf_lsm_h%end_index(j,i)
	1817	IF ( ALLOCATED( surf_lsm_h%us ) ) &
	1818	surf_h(0)%us(mm(0)) = surf_lsm_h%us(m)
	1819	IF ( ALLOCATED( surf_lsm_h%ts ) ) &
	1820	surf_h(0)%ts(mm(0)) = surf_lsm_h%ts(m)
	1821	IF ( ALLOCATED( surf_lsm_h%qs ) ) &
	1822	surf_h(0)%qs(mm(0)) = surf_lsm_h%qs(m)
	1823	IF ( ALLOCATED( surf_lsm_h%ss ) ) &
	1824	surf_h(0)%ss(mm(0)) = surf_lsm_h%ss(m)
[2292]	1825	IF ( ALLOCATED( surf_lsm_h%qcs ) ) &
	1826	surf_h(0)%qcs(mm(0)) = surf_lsm_h%qcs(m)
	1827	IF ( ALLOCATED( surf_lsm_h%ncs ) ) &
	1828	surf_h(0)%ncs(mm(0)) = surf_lsm_h%ncs(m)
[2232]	1829	IF ( ALLOCATED( surf_lsm_h%qrs ) ) &
	1830	surf_h(0)%qrs(mm(0)) = surf_lsm_h%qrs(m)
	1831	IF ( ALLOCATED( surf_lsm_h%nrs ) ) &
	1832	surf_h(0)%nrs(mm(0)) = surf_lsm_h%nrs(m)
	1833	IF ( ALLOCATED( surf_lsm_h%ol ) ) &
	1834	surf_h(0)%ol(mm(0)) = surf_lsm_h%ol(m)
	1835	IF ( ALLOCATED( surf_lsm_h%rib ) ) &
	1836	surf_h(0)%rib(mm(0)) = surf_lsm_h%rib(m)
	1837	IF ( ALLOCATED( surf_lsm_h%usws ) ) &
	1838	surf_h(0)%usws(mm(0)) = surf_lsm_h%usws(m)
	1839	IF ( ALLOCATED( surf_lsm_h%vsws ) ) &
	1840	surf_h(0)%vsws(mm(0)) = surf_lsm_h%vsws(m)
	1841	IF ( ALLOCATED( surf_lsm_h%shf ) ) &
	1842	surf_h(0)%shf(mm(0)) = surf_lsm_h%shf(m)
	1843	IF ( ALLOCATED( surf_lsm_h%qsws ) ) &
	1844	surf_h(0)%qsws(mm(0)) = surf_lsm_h%qsws(m)
	1845	IF ( ALLOCATED( surf_lsm_h%ssws ) ) &
[2378]	1846	surf_h(0)%ssws(mm(0)) = surf_lsm_h%ssws(m)
[2292]	1847	IF ( ALLOCATED( surf_lsm_h%ncsws ) ) &
	1848	surf_h(0)%ncsws(mm(0)) = surf_lsm_h%ncsws(m)
[2232]	1849	IF ( ALLOCATED( surf_lsm_h%nrsws ) ) &
	1850	surf_h(0)%nrsws(mm(0)) = surf_lsm_h%nrsws(m)
	1851	IF ( ALLOCATED( surf_lsm_h%sasws ) ) &
	1852	surf_h(0)%sasws(mm(0)) = surf_lsm_h%sasws(m)
	1853
	1854	mm(0) = mm(0) + 1
	1855
	1856	ENDDO
	1857
	1858	DO m = surf_usm_h%start_index(j,i), &
	1859	surf_usm_h%end_index(j,i)
	1860	IF ( ALLOCATED( surf_usm_h%us ) ) &
	1861	surf_h(0)%us(mm(0)) = surf_usm_h%us(m)
	1862	IF ( ALLOCATED( surf_usm_h%ts ) ) &
	1863	surf_h(0)%ts(mm(0)) = surf_usm_h%ts(m)
	1864	IF ( ALLOCATED( surf_usm_h%qs ) ) &
	1865	surf_h(0)%qs(mm(0)) = surf_usm_h%qs(m)
	1866	IF ( ALLOCATED( surf_usm_h%ss ) ) &
	1867	surf_h(0)%ss(mm(0)) = surf_usm_h%ss(m)
[2292]	1868	IF ( ALLOCATED( surf_usm_h%qcs ) ) &
	1869	surf_h(0)%qcs(mm(0)) = surf_usm_h%qcs(m)
	1870	IF ( ALLOCATED( surf_usm_h%ncs ) ) &
	1871	surf_h(0)%ncs(mm(0)) = surf_usm_h%ncs(m)
[2232]	1872	IF ( ALLOCATED( surf_usm_h%qrs ) ) &
	1873	surf_h(0)%qrs(mm(0)) = surf_usm_h%qrs(m)
	1874	IF ( ALLOCATED( surf_usm_h%nrs ) ) &
	1875	surf_h(0)%nrs(mm(0)) = surf_usm_h%nrs(m)
	1876	IF ( ALLOCATED( surf_usm_h%ol ) ) &
	1877	surf_h(0)%ol(mm(0)) = surf_usm_h%ol(m)
	1878	IF ( ALLOCATED( surf_usm_h%rib ) ) &
	1879	surf_h(0)%rib(mm(0)) = surf_usm_h%rib(m)
	1880	IF ( ALLOCATED( surf_usm_h%usws ) ) &
	1881	surf_h(0)%usws(mm(0)) = surf_usm_h%usws(m)
	1882	IF ( ALLOCATED( surf_usm_h%vsws ) ) &
	1883	surf_h(0)%vsws(mm(0)) = surf_usm_h%vsws(m)
	1884	IF ( ALLOCATED( surf_usm_h%shf ) ) &
	1885	surf_h(0)%shf(mm(0)) = surf_usm_h%shf(m)
	1886	IF ( ALLOCATED( surf_usm_h%qsws ) ) &
	1887	surf_h(0)%qsws(mm(0)) = surf_usm_h%qsws(m)
	1888	IF ( ALLOCATED( surf_usm_h%ssws ) ) &
[2378]	1889	surf_h(0)%ssws(mm(0)) = surf_usm_h%ssws(m)
[2292]	1890	IF ( ALLOCATED( surf_usm_h%ncsws ) ) &
	1891	surf_h(0)%ncsws(mm(0)) = surf_usm_h%ncsws(m)
[2232]	1892	IF ( ALLOCATED( surf_usm_h%nrsws ) ) &
	1893	surf_h(0)%nrsws(mm(0)) = surf_usm_h%nrsws(m)
	1894	IF ( ALLOCATED( surf_usm_h%sasws ) ) &
	1895	surf_h(0)%sasws(mm(0)) = surf_usm_h%sasws(m)
	1896
	1897	mm(0) = mm(0) + 1
	1898
	1899	ENDDO
	1900
	1901
	1902	ENDIF
	1903
	1904	ENDDO
	1905
	1906	ENDDO
	1907	IF ( l == 0 ) THEN
	1908	surf_h(l)%start_index = MAX( surf_def_h(l)%start_index, &
	1909	surf_lsm_h%start_index, &
	1910	surf_usm_h%start_index )
	1911	surf_h(l)%end_index = MAX( surf_def_h(l)%end_index, &
	1912	surf_lsm_h%end_index, &
	1913	surf_usm_h%end_index )
	1914	ELSE
	1915	surf_h(l)%start_index = surf_def_h(l)%start_index
	1916	surf_h(l)%end_index = surf_def_h(l)%end_index
	1917	ENDIF
	1918	ENDDO
	1919
	1920
	1921	mm(0:3) = 1
	1922	DO l = 0, 3
	1923	DO i = nxl, nxr
	1924	DO j = nys, nyn
	1925	DO m = surf_def_v(l)%start_index(j,i), &
	1926	surf_def_v(l)%end_index(j,i)
	1927	IF ( ALLOCATED( surf_def_v(l)%us ) ) &
	1928	surf_v(l)%us(mm(l)) = surf_def_v(l)%us(m)
	1929	IF ( ALLOCATED( surf_def_v(l)%ts ) ) &
	1930	surf_v(l)%ts(mm(l)) = surf_def_v(l)%ts(m)
	1931	IF ( ALLOCATED( surf_def_v(l)%qs ) ) &
	1932	surf_v(l)%qs(mm(l)) = surf_def_v(l)%qs(m)
	1933	IF ( ALLOCATED( surf_def_v(l)%ss ) ) &
	1934	surf_v(l)%ss(mm(l)) = surf_def_v(l)%ss(m)
[2292]	1935	IF ( ALLOCATED( surf_def_v(l)%qcs ) ) &
	1936	surf_v(l)%qcs(mm(l)) = surf_def_v(l)%qcs(m)
	1937	IF ( ALLOCATED( surf_def_v(l)%ncs ) ) &
	1938	surf_v(l)%ncs(mm(l)) = surf_def_v(l)%ncs(m)
[2232]	1939	IF ( ALLOCATED( surf_def_v(l)%qrs ) ) &
	1940	surf_v(l)%qrs(mm(l)) = surf_def_v(l)%qrs(m)
	1941	IF ( ALLOCATED( surf_def_v(l)%nrs ) ) &
	1942	surf_v(l)%nrs(mm(l)) = surf_def_v(l)%nrs(m)
	1943	IF ( ALLOCATED( surf_def_v(l)%ol ) ) &
	1944	surf_v(l)%ol(mm(l)) = surf_def_v(l)%ol(m)
	1945	IF ( ALLOCATED( surf_def_v(l)%rib ) ) &
	1946	surf_v(l)%rib(mm(l)) = surf_def_v(l)%rib(m)
	1947	IF ( ALLOCATED( surf_def_v(l)%shf ) ) &
	1948	surf_v(l)%shf(mm(l)) = surf_def_v(l)%shf(m)
	1949	IF ( ALLOCATED( surf_def_v(l)%qsws ) ) &
	1950	surf_v(l)%qsws(mm(l)) = surf_def_v(l)%qsws(m)
	1951	IF ( ALLOCATED( surf_def_v(l)%ssws ) ) &
[2378]	1952	surf_v(l)%ssws(mm(l)) = surf_def_v(l)%ssws(m)
[2292]	1953	IF ( ALLOCATED( surf_def_v(l)%ncsws ) ) &
	1954	surf_v(l)%ncsws(mm(l)) = surf_def_v(l)%ncsws(m)
[2232]	1955	IF ( ALLOCATED( surf_def_v(l)%nrsws ) ) &
	1956	surf_v(l)%nrsws(mm(l)) = surf_def_v(l)%nrsws(m)
	1957	IF ( ALLOCATED( surf_def_v(l)%sasws ) ) &
	1958	surf_v(l)%sasws(mm(l)) = surf_def_v(l)%sasws(m)
	1959	IF ( ALLOCATED( surf_def_v(l)%mom_flux_uv) ) &
	1960	surf_v(l)%mom_flux_uv(mm(l)) = surf_def_v(l)%mom_flux_uv(m)
	1961	IF ( ALLOCATED( surf_def_v(l)%mom_flux_w) ) &
	1962	surf_v(l)%mom_flux_w(mm(l)) = surf_def_v(l)%mom_flux_w(m)
	1963	IF ( ALLOCATED( surf_def_v(l)%mom_flux_tke) ) &
	1964	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_def_v(l)%mom_flux_tke(0:1,m)
	1965
	1966	mm(l) = mm(l) + 1
	1967	ENDDO
	1968
	1969	DO m = surf_lsm_v(l)%start_index(j,i), &
	1970	surf_lsm_v(l)%end_index(j,i)
	1971	IF ( ALLOCATED( surf_lsm_v(l)%us ) ) &
	1972	surf_v(l)%us(mm(l)) = surf_lsm_v(l)%us(m)
	1973	IF ( ALLOCATED( surf_lsm_v(l)%ts ) ) &
	1974	surf_v(l)%ts(mm(l)) = surf_lsm_v(l)%ts(m)
	1975	IF ( ALLOCATED( surf_lsm_v(l)%qs ) ) &
	1976	surf_v(l)%qs(mm(l)) = surf_lsm_v(l)%qs(m)
	1977	IF ( ALLOCATED( surf_lsm_v(l)%ss ) ) &
	1978	surf_v(l)%ss(mm(l)) = surf_lsm_v(l)%ss(m)
[2292]	1979	IF ( ALLOCATED( surf_lsm_v(l)%qcs ) ) &
	1980	surf_v(l)%qcs(mm(l)) = surf_lsm_v(l)%qcs(m)
	1981	IF ( ALLOCATED( surf_lsm_v(l)%ncs ) ) &
	1982	surf_v(l)%ncs(mm(l)) = surf_lsm_v(l)%ncs(m)
[2232]	1983	IF ( ALLOCATED( surf_lsm_v(l)%qrs ) ) &
	1984	surf_v(l)%qrs(mm(l)) = surf_lsm_v(l)%qrs(m)
	1985	IF ( ALLOCATED( surf_lsm_v(l)%nrs ) ) &
	1986	surf_v(l)%nrs(mm(l)) = surf_lsm_v(l)%nrs(m)
	1987	IF ( ALLOCATED( surf_lsm_v(l)%ol ) ) &
	1988	surf_v(l)%ol(mm(l)) = surf_lsm_v(l)%ol(m)
	1989	IF ( ALLOCATED( surf_lsm_v(l)%rib ) ) &
	1990	surf_v(l)%rib(mm(l)) = surf_lsm_v(l)%rib(m)
	1991	IF ( ALLOCATED( surf_lsm_v(l)%usws ) ) &
	1992	surf_v(l)%usws(mm(l)) = surf_lsm_v(l)%usws(m)
	1993	IF ( ALLOCATED( surf_lsm_v(l)%vsws ) ) &
	1994	surf_v(l)%vsws(mm(l)) = surf_lsm_v(l)%vsws(m)
	1995	IF ( ALLOCATED( surf_lsm_v(l)%shf ) ) &
	1996	surf_v(l)%shf(mm(l)) = surf_lsm_v(l)%shf(m)
	1997	IF ( ALLOCATED( surf_lsm_v(l)%qsws ) ) &
	1998	surf_v(l)%qsws(mm(l)) = surf_lsm_v(l)%qsws(m)
	1999	IF ( ALLOCATED( surf_lsm_v(l)%ssws ) ) &
[2378]	2000	surf_v(l)%ssws(mm(l)) = surf_lsm_v(l)%ssws(m)
[2292]	2001	IF ( ALLOCATED( surf_lsm_v(l)%ncsws ) ) &
	2002	surf_v(l)%ncsws(mm(l)) = surf_lsm_v(l)%ncsws(m)
[2232]	2003	IF ( ALLOCATED( surf_lsm_v(l)%nrsws ) ) &
	2004	surf_v(l)%nrsws(mm(l)) = surf_lsm_v(l)%nrsws(m)
	2005	IF ( ALLOCATED( surf_lsm_v(l)%sasws ) ) &
	2006	surf_v(l)%sasws(mm(l)) = surf_lsm_v(l)%sasws(m)
	2007	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_uv) ) &
	2008	surf_v(l)%mom_flux_uv(mm(l)) = surf_lsm_v(l)%mom_flux_uv(m)
	2009	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_w) ) &
	2010	surf_v(l)%mom_flux_w(mm(l)) = surf_lsm_v(l)%mom_flux_w(m)
	2011	IF ( ALLOCATED( surf_lsm_v(l)%mom_flux_tke) ) &
	2012	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_lsm_v(l)%mom_flux_tke(0:1,m)
	2013
	2014	mm(l) = mm(l) + 1
	2015	ENDDO
	2016
	2017	DO m = surf_usm_v(l)%start_index(j,i), &
	2018	surf_usm_v(l)%end_index(j,i)
	2019	IF ( ALLOCATED( surf_usm_v(l)%us ) ) &
	2020	surf_v(l)%us(mm(l)) = surf_usm_v(l)%us(m)
	2021	IF ( ALLOCATED( surf_usm_v(l)%ts ) ) &
	2022	surf_v(l)%ts(mm(l)) = surf_usm_v(l)%ts(m)
	2023	IF ( ALLOCATED( surf_usm_v(l)%qs ) ) &
	2024	surf_v(l)%qs(mm(l)) = surf_usm_v(l)%qs(m)
	2025	IF ( ALLOCATED( surf_usm_v(l)%ss ) ) &
	2026	surf_v(l)%ss(mm(l)) = surf_usm_v(l)%ss(m)
[2292]	2027	IF ( ALLOCATED( surf_usm_v(l)%qcs ) ) &
	2028	surf_v(l)%qcs(mm(l)) = surf_usm_v(l)%qcs(m)
	2029	IF ( ALLOCATED( surf_usm_v(l)%ncs ) ) &
	2030	surf_v(l)%ncs(mm(l)) = surf_usm_v(l)%ncs(m)
[2232]	2031	IF ( ALLOCATED( surf_usm_v(l)%qrs ) ) &
	2032	surf_v(l)%qrs(mm(l)) = surf_usm_v(l)%qrs(m)
	2033	IF ( ALLOCATED( surf_usm_v(l)%nrs ) ) &
	2034	surf_v(l)%nrs(mm(l)) = surf_usm_v(l)%nrs(m)
	2035	IF ( ALLOCATED( surf_usm_v(l)%ol ) ) &
	2036	surf_v(l)%ol(mm(l)) = surf_usm_v(l)%ol(m)
	2037	IF ( ALLOCATED( surf_usm_v(l)%rib ) ) &
	2038	surf_v(l)%rib(mm(l)) = surf_usm_v(l)%rib(m)
	2039	IF ( ALLOCATED( surf_usm_v(l)%usws ) ) &
	2040	surf_v(l)%usws(mm(l)) = surf_usm_v(l)%usws(m)
	2041	IF ( ALLOCATED( surf_usm_v(l)%vsws ) ) &
	2042	surf_v(l)%vsws(mm(l)) = surf_usm_v(l)%vsws(m)
	2043	IF ( ALLOCATED( surf_usm_v(l)%shf ) ) &
	2044	surf_v(l)%shf(mm(l)) = surf_usm_v(l)%shf(m)
	2045	IF ( ALLOCATED( surf_usm_v(l)%qsws ) ) &
	2046	surf_v(l)%qsws(mm(l)) = surf_usm_v(l)%qsws(m)
	2047	IF ( ALLOCATED( surf_usm_v(l)%ssws ) ) &
[2378]	2048	surf_v(l)%ssws(mm(l)) = surf_usm_v(l)%ssws(m)
[2292]	2049	IF ( ALLOCATED( surf_usm_v(l)%ncsws ) ) &
	2050	surf_v(l)%ncsws(mm(l)) = surf_usm_v(l)%ncsws(m)
[2232]	2051	IF ( ALLOCATED( surf_usm_v(l)%nrsws ) ) &
	2052	surf_v(l)%nrsws(mm(l)) = surf_usm_v(l)%nrsws(m)
	2053	IF ( ALLOCATED( surf_usm_v(l)%sasws ) ) &
	2054	surf_v(l)%sasws(mm(l)) = surf_usm_v(l)%sasws(m)
	2055	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_uv) ) &
	2056	surf_v(l)%mom_flux_uv(mm(l)) = surf_usm_v(l)%mom_flux_uv(m)
	2057	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_w) ) &
	2058	surf_v(l)%mom_flux_w(mm(l)) = surf_usm_v(l)%mom_flux_w(m)
	2059	IF ( ALLOCATED( surf_usm_v(l)%mom_flux_tke) ) &
	2060	surf_v(l)%mom_flux_tke(0:1,mm(l)) = surf_usm_v(l)%mom_flux_tke(0:1,m)
	2061
	2062	mm(l) = mm(l) + 1
	2063	ENDDO
	2064
	2065	ENDDO
	2066	ENDDO
	2067	!
	2068	!-- Finally, determine start- and end-index for the respective surface
	2069	surf_v(l)%start_index = MAX( surf_def_v(l)%start_index, &
	2070	surf_lsm_v(l)%start_index, &
	2071	surf_usm_v(l)%start_index )
	2072	surf_v(l)%end_index = MAX( surf_def_v(l)%end_index, &
	2073	surf_lsm_v(l)%end_index, &
	2074	surf_usm_v(l)%end_index )
	2075	ENDDO
	2076
	2077	WRITE ( 14 ) 'ns_h_on_file '
	2078	WRITE ( 14 ) ns_h_on_file
	2079	WRITE ( 14 ) 'ns_v_on_file '
	2080	WRITE ( 14 ) ns_v_on_file
	2081	!
	2082	!-- Write required restart data.
	2083	!-- Start with horizontal surfaces (upward-, downward-facing, and model top)
	2084	DO l = 0, 2
	2085	WRITE( dum, '(I1)') l
	2086
	2087	WRITE ( 14 ) 'surf_h(' // dum // ')%start_index '
	2088	WRITE ( 14 ) surf_h(l)%start_index
	2089	WRITE ( 14 ) 'surf_h(' // dum // ')%end_index '
	2090	WRITE ( 14 ) surf_h(l)%end_index
	2091
	2092	WRITE ( 14 ) 'surf_h(' // dum // ')%us '
	2093	IF ( ALLOCATED ( surf_h(l)%us ) ) THEN
	2094	WRITE ( 14 ) surf_h(l)%us
	2095	ENDIF
	2096	WRITE ( 14 ) 'surf_h(' // dum // ')%ts '
	2097	IF ( ALLOCATED ( surf_h(l)%ts ) ) THEN
	2098	WRITE ( 14 ) surf_h(l)%ts
	2099	ENDIF
	2100	WRITE ( 14 ) 'surf_h(' // dum // ')%qs '
	2101	IF ( ALLOCATED ( surf_h(l)%qs ) ) THEN
	2102	WRITE ( 14 ) surf_h(l)%qs
	2103	ENDIF
	2104	WRITE ( 14 ) 'surf_h(' // dum // ')%ss '
	2105	IF ( ALLOCATED ( surf_h(l)%ss ) ) THEN
	2106	WRITE ( 14 ) surf_h(l)%ss
	2107	ENDIF
[2292]	2108	WRITE ( 14 ) 'surf_h(' // dum // ')%qcs '
	2109	IF ( ALLOCATED ( surf_h(l)%qcs ) ) THEN
	2110	WRITE ( 14 ) surf_h(l)%qcs
	2111	ENDIF
	2112	WRITE ( 14 ) 'surf_h(' // dum // ')%ncs '
	2113	IF ( ALLOCATED ( surf_h(l)%ncs ) ) THEN
	2114	WRITE ( 14 ) surf_h(l)%ncs
	2115	ENDIF
[2232]	2116	WRITE ( 14 ) 'surf_h(' // dum // ')%qrs '
	2117	IF ( ALLOCATED ( surf_h(l)%qrs ) ) THEN
	2118	WRITE ( 14 ) surf_h(l)%qrs
	2119	ENDIF
	2120	WRITE ( 14 ) 'surf_h(' // dum // ')%nrs '
	2121	IF ( ALLOCATED ( surf_h(l)%nrs ) ) THEN
	2122	WRITE ( 14 ) surf_h(l)%nrs
	2123	ENDIF
	2124	WRITE ( 14 ) 'surf_h(' // dum // ')%ol '
	2125	IF ( ALLOCATED ( surf_h(l)%ol ) ) THEN
	2126	WRITE ( 14 ) surf_h(l)%ol
	2127	ENDIF
	2128	WRITE ( 14 ) 'surf_h(' // dum // ')%rib '
	2129	IF ( ALLOCATED ( surf_h(l)%rib ) ) THEN
	2130	WRITE ( 14 ) surf_h(l)%rib
	2131	ENDIF
	2132	WRITE ( 14 ) 'surf_h(' // dum // ')%usws '
	2133	IF ( ALLOCATED ( surf_h(l)%usws ) ) THEN
	2134	WRITE ( 14 ) surf_h(l)%usws
	2135	ENDIF
	2136	WRITE ( 14 ) 'surf_h(' // dum // ')%vsws '
	2137	IF ( ALLOCATED ( surf_h(l)%vsws ) ) THEN
	2138	WRITE ( 14 ) surf_h(l)%vsws
	2139	ENDIF
	2140	WRITE ( 14 ) 'surf_h(' // dum // ')%shf '
	2141	IF ( ALLOCATED ( surf_h(l)%shf ) ) THEN
	2142	WRITE ( 14 ) surf_h(l)%shf
	2143	ENDIF
	2144	WRITE ( 14 ) 'surf_h(' // dum // ')%qsws '
	2145	IF ( ALLOCATED ( surf_h(l)%qsws ) ) THEN
	2146	WRITE ( 14 ) surf_h(l)%qsws
	2147	ENDIF
	2148	WRITE ( 14 ) 'surf_h(' // dum // ')%ssws '
	2149	IF ( ALLOCATED ( surf_h(l)%ssws ) ) THEN
	2150	WRITE ( 14 ) surf_h(l)%ssws
	2151	ENDIF
[2292]	2152	WRITE ( 14 ) 'surf_h(' // dum // ')%qcsws '
	2153	IF ( ALLOCATED ( surf_h(l)%qcsws ) ) THEN
	2154	WRITE ( 14 ) surf_h(l)%qcsws
	2155	ENDIF
	2156	WRITE ( 14 ) 'surf_h(' // dum // ')%ncsws '
	2157	IF ( ALLOCATED ( surf_h(l)%ncsws ) ) THEN
	2158	WRITE ( 14 ) surf_h(l)%ncsws
	2159	ENDIF
[2232]	2160	WRITE ( 14 ) 'surf_h(' // dum // ')%qrsws '
	2161	IF ( ALLOCATED ( surf_h(l)%qrsws ) ) THEN
	2162	WRITE ( 14 ) surf_h(l)%qrsws
	2163	ENDIF
	2164	WRITE ( 14 ) 'surf_h(' // dum // ')%nrsws '
	2165	IF ( ALLOCATED ( surf_h(l)%nrsws ) ) THEN
	2166	WRITE ( 14 ) surf_h(l)%nrsws
	2167	ENDIF
	2168	WRITE ( 14 ) 'surf_h(' // dum // ')%sasws '
	2169	IF ( ALLOCATED ( surf_h(l)%sasws ) ) THEN
	2170	WRITE ( 14 ) surf_h(l)%sasws
	2171	ENDIF
	2172	ENDDO
	2173	!
	2174	!-- Write vertical surfaces
	2175	DO l = 0, 3
	2176	WRITE( dum, '(I1)') l
	2177
	2178	WRITE ( 14 ) 'surf_v(' // dum // ')%start_index '
	2179	WRITE ( 14 ) surf_v(l)%start_index
	2180	WRITE ( 14 ) 'surf_v(' // dum // ')%end_index '
	2181	WRITE ( 14 ) surf_v(l)%end_index
	2182
	2183	WRITE ( 14 ) 'surf_v(' // dum // ')%us '
	2184	IF ( ALLOCATED ( surf_v(l)%us ) ) THEN
	2185	WRITE ( 14 ) surf_v(l)%us
	2186	ENDIF
	2187	WRITE ( 14 ) 'surf_v(' // dum // ')%ts '
	2188	IF ( ALLOCATED ( surf_v(l)%ts ) ) THEN
	2189	WRITE ( 14 ) surf_v(l)%ts
	2190	ENDIF
	2191	WRITE ( 14 ) 'surf_v(' // dum // ')%qs '
	2192	IF ( ALLOCATED ( surf_v(l)%qs ) ) THEN
	2193	WRITE ( 14 ) surf_v(l)%qs
	2194	ENDIF
	2195	WRITE ( 14 ) 'surf_v(' // dum // ')%ss '
	2196	IF ( ALLOCATED ( surf_v(l)%ss ) ) THEN
	2197	WRITE ( 14 ) surf_v(l)%ss
	2198	ENDIF
[2292]	2199	WRITE ( 14 ) 'surf_v(' // dum // ')%qcs '
	2200	IF ( ALLOCATED ( surf_v(l)%qcs ) ) THEN
	2201	WRITE ( 14 ) surf_v(l)%qcs
	2202	ENDIF
	2203	WRITE ( 14 ) 'surf_v(' // dum // ')%ncs '
	2204	IF ( ALLOCATED ( surf_v(l)%ncs ) ) THEN
	2205	WRITE ( 14 ) surf_v(l)%ncs
	2206	ENDIF
[2232]	2207	WRITE ( 14 ) 'surf_v(' // dum // ')%qrs '
	2208	IF ( ALLOCATED ( surf_v(l)%qrs ) ) THEN
	2209	WRITE ( 14 ) surf_v(l)%qrs
	2210	ENDIF
	2211	WRITE ( 14 ) 'surf_v(' // dum // ')%nrs '
	2212	IF ( ALLOCATED ( surf_v(l)%nrs ) ) THEN
	2213	WRITE ( 14 ) surf_v(l)%nrs
	2214	ENDIF
	2215	WRITE ( 14 ) 'surf_v(' // dum // ')%ol '
	2216	IF ( ALLOCATED ( surf_v(l)%ol ) ) THEN
	2217	WRITE ( 14 ) surf_v(l)%ol
	2218	ENDIF
	2219	WRITE ( 14 ) 'surf_v(' // dum // ')%rib '
	2220	IF ( ALLOCATED ( surf_v(l)%rib ) ) THEN
	2221	WRITE ( 14 ) surf_v(l)%rib
	2222	ENDIF
	2223	WRITE ( 14 ) 'surf_v(' // dum // ')%shf '
	2224	IF ( ALLOCATED ( surf_v(l)%shf ) ) THEN
	2225	WRITE ( 14 ) surf_v(l)%shf
	2226	ENDIF
	2227	WRITE ( 14 ) 'surf_v(' // dum // ')%qsws '
	2228	IF ( ALLOCATED ( surf_v(l)%qsws ) ) THEN
	2229	WRITE ( 14 ) surf_v(l)%qsws
	2230	ENDIF
	2231	WRITE ( 14 ) 'surf_v(' // dum // ')%ssws '
	2232	IF ( ALLOCATED ( surf_v(l)%ssws ) ) THEN
	2233	WRITE ( 14 ) surf_v(l)%ssws
	2234	ENDIF
[2292]	2235	WRITE ( 14 ) 'surf_v(' // dum // ')%qcsws '
	2236	IF ( ALLOCATED ( surf_v(l)%qcsws ) ) THEN
	2237	WRITE ( 14 ) surf_v(l)%qcsws
	2238	ENDIF
	2239	WRITE ( 14 ) 'surf_v(' // dum // ')%ncsws '
	2240	IF ( ALLOCATED ( surf_v(l)%ncsws ) ) THEN
	2241	WRITE ( 14 ) surf_v(l)%ncsws
	2242	ENDIF
[2232]	2243	WRITE ( 14 ) 'surf_v(' // dum // ')%qrsws '
	2244	IF ( ALLOCATED ( surf_v(l)%qrsws ) ) THEN
	2245	WRITE ( 14 ) surf_v(l)%qrsws
	2246	ENDIF
	2247	WRITE ( 14 ) 'surf_v(' // dum // ')%nrsws '
	2248	IF ( ALLOCATED ( surf_v(l)%nrsws ) ) THEN
	2249	WRITE ( 14 ) surf_v(l)%nrsws
	2250	ENDIF
	2251	WRITE ( 14 ) 'surf_v(' // dum // ')%sasws '
	2252	IF ( ALLOCATED ( surf_v(l)%sasws ) ) THEN
	2253	WRITE ( 14 ) surf_v(l)%sasws
	2254	ENDIF
	2255	WRITE ( 14 ) 'surf_v(' // dum // ')%mom_uv '
	2256	IF ( ALLOCATED ( surf_v(l)%mom_flux_uv ) ) THEN
	2257	WRITE ( 14 ) surf_v(l)%mom_flux_uv
	2258	ENDIF
	2259	WRITE ( 14 ) 'surf_v(' // dum // ')%mom_w '
	2260	IF ( ALLOCATED ( surf_v(l)%mom_flux_w ) ) THEN
	2261	WRITE ( 14 ) surf_v(l)%mom_flux_w
	2262	ENDIF
	2263	WRITE ( 14 ) 'surf_v(' // dum // ')%mom_tke '
	2264	IF ( ALLOCATED ( surf_v(l)%mom_flux_tke ) ) THEN
	2265	WRITE ( 14 ) surf_v(l)%mom_flux_tke
	2266	ENDIF
	2267
	2268	ENDDO
	2269
	2270	WRITE ( 14 ) '* end surf * '
	2271
	2272	END SUBROUTINE surface_write_restart_data
	2273
	2274
	2275	!------------------------------------------------------------------------------!
	2276	! Description:
	2277	! ------------
	2278	!> Reads surface-related restart data. Please note, restart data for a certain
	2279	!> surface orientation (e.g. horizontal upward-facing) is stored in one
	2280	!> array, even if surface elements may belong to different surface types
	2281	!> natural or urban for example). Surface elements are redistributed into its
	2282	!> respective surface types within this routine. This allows e.g. changing the
	2283	!> surface type after reading the restart data, which might be required in case
[2269]	2284	!> of cyclic_fill mode.
[2232]	2285	!------------------------------------------------------------------------------!
	2286	SUBROUTINE surface_read_restart_data( ii, &
	2287	nxlfa, nxl_on_file, nxrfa, nxr_on_file, &
	2288	nynfa, nyn_on_file, nysfa, nys_on_file, &
	2289	offset_xa, offset_ya, overlap_count )
	2290
	2291	USE pegrid, &
	2292	ONLY: numprocs_previous_run
	2293
	2294	CHARACTER (LEN=1) :: dum !< dummy to create correct string for reading input variable
	2295	CHARACTER (LEN=30) :: field_chr !< input variable
	2296
	2297	INTEGER(iwp) :: i !< running index along x-direction, refers to former domain size
	2298	INTEGER(iwp) :: ic !< running index along x-direction, refers to current domain size
	2299	INTEGER(iwp) :: j !< running index along y-direction, refers to former domain size
	2300	INTEGER(iwp) :: jc !< running index along y-direction, refers to former domain size
	2301	INTEGER(iwp) :: k !< running index along z-direction
	2302	INTEGER(iwp) :: l !< index variable for surface type
	2303	INTEGER(iwp) :: m !< running index for surface elements, refers to gathered array encompassing all surface types
	2304	INTEGER(iwp) :: mm !< running index for surface elements, refers to individual surface types
	2305
	2306	INTEGER(iwp) :: ii !< running index over input files
	2307	INTEGER(iwp) :: kk !< running index over previous input files covering current local domain
	2308	INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain
	2309	INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain
	2310	INTEGER(iwp) :: nxl_on_file !< index of left boundary on former local domain
	2311	INTEGER(iwp) :: nxrc !< index of right boundary on current subdomain
	2312	INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain
	2313	INTEGER(iwp) :: nxr_on_file !< index of right boundary on former local domain
	2314	INTEGER(iwp) :: nync !< index of north boundary on current subdomain
	2315	INTEGER(iwp) :: nynf !< index of north boundary on former subdomain
	2316	INTEGER(iwp) :: nyn_on_file !< index of norht boundary on former local domain
	2317	INTEGER(iwp) :: nysc !< index of south boundary on current subdomain
	2318	INTEGER(iwp) :: nysf !< index of south boundary on former subdomain
	2319	INTEGER(iwp) :: nys_on_file !< index of south boundary on former local domain
	2320	INTEGER(iwp) :: overlap_count !< number of overlaps
	2321
	2322	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nxlfa !<
	2323	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nxrfa !<
	2324	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nynfa !<
	2325	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nysfa !<
	2326	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: offset_xa !<
	2327	INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: offset_ya !<
	2328
	2329
	2330	LOGICAL :: horizontal_surface !< flag indicating horizontal surfaces
	2331	LOGICAL :: surf_match_def !< flag indicating that surface element is of default type
	2332	LOGICAL :: surf_match_lsm !< flag indicating that surface element is of natural type
	2333	LOGICAL :: surf_match_usm !< flag indicating that surface element is of urban type
[2269]	2334	LOGICAL :: vertical_surface !< flag indicating vertical surfaces
[2232]	2335
	2336	TYPE(surf_type), DIMENSION(0:2) :: surf_h !< horizontal surface type on file
	2337	TYPE(surf_type), DIMENSION(0:3) :: surf_v !< vertical surface type on file
	2338
	2339	!
	2340	!-- Read number of respective surface elements on file
	2341	READ ( 13 ) field_chr
	2342	IF ( TRIM( field_chr ) /= 'ns_h_on_file' ) THEN
	2343	!
	2344	!-- Add a proper error message
	2345	ENDIF
	2346	READ ( 13 ) ns_h_on_file
	2347
	2348	READ ( 13 ) field_chr
	2349	IF ( TRIM( field_chr ) /= 'ns_v_on_file' ) THEN
	2350	!
	2351	!-- Add a proper error message
	2352	ENDIF
	2353	READ ( 13 ) ns_v_on_file
	2354	!
	2355	!-- Allocate memory for number of surface elements on file. Please note,
	2356	!-- these number is not necessarily the same as the final number of surface
	2357	!-- elements on local domain, which is the case if processor topology changes
	2358	!-- during restart runs.
	2359	!-- Horizontal upward facing
	2360	surf_h(0)%ns = ns_h_on_file(0)
	2361	CALL allocate_surface_attributes_h( surf_h(0), &
	2362	nys_on_file, nyn_on_file, &
	2363	nxl_on_file, nxr_on_file )
	2364	!
	2365	!-- Horizontal downward facing
	2366	surf_h(1)%ns = ns_h_on_file(1)
	2367	CALL allocate_surface_attributes_h( surf_h(1), &
	2368	nys_on_file, nyn_on_file, &
	2369	nxl_on_file, nxr_on_file )
	2370	!
	2371	!-- Model top
	2372	surf_h(2)%ns = ns_h_on_file(2)
	2373	CALL allocate_surface_attributes_h_top( surf_h(2), &
	2374	nys_on_file, nyn_on_file, &
	2375	nxl_on_file, nxr_on_file )
	2376	!
	2377	!-- Vertical surfaces
	2378	DO l = 0, 3
	2379	surf_v(l)%ns = ns_v_on_file(l)
	2380	CALL allocate_surface_attributes_v( surf_v(l), .FALSE., &
	2381	nys_on_file, nyn_on_file, &
	2382	nxl_on_file, nxr_on_file )
	2383	ENDDO
	2384
	2385	IF ( initializing_actions == 'read_restart_data' .OR. &
	2386	initializing_actions == 'cyclic_fill' ) THEN
	2387	!
	2388	!-- Initial setting of flags for horizontal and vertical surfaces, will
	2389	!-- be set after start- and end-indices are read.
	2390	horizontal_surface = .FALSE.
	2391	vertical_surface = .FALSE.
	2392
	2393	READ ( 13 ) field_chr
	2394
	2395	DO WHILE ( TRIM( field_chr ) /= '* end surf *' )
	2396	!
	2397	!-- Map data on file as often as needed (data are read only for k=1)
	2398	DO kk = 1, overlap_count
	2399	!
	2400	!-- Get the index range of the subdomain on file which overlap with the
	2401	!-- current subdomain
	2402	nxlf = nxlfa(ii,kk)
	2403	nxlc = nxlfa(ii,kk) + offset_xa(ii,kk)
	2404	nxrf = nxrfa(ii,kk)
	2405	nxrc = nxrfa(ii,kk) + offset_xa(ii,kk)
	2406	nysf = nysfa(ii,kk)
	2407	nysc = nysfa(ii,kk) + offset_ya(ii,kk)
	2408	nynf = nynfa(ii,kk)
	2409	nync = nynfa(ii,kk) + offset_ya(ii,kk)
	2410
	2411	SELECT CASE ( TRIM( field_chr ) )
	2412
	2413	CASE ( 'surf_h(0)%start_index' )
	2414	IF ( kk == 1 ) &
	2415	READ ( 13 ) surf_h(0)%start_index
	2416	l = 0
	2417	CASE ( 'surf_h(0)%end_index' )
	2418	IF ( kk == 1 ) &
	2419	READ ( 13 ) surf_h(0)%end_index
	2420	horizontal_surface = .TRUE.
	2421	vertical_surface = .FALSE.
	2422	CASE ( 'surf_h(0)%us' )
	2423	IF ( ALLOCATED( surf_h(0)%us ) .AND. kk == 1 ) &
	2424	READ ( 13 ) surf_h(0)%us
	2425	CASE ( 'surf_h(0)%ts' )
	2426	IF ( ALLOCATED( surf_h(0)%ts ) .AND. kk == 1 ) &
	2427	READ ( 13 ) surf_h(0)%ts
	2428	CASE ( 'surf_h(0)%qs' )
	2429	IF ( ALLOCATED( surf_h(0)%qs ) .AND. kk == 1 ) &
	2430	READ ( 13 ) surf_h(0)%qs
	2431	CASE ( 'surf_h(0)%ss' )
	2432	IF ( ALLOCATED( surf_h(0)%ss ) .AND. kk == 1 ) &
	2433	READ ( 13 ) surf_h(0)%ss
[2292]	2434	CASE ( 'surf_h(0)%qcs' )
	2435	IF ( ALLOCATED( surf_h(0)%qcs ) .AND. kk == 1 ) &
	2436	READ ( 13 ) surf_h(0)%qcs
	2437	CASE ( 'surf_h(0)%ncs' )
	2438	IF ( ALLOCATED( surf_h(0)%ncs ) .AND. kk == 1 ) &
	2439	READ ( 13 ) surf_h(0)%ncs
[2232]	2440	CASE ( 'surf_h(0)%qrs' )
	2441	IF ( ALLOCATED( surf_h(0)%qrs ) .AND. kk == 1 ) &
	2442	READ ( 13 ) surf_h(0)%qrs
	2443	CASE ( 'surf_h(0)%nrs' )
	2444	IF ( ALLOCATED( surf_h(0)%nrs ) .AND. kk == 1 ) &
	2445	READ ( 13 ) surf_h(0)%nrs
	2446	CASE ( 'surf_h(0)%ol' )
	2447	IF ( ALLOCATED( surf_h(0)%ol ) .AND. kk == 1 ) &
	2448	READ ( 13 ) surf_h(0)%ol
	2449	CASE ( 'surf_h(0)%rib' )
	2450	IF ( ALLOCATED( surf_h(0)%rib ) .AND. kk == 1 ) &
	2451	READ ( 13 ) surf_h(0)%rib
	2452	CASE ( 'surf_h(0)%usws' )
	2453	IF ( ALLOCATED( surf_h(0)%usws ) .AND. kk == 1 ) &
	2454	READ ( 13 ) surf_h(0)%usws
	2455	CASE ( 'surf_h(0)%vsws' )
	2456	IF ( ALLOCATED( surf_h(0)%vsws ) .AND. kk == 1 ) &
	2457	READ ( 13 ) surf_h(0)%vsws
	2458	CASE ( 'surf_h(0)%shf' )
	2459	IF ( ALLOCATED( surf_h(0)%shf ) .AND. kk == 1 ) &
	2460	READ ( 13 ) surf_h(0)%shf
	2461	CASE ( 'surf_h(0)%qsws' )
	2462	IF ( ALLOCATED( surf_h(0)%qsws ) .AND. kk == 1 ) &
	2463	READ ( 13 ) surf_h(0)%qsws
	2464	CASE ( 'surf_h(0)%ssws' )
	2465	IF ( ALLOCATED( surf_h(0)%ssws ) .AND. kk == 1 ) &
	2466	READ ( 13 ) surf_h(0)%ssws
[2292]	2467	CASE ( 'surf_h(0)%qcsws' )
	2468	IF ( ALLOCATED( surf_h(0)%qcsws ) .AND. kk == 1 ) &
	2469	READ ( 13 ) surf_h(0)%qcsws
	2470	CASE ( 'surf_h(0)%ncsws' )
	2471	IF ( ALLOCATED( surf_h(0)%ncsws ) .AND. kk == 1 ) &
	2472	READ ( 13 ) surf_h(0)%ncsws
[2232]	2473	CASE ( 'surf_h(0)%qrsws' )
	2474	IF ( ALLOCATED( surf_h(0)%qrsws ) .AND. kk == 1 ) &
	2475	READ ( 13 ) surf_h(0)%qrsws
	2476	CASE ( 'surf_h(0)%nrsws' )
	2477	IF ( ALLOCATED( surf_h(0)%nrsws ) .AND. kk == 1 ) &
	2478	READ ( 13 ) surf_h(0)%nrsws
	2479	CASE ( 'surf_h(0)%sasws' )
	2480	IF ( ALLOCATED( surf_h(0)%sasws ) .AND. kk == 1 ) &
	2481	READ ( 13 ) surf_h(0)%sasws
	2482
	2483	CASE ( 'surf_h(1)%start_index' )
	2484	IF ( kk == 1 ) &
	2485	READ ( 13 ) surf_h(1)%start_index
	2486	l = 1
	2487	CASE ( 'surf_h(1)%end_index' )
	2488	IF ( kk == 1 ) &
	2489	READ ( 13 ) surf_h(1)%end_index
	2490	CASE ( 'surf_h(1)%us' )
	2491	IF ( ALLOCATED( surf_h(1)%us ) .AND. kk == 1 ) &
	2492	READ ( 13 ) surf_h(1)%us
	2493	CASE ( 'surf_h(1)%ts' )
	2494	IF ( ALLOCATED( surf_h(1)%ts ) .AND. kk == 1 ) &
	2495	READ ( 13 ) surf_h(1)%ts
	2496	CASE ( 'surf_h(1)%qs' )
	2497	IF ( ALLOCATED( surf_h(1)%qs ) .AND. kk == 1 ) &
	2498	READ ( 13 ) surf_h(1)%qs
	2499	CASE ( 'surf_h(1)%ss' )
	2500	IF ( ALLOCATED( surf_h(1)%ss ) .AND. kk == 1 ) &
	2501	READ ( 13 ) surf_h(1)%ss
[2292]	2502	CASE ( 'surf_h(1)%qcs' )
	2503	IF ( ALLOCATED( surf_h(1)%qcs ) .AND. kk == 1 ) &
	2504	READ ( 13 ) surf_h(1)%qcs
	2505	CASE ( 'surf_h(1)%ncs' )
	2506	IF ( ALLOCATED( surf_h(1)%ncs ) .AND. kk == 1 ) &
	2507	READ ( 13 ) surf_h(1)%ncs
[2232]	2508	CASE ( 'surf_h(1)%qrs' )
	2509	IF ( ALLOCATED( surf_h(1)%qrs ) .AND. kk == 1 ) &
	2510	READ ( 13 ) surf_h(1)%qrs
	2511	CASE ( 'surf_h(1)%nrs' )
	2512	IF ( ALLOCATED( surf_h(1)%nrs ) .AND. kk == 1 ) &
	2513	READ ( 13 ) surf_h(1)%nrs
	2514	CASE ( 'surf_h(1)%ol' )
	2515	IF ( ALLOCATED( surf_h(1)%ol ) .AND. kk == 1 ) &
	2516	READ ( 13 ) surf_h(1)%ol
	2517	CASE ( 'surf_h(1)%rib' )
	2518	IF ( ALLOCATED( surf_h(1)%rib ) .AND. kk == 1 ) &
	2519	READ ( 13 ) surf_h(1)%rib
	2520	CASE ( 'surf_h(1)%usws' )
	2521	IF ( ALLOCATED( surf_h(1)%usws ) .AND. kk == 1 ) &
	2522	READ ( 13 ) surf_h(1)%usws
	2523	CASE ( 'surf_h(1)%vsws' )
	2524	IF ( ALLOCATED( surf_h(1)%vsws ) .AND. kk == 1 ) &
	2525	READ ( 13 ) surf_h(1)%vsws
	2526	CASE ( 'surf_h(1)%shf' )
	2527	IF ( ALLOCATED( surf_h(1)%shf ) .AND. kk == 1 ) &
	2528	READ ( 13 ) surf_h(1)%shf
	2529	CASE ( 'surf_h(1)%qsws' )
	2530	IF ( ALLOCATED( surf_h(1)%qsws ) .AND. kk == 1 ) &
	2531	READ ( 13 ) surf_h(1)%qsws
	2532	CASE ( 'surf_h(1)%ssws' )
	2533	IF ( ALLOCATED( surf_h(1)%ssws ) .AND. kk == 1 ) &
	2534	READ ( 13 ) surf_h(1)%ssws
[2292]	2535	CASE ( 'surf_h(1)%qcsws' )
	2536	IF ( ALLOCATED( surf_h(1)%qcsws ) .AND. kk == 1 ) &
	2537	READ ( 13 ) surf_h(1)%qcsws
	2538	CASE ( 'surf_h(1)%ncsws' )
	2539	IF ( ALLOCATED( surf_h(1)%ncsws ) .AND. kk == 1 ) &
	2540	READ ( 13 ) surf_h(1)%ncsws
[2232]	2541	CASE ( 'surf_h(1)%qrsws' )
	2542	IF ( ALLOCATED( surf_h(1)%qrsws ) .AND. kk == 1 ) &
	2543	READ ( 13 ) surf_h(1)%qrsws
	2544	CASE ( 'surf_h(1)%nrsws' )
	2545	IF ( ALLOCATED( surf_h(1)%nrsws ) .AND. kk == 1 ) &
	2546	READ ( 13 ) surf_h(1)%nrsws
	2547	CASE ( 'surf_h(1)%sasws' )
	2548	IF ( ALLOCATED( surf_h(1)%sasws ) .AND. kk == 1 ) &
	2549	READ ( 13 ) surf_h(1)%sasws
	2550
	2551	CASE ( 'surf_h(2)%start_index' )
	2552	IF ( kk == 1 ) &
	2553	READ ( 13 ) surf_h(2)%start_index
	2554	l = 2
	2555	CASE ( 'surf_h(2)%end_index' )
	2556	IF ( kk == 1 ) &
	2557	READ ( 13 ) surf_h(2)%end_index
	2558	CASE ( 'surf_h(2)%us' )
	2559	IF ( ALLOCATED( surf_h(2)%us ) .AND. kk == 1 ) &
	2560	READ ( 13 ) surf_h(2)%us
	2561	CASE ( 'surf_h(2)%ts' )
	2562	IF ( ALLOCATED( surf_h(2)%ts ) .AND. kk == 1 ) &
	2563	READ ( 13 ) surf_h(2)%ts
	2564	CASE ( 'surf_h(2)%qs' )
	2565	IF ( ALLOCATED( surf_h(2)%qs ) .AND. kk == 1 ) &
	2566	READ ( 13 ) surf_h(2)%qs
	2567	CASE ( 'surf_h(2)%ss' )
	2568	IF ( ALLOCATED( surf_h(2)%ss ) .AND. kk == 1 ) &
	2569	READ ( 13 ) surf_h(2)%ss
[2292]	2570	CASE ( 'surf_h(2)%qcs' )
	2571	IF ( ALLOCATED( surf_h(2)%qcs ) .AND. kk == 1 ) &
	2572	READ ( 13 ) surf_h(2)%qcs
	2573	CASE ( 'surf_h(2)%ncs' )
	2574	IF ( ALLOCATED( surf_h(2)%ncs ) .AND. kk == 1 ) &
	2575	READ ( 13 ) surf_h(2)%ncs
[2232]	2576	CASE ( 'surf_h(2)%qrs' )
	2577	IF ( ALLOCATED( surf_h(2)%qrs ) .AND. kk == 1 ) &
	2578	READ ( 13 ) surf_h(2)%qrs
	2579	CASE ( 'surf_h(2)%nrs' )
	2580	IF ( ALLOCATED( surf_h(2)%nrs ) .AND. kk == 1 ) &
	2581	READ ( 13 ) surf_h(2)%nrs
	2582	CASE ( 'surf_h(2)%ol' )
	2583	IF ( ALLOCATED( surf_h(2)%ol ) .AND. kk == 1 ) &
	2584	READ ( 13 ) surf_h(2)%ol
	2585	CASE ( 'surf_h(2)%rib' )
	2586	IF ( ALLOCATED( surf_h(2)%rib ) .AND. kk == 1 ) &
	2587	READ ( 13 ) surf_h(2)%rib
	2588	CASE ( 'surf_h(2)%usws' )
	2589	IF ( ALLOCATED( surf_h(2)%usws ) .AND. kk == 1 ) &
	2590	READ ( 13 ) surf_h(2)%usws
	2591	CASE ( 'surf_h(2)%vsws' )
	2592	IF ( ALLOCATED( surf_h(2)%vsws ) .AND. kk == 1 ) &
	2593	READ ( 13 ) surf_h(2)%vsws
	2594	CASE ( 'surf_h(2)%shf' )
	2595	IF ( ALLOCATED( surf_h(2)%shf ) .AND. kk == 1 ) &
	2596	READ ( 13 ) surf_h(2)%shf
	2597	CASE ( 'surf_h(2)%qsws' )
	2598	IF ( ALLOCATED( surf_h(2)%qsws ) .AND. kk == 1 ) &
	2599	READ ( 13 ) surf_h(2)%qsws
	2600	CASE ( 'surf_h(2)%ssws' )
	2601	IF ( ALLOCATED( surf_h(2)%ssws ) .AND. kk == 1 ) &
	2602	READ ( 13 ) surf_h(2)%ssws
[2292]	2603	CASE ( 'surf_h(2)%qcsws' )
	2604	IF ( ALLOCATED( surf_h(2)%qcsws ) .AND. kk == 1 ) &
	2605	READ ( 13 ) surf_h(2)%qcsws
	2606	CASE ( 'surf_h(2)%ncsws' )
	2607	IF ( ALLOCATED( surf_h(2)%ncsws ) .AND. kk == 1 ) &
	2608	READ ( 13 ) surf_h(2)%ncsws
[2232]	2609	CASE ( 'surf_h(2)%qrsws' )
	2610	IF ( ALLOCATED( surf_h(2)%qrsws ) .AND. kk == 1 ) &
	2611	READ ( 13 ) surf_h(2)%qrsws
	2612	CASE ( 'surf_h(2)%nrsws' )
	2613	IF ( ALLOCATED( surf_h(2)%nrsws ) .AND. kk == 1 ) &
	2614	READ ( 13 ) surf_h(2)%nrsws
	2615	CASE ( 'surf_h(2)%sasws' )
	2616	IF ( ALLOCATED( surf_h(2)%sasws ) .AND. kk == 1 ) &
	2617	READ ( 13 ) surf_h(2)%sasws
	2618
	2619	CASE ( 'surf_v(0)%start_index' )
	2620	IF ( kk == 1 ) &
	2621	READ ( 13 ) surf_v(0)%start_index
	2622	l = 0
	2623	horizontal_surface = .FALSE.
	2624	vertical_surface = .TRUE.
	2625	CASE ( 'surf_v(0)%end_index' )
	2626	IF ( kk == 1 ) &
	2627	READ ( 13 ) surf_v(0)%end_index
	2628	CASE ( 'surf_v(0)%us' )
	2629	IF ( ALLOCATED( surf_v(0)%us ) .AND. kk == 1 ) &
	2630	READ ( 13 ) surf_v(0)%us
	2631	CASE ( 'surf_v(0)%ts' )
	2632	IF ( ALLOCATED( surf_v(0)%ts ) .AND. kk == 1 ) &
	2633	READ ( 13 ) surf_v(0)%ts
	2634	CASE ( 'surf_v(0)%qs' )
	2635	IF ( ALLOCATED( surf_v(0)%qs ) .AND. kk == 1 ) &
	2636	READ ( 13 ) surf_v(0)%qs
	2637	CASE ( 'surf_v(0)%ss' )
	2638	IF ( ALLOCATED( surf_v(0)%ss ) .AND. kk == 1 ) &
	2639	READ ( 13 ) surf_v(0)%ss
[2292]	2640	CASE ( 'surf_v(0)%qcs' )
	2641	IF ( ALLOCATED( surf_v(0)%qcs ) .AND. kk == 1 ) &
	2642	READ ( 13 ) surf_v(0)%qcs
	2643	CASE ( 'surf_v(0)%ncs' )
	2644	IF ( ALLOCATED( surf_v(0)%ncs ) .AND. kk == 1 ) &
	2645	READ ( 13 ) surf_v(0)%ncs
[2232]	2646	CASE ( 'surf_v(0)%qrs' )
	2647	IF ( ALLOCATED( surf_v(0)%qrs ) .AND. kk == 1 ) &
	2648	READ ( 13 ) surf_v(0)%qrs
	2649	CASE ( 'surf_v(0)%nrs' )
	2650	IF ( ALLOCATED( surf_v(0)%nrs ) .AND. kk == 1 ) &
	2651	READ ( 13 ) surf_v(0)%nrs
	2652	CASE ( 'surf_v(0)%ol' )
	2653	IF ( ALLOCATED( surf_v(0)%ol ) .AND. kk == 1 ) &
	2654	READ ( 13 ) surf_v(0)%ol
	2655	CASE ( 'surf_v(0)%rib' )
	2656	IF ( ALLOCATED( surf_v(0)%rib ) .AND. kk == 1 ) &
	2657	READ ( 13 ) surf_v(0)%rib
	2658	CASE ( 'surf_v(0)%shf' )
	2659	IF ( ALLOCATED( surf_v(0)%shf ) .AND. kk == 1 ) &
	2660	READ ( 13 ) surf_v(0)%shf
	2661	CASE ( 'surf_v(0)%qsws' )
	2662	IF ( ALLOCATED( surf_v(0)%qsws ) .AND. kk == 1 ) &
	2663	READ ( 13 ) surf_v(0)%qsws
	2664	CASE ( 'surf_v(0)%ssws' )
	2665	IF ( ALLOCATED( surf_v(0)%ssws ) .AND. kk == 1 ) &
	2666	READ ( 13 ) surf_v(0)%ssws
[2292]	2667	CASE ( 'surf_v(0)%qcsws' )
	2668	IF ( ALLOCATED( surf_v(0)%qcsws ) .AND. kk == 1 ) &
	2669	READ ( 13 ) surf_v(0)%qcsws
	2670	CASE ( 'surf_v(0)%ncsws' )
	2671	IF ( ALLOCATED( surf_v(0)%ncsws ) .AND. kk == 1 ) &
	2672	READ ( 13 ) surf_v(0)%ncsws
[2232]	2673	CASE ( 'surf_v(0)%qrsws' )
	2674	IF ( ALLOCATED( surf_v(0)%qrsws ) .AND. kk == 1 ) &
	2675	READ ( 13 ) surf_v(0)%qrsws
	2676	CASE ( 'surf_v(0)%nrsws' )
	2677	IF ( ALLOCATED( surf_v(0)%nrsws ) .AND. kk == 1 ) &
	2678	READ ( 13 ) surf_v(0)%nrsws
	2679	CASE ( 'surf_v(0)%sasws' )
	2680	IF ( ALLOCATED( surf_v(0)%sasws ) .AND. kk == 1 ) &
	2681	READ ( 13 ) surf_v(0)%sasws
	2682	CASE ( 'surf_v(0)%mom_uv' )
	2683	IF ( ALLOCATED( surf_v(0)%mom_flux_uv ) .AND. kk == 1 )&
	2684	READ ( 13 ) surf_v(0)%mom_flux_uv
	2685	CASE ( 'surf_v(0)%mom_w' )
	2686	IF ( ALLOCATED( surf_v(0)%mom_flux_w ) .AND. kk == 1 ) &
	2687	READ ( 13 ) surf_v(0)%mom_flux_w
	2688	CASE ( 'surf_v(0)%mom_tke' )
	2689	IF ( ALLOCATED( surf_v(0)%mom_flux_tke ) .AND. kk == 1 )&
	2690	READ ( 13 ) surf_v(0)%mom_flux_tke
	2691
	2692	CASE ( 'surf_v(1)%start_index' )
	2693	IF ( kk == 1 ) &
	2694	READ ( 13 ) surf_v(1)%start_index
	2695	l = 1
	2696	CASE ( 'surf_v(1)%end_index' )
	2697	IF ( kk == 1 ) &
	2698	READ ( 13 ) surf_v(1)%end_index
	2699	CASE ( 'surf_v(1)%us' )
	2700	IF ( ALLOCATED( surf_v(1)%us ) .AND. kk == 1 ) &
	2701	READ ( 13 ) surf_v(1)%us
	2702	CASE ( 'surf_v(1)%ts' )
	2703	IF ( ALLOCATED( surf_v(1)%ts ) .AND. kk == 1 ) &
	2704	READ ( 13 ) surf_v(1)%ts
	2705	CASE ( 'surf_v(1)%qs' )
	2706	IF ( ALLOCATED( surf_v(1)%qs ) .AND. kk == 1 ) &
	2707	READ ( 13 ) surf_v(1)%qs
	2708	CASE ( 'surf_v(1)%ss' )
	2709	IF ( ALLOCATED( surf_v(1)%ss ) .AND. kk == 1 ) &
	2710	READ ( 13 ) surf_v(1)%ss
[2292]	2711	CASE ( 'surf_v(1)%qcs' )
	2712	IF ( ALLOCATED( surf_v(1)%qcs ) .AND. kk == 1 ) &
	2713	READ ( 13 ) surf_v(1)%qcs
	2714	CASE ( 'surf_v(1)%ncs' )
	2715	IF ( ALLOCATED( surf_v(1)%ncs ) .AND. kk == 1 ) &
	2716	READ ( 13 ) surf_v(1)%ncs
[2232]	2717	CASE ( 'surf_v(1)%qrs' )
	2718	IF ( ALLOCATED( surf_v(1)%qrs ) .AND. kk == 1 ) &
	2719	READ ( 13 ) surf_v(1)%qrs
	2720	CASE ( 'surf_v(1)%nrs' )
	2721	IF ( ALLOCATED( surf_v(1)%nrs ) .AND. kk == 1 ) &
	2722	READ ( 13 ) surf_v(1)%nrs
	2723	CASE ( 'surf_v(1)%ol' )
	2724	IF ( ALLOCATED( surf_v(1)%ol ) .AND. kk == 1 ) &
	2725	READ ( 13 ) surf_v(1)%ol
	2726	CASE ( 'surf_v(1)%rib' )
	2727	IF ( ALLOCATED( surf_v(1)%rib ) .AND. kk == 1 ) &
	2728	READ ( 13 ) surf_v(1)%rib
	2729	CASE ( 'surf_v(1)%shf' )
	2730	IF ( ALLOCATED( surf_v(1)%shf ) .AND. kk == 1 ) &
	2731	READ ( 13 ) surf_v(1)%shf
	2732	CASE ( 'surf_v(1)%qsws' )
	2733	IF ( ALLOCATED( surf_v(1)%qsws ) .AND. kk == 1 ) &
	2734	READ ( 13 ) surf_v(1)%qsws
	2735	CASE ( 'surf_v(1)%ssws' )
	2736	IF ( ALLOCATED( surf_v(1)%ssws ) .AND. kk == 1 ) &
	2737	READ ( 13 ) surf_v(1)%ssws
[2292]	2738	CASE ( 'surf_v(1)%qcsws' )
	2739	IF ( ALLOCATED( surf_v(1)%qcsws ) .AND. kk == 1 ) &
	2740	READ ( 13 ) surf_v(1)%qcsws
	2741	CASE ( 'surf_v(1)%ncsws' )
	2742	IF ( ALLOCATED( surf_v(1)%ncsws ) .AND. kk == 1 ) &
	2743	READ ( 13 ) surf_v(1)%ncsws
[2232]	2744	CASE ( 'surf_v(1)%qrsws' )
	2745	IF ( ALLOCATED( surf_v(1)%qrsws ) .AND. kk == 1 ) &
	2746	READ ( 13 ) surf_v(1)%qrsws
	2747	CASE ( 'surf_v(1)%nrsws' )
	2748	IF ( ALLOCATED( surf_v(1)%nrsws ) .AND. kk == 1 ) &
	2749	READ ( 13 ) surf_v(1)%nrsws
	2750	CASE ( 'surf_v(1)%sasws' )
	2751	IF ( ALLOCATED( surf_v(1)%sasws ) .AND. kk == 1 ) &
	2752	READ ( 13 ) surf_v(1)%sasws
	2753	CASE ( 'surf_v(1)%mom_uv' )
	2754	IF ( ALLOCATED( surf_v(1)%mom_flux_uv ) .AND. kk == 1 )&
	2755	READ ( 13 ) surf_v(1)%mom_flux_uv
	2756	CASE ( 'surf_v(1)%mom_w' )
	2757	IF ( ALLOCATED( surf_v(1)%mom_flux_w ) .AND. kk == 1 ) &
	2758	READ ( 13 ) surf_v(1)%mom_flux_w
	2759	CASE ( 'surf_v(1)%mom_tke' )
	2760	IF ( ALLOCATED( surf_v(1)%mom_flux_tke ) .AND. kk == 1 )&
	2761	READ ( 13 ) surf_v(1)%mom_flux_tke
	2762
	2763	CASE ( 'surf_v(2)%start_index' )
	2764	IF ( kk == 1 ) &
	2765	READ ( 13 ) surf_v(2)%start_index
	2766	l = 2
	2767	CASE ( 'surf_v(2)%end_index' )
	2768	IF ( kk == 1 ) &
	2769	READ ( 13 ) surf_v(2)%end_index
	2770	CASE ( 'surf_v(2)%us' )
	2771	IF ( ALLOCATED( surf_v(2)%us ) .AND. kk == 1 ) &
	2772	READ ( 13 ) surf_v(2)%us
	2773	CASE ( 'surf_v(2)%ts' )
	2774	IF ( ALLOCATED( surf_v(2)%ts ) .AND. kk == 1 ) &
	2775	READ ( 13 ) surf_v(2)%ts
	2776	CASE ( 'surf_v(2)%qs' )
	2777	IF ( ALLOCATED( surf_v(2)%qs ) .AND. kk == 1 ) &
	2778	READ ( 13 ) surf_v(2)%qs
	2779	CASE ( 'surf_v(2)%ss' )
	2780	IF ( ALLOCATED( surf_v(2)%ss ) .AND. kk == 1 ) &
	2781	READ ( 13 ) surf_v(2)%ss
[2292]	2782	CASE ( 'surf_v(2)%qcs' )
	2783	IF ( ALLOCATED( surf_v(2)%qcs ) .AND. kk == 1 ) &
	2784	READ ( 13 ) surf_v(2)%qcs
	2785	CASE ( 'surf_v(2)%ncs' )
	2786	IF ( ALLOCATED( surf_v(2)%ncs ) .AND. kk == 1 ) &
	2787	READ ( 13 ) surf_v(2)%ncs
[2232]	2788	CASE ( 'surf_v(2)%qrs' )
	2789	IF ( ALLOCATED( surf_v(2)%qrs ) .AND. kk == 1 ) &
	2790	READ ( 13 ) surf_v(2)%qrs
	2791	CASE ( 'surf_v(2)%nrs' )
	2792	IF ( ALLOCATED( surf_v(2)%nrs ) .AND. kk == 1 ) &
	2793	READ ( 13 ) surf_v(2)%nrs
	2794	CASE ( 'surf_v(2)%ol' )
	2795	IF ( ALLOCATED( surf_v(2)%ol ) .AND. kk == 1 ) &
	2796	READ ( 13 ) surf_v(2)%ol
	2797	CASE ( 'surf_v(2)%rib' )
	2798	IF ( ALLOCATED( surf_v(2)%rib ) .AND. kk == 1 ) &
	2799	READ ( 13 ) surf_v(2)%rib
	2800	CASE ( 'surf_v(2)%shf' )
	2801	IF ( ALLOCATED( surf_v(2)%shf ) .AND. kk == 1 ) &
	2802	READ ( 13 ) surf_v(2)%shf
	2803	CASE ( 'surf_v(2)%qsws' )
	2804	IF ( ALLOCATED( surf_v(2)%qsws ) .AND. kk == 1 ) &
	2805	READ ( 13 ) surf_v(2)%qsws
	2806	CASE ( 'surf_v(2)%ssws' )
	2807	IF ( ALLOCATED( surf_v(2)%ssws ) .AND. kk == 1 ) &
	2808	READ ( 13 ) surf_v(2)%ssws
[2292]	2809	CASE ( 'surf_v(2)%qcsws' )
	2810	IF ( ALLOCATED( surf_v(2)%qcsws ) .AND. kk == 1 ) &
	2811	READ ( 13 ) surf_v(2)%qcsws
	2812	CASE ( 'surf_v(2)%ncsws' )
	2813	IF ( ALLOCATED( surf_v(2)%ncsws ) .AND. kk == 1 ) &
	2814	READ ( 13 ) surf_v(2)%ncsws
[2232]	2815	CASE ( 'surf_v(2)%qrsws' )
	2816	IF ( ALLOCATED( surf_v(2)%qrsws ) .AND. kk == 1 ) &
	2817	READ ( 13 ) surf_v(2)%qrsws
	2818	CASE ( 'surf_v(2)%nrsws' )
	2819	IF ( ALLOCATED( surf_v(2)%nrsws ) .AND. kk == 1 ) &
	2820	READ ( 13 ) surf_v(2)%nrsws
	2821	CASE ( 'surf_v(2)%sasws' )
	2822	IF ( ALLOCATED( surf_v(2)%sasws ) .AND. kk == 1 ) &
	2823	READ ( 13 ) surf_v(2)%sasws
	2824	CASE ( 'surf_v(2)%mom_uv' )
	2825	IF ( ALLOCATED( surf_v(2)%mom_flux_uv ) .AND. kk == 1 )&
	2826	READ ( 13 ) surf_v(2)%mom_flux_uv
	2827	CASE ( 'surf_v(2)%mom_w' )
	2828	IF ( ALLOCATED( surf_v(2)%mom_flux_w ) .AND. kk == 1 ) &
	2829	READ ( 13 ) surf_v(2)%mom_flux_w
	2830	CASE ( 'surf_v(2)%mom_tke' )
	2831	IF ( ALLOCATED( surf_v(2)%mom_flux_tke ) .AND. kk == 1 )&
	2832	READ ( 13 ) surf_v(2)%mom_flux_tke
	2833
	2834	CASE ( 'surf_v(3)%start_index' )
	2835	IF ( kk == 1 ) &
	2836	READ ( 13 ) surf_v(3)%start_index
	2837	l = 3
	2838	CASE ( 'surf_v(3)%end_index' )
	2839	IF ( kk == 1 ) &
	2840	READ ( 13 ) surf_v(3)%end_index
	2841	CASE ( 'surf_v(3)%us' )
	2842	IF ( ALLOCATED( surf_v(3)%us ) .AND. kk == 1 ) &
	2843	READ ( 13 ) surf_v(3)%us
	2844	CASE ( 'surf_v(3)%ts' )
	2845	IF ( ALLOCATED( surf_v(3)%ts ) .AND. kk == 1 ) &
	2846	READ ( 13 ) surf_v(3)%ts
	2847	CASE ( 'surf_v(3)%qs' )
	2848	IF ( ALLOCATED( surf_v(3)%qs ) .AND. kk == 1 ) &
	2849	READ ( 13 ) surf_v(3)%qs
	2850	CASE ( 'surf_v(3)%ss' )
	2851	IF ( ALLOCATED( surf_v(3)%ss ) .AND. kk == 1 ) &
	2852	READ ( 13 ) surf_v(3)%ss
[2292]	2853	CASE ( 'surf_v(3)%qcs' )
	2854	IF ( ALLOCATED( surf_v(3)%qcs ) .AND. kk == 1 ) &
	2855	READ ( 13 ) surf_v(3)%qcs
	2856	CASE ( 'surf_v(3)%ncs' )
	2857	IF ( ALLOCATED( surf_v(3)%ncs ) .AND. kk == 1 ) &
	2858	READ ( 13 ) surf_v(3)%ncs
[2232]	2859	CASE ( 'surf_v(3)%qrs' )
	2860	IF ( ALLOCATED( surf_v(3)%qrs ) .AND. kk == 1 ) &
	2861	READ ( 13 ) surf_v(3)%qrs
	2862	CASE ( 'surf_v(3)%nrs' )
	2863	IF ( ALLOCATED( surf_v(3)%nrs ) .AND. kk == 1 ) &
	2864	READ ( 13 ) surf_v(3)%nrs
	2865	CASE ( 'surf_v(3)%ol' )
	2866	IF ( ALLOCATED( surf_v(3)%ol ) .AND. kk == 1 ) &
	2867	READ ( 13 ) surf_v(3)%ol
	2868	CASE ( 'surf_v(3)%rib' )
	2869	IF ( ALLOCATED( surf_v(3)%rib ) .AND. kk == 1 ) &
	2870	READ ( 13 ) surf_v(3)%rib
	2871	CASE ( 'surf_v(3)%shf' )
	2872	IF ( ALLOCATED( surf_v(3)%shf ) .AND. kk == 1 ) &
	2873	READ ( 13 ) surf_v(3)%shf
	2874	CASE ( 'surf_v(3)%qsws' )
	2875	IF ( ALLOCATED( surf_v(3)%qsws ) .AND. kk == 1 ) &
	2876	READ ( 13 ) surf_v(3)%qsws
	2877	CASE ( 'surf_v(3)%ssws' )
	2878	IF ( ALLOCATED( surf_v(3)%ssws ) .AND. kk == 1 ) &
	2879	READ ( 13 ) surf_v(3)%ssws
[2292]	2880	CASE ( 'surf_v(3)%qcsws' )
	2881	IF ( ALLOCATED( surf_v(3)%qcsws ) .AND. kk == 1 ) &
	2882	READ ( 13 ) surf_v(3)%qcsws
	2883	CASE ( 'surf_v(3)%ncsws' )
	2884	IF ( ALLOCATED( surf_v(3)%ncsws ) .AND. kk == 1 ) &
	2885	READ ( 13 ) surf_v(3)%ncsws
[2232]	2886	CASE ( 'surf_v(3)%qrsws' )
	2887	IF ( ALLOCATED( surf_v(3)%qrsws ) .AND. kk == 1 ) &
	2888	READ ( 13 ) surf_v(3)%qrsws
	2889	CASE ( 'surf_v(3)%nrsws' )
	2890	IF ( ALLOCATED( surf_v(3)%nrsws ) .AND. kk == 1 ) &
	2891	READ ( 13 ) surf_v(3)%nrsws
	2892	CASE ( 'surf_v(3)%sasws' )
	2893	IF ( ALLOCATED( surf_v(3)%sasws ) .AND. kk == 1 ) &
	2894	READ ( 13 ) surf_v(3)%sasws
	2895	CASE ( 'surf_v(3)%mom_uv' )
	2896	IF ( ALLOCATED( surf_v(3)%mom_flux_uv ) .AND. kk == 1 )&
	2897	READ ( 13 ) surf_v(3)%mom_flux_uv
	2898	CASE ( 'surf_v(3)%mom_w' )
	2899	IF ( ALLOCATED( surf_v(3)%mom_flux_w ) .AND. kk == 1 ) &
	2900	READ ( 13 ) surf_v(3)%mom_flux_w
	2901	CASE ( 'surf_v(3)%mom_tke' )
	2902	IF ( ALLOCATED( surf_v(3)%mom_flux_tke ) .AND. kk == 1 )&
	2903	READ ( 13 ) surf_v(3)%mom_flux_tke
	2904
	2905	END SELECT
	2906	!
	2907	!-- Redistribute surface elements on its respective type.
	2908	IF ( horizontal_surface ) THEN
	2909	ic = nxlc
	2910	DO i = nxlf, nxrf
	2911	jc = nysc
	2912	DO j = nysf, nynf
	2913
	2914	surf_match_def = surf_def_h(l)%end_index(jc,ic) >= &
	2915	surf_def_h(l)%start_index(jc,ic)
	2916	surf_match_lsm = surf_lsm_h%end_index(jc,ic) >= &
	2917	surf_lsm_h%start_index(jc,ic)
	2918	surf_match_usm = surf_usm_h%end_index(jc,ic) >= &
	2919	surf_usm_h%start_index(jc,ic)
	2920
	2921	IF ( surf_match_def ) THEN
	2922	mm = surf_def_h(l)%start_index(jc,ic)
	2923	DO m = surf_h(l)%start_index(j,i), &
	2924	surf_h(l)%end_index(j,i)
	2925	CALL restore_surface_elements( surf_def_h(l), &
	2926	mm, surf_h(l), m )
	2927	mm = mm + 1
	2928	ENDDO
	2929	ENDIF
	2930
	2931	IF ( surf_match_lsm ) THEN
	2932	mm = surf_lsm_h%start_index(jc,ic)
	2933	DO m = surf_h(l)%start_index(j,i), &
	2934	surf_h(l)%end_index(j,i)
	2935	CALL restore_surface_elements( surf_lsm_h, &
	2936	mm, surf_h(l), m )
	2937	mm = mm + 1
	2938	ENDDO
	2939	ENDIF
	2940
	2941	IF ( surf_match_usm ) THEN
	2942	mm = surf_usm_h%start_index(jc,ic)
	2943	DO m = surf_h(l)%start_index(j,i), &
	2944	surf_h(l)%end_index(j,i)
	2945	CALL restore_surface_elements( surf_usm_h, &
	2946	mm, surf_h(l), m )
	2947	mm = mm + 1
	2948	ENDDO
	2949	ENDIF
	2950
	2951	jc = jc + 1
	2952	ENDDO
	2953	ic = ic + 1
	2954	ENDDO
	2955	ELSEIF ( vertical_surface ) THEN
	2956	ic = nxlc
	2957	DO i = nxlf, nxrf
	2958	jc = nysc
	2959	DO j = nysf, nynf
	2960
	2961	surf_match_def = surf_def_v(l)%end_index(jc,ic) >= &
	2962	surf_def_v(l)%start_index(jc,ic)
	2963	surf_match_lsm = surf_lsm_v(l)%end_index(jc,ic) >= &
	2964	surf_lsm_v(l)%start_index(jc,ic)
	2965	surf_match_usm = surf_usm_v(l)%end_index(jc,ic) >= &
	2966	surf_usm_v(l)%start_index(jc,ic)
	2967
	2968
	2969
	2970	IF ( surf_match_def ) THEN
	2971	mm = surf_def_v(l)%start_index(jc,ic)
	2972	DO m = surf_v(l)%start_index(j,i), &
	2973	surf_v(l)%end_index(j,i)
	2974	CALL restore_surface_elements( surf_def_v(l), &
	2975	mm, surf_v(l), m )
	2976	mm = mm + 1
	2977	ENDDO
	2978	ENDIF
	2979
	2980	IF ( surf_match_lsm ) THEN
	2981	mm = surf_lsm_v(l)%start_index(jc,ic)
	2982	DO m = surf_v(l)%start_index(j,i), &
	2983	surf_v(l)%end_index(j,i)
	2984	CALL restore_surface_elements( surf_lsm_v(l), &
	2985	mm, surf_v(l), m )
	2986	mm = mm + 1
	2987	ENDDO
	2988	ENDIF
	2989
	2990	IF ( surf_match_usm ) THEN
	2991	mm = surf_usm_v(l)%start_index(jc,ic)
	2992	DO m = surf_v(l)%start_index(j,i), &
	2993	surf_v(l)%end_index(j,i)
	2994	CALL restore_surface_elements( surf_usm_v(l), &
	2995	mm, surf_v(l), m )
	2996	mm = mm + 1
	2997	ENDDO
	2998	ENDIF
	2999
	3000	jc = jc + 1
	3001	ENDDO
	3002	ic = ic + 1
	3003	ENDDO
	3004	ENDIF
	3005
	3006	ENDDO
	3007
	3008	READ ( 13 ) field_chr
	3009
	3010	ENDDO
	3011
	3012	ENDIF
	3013
	3014
	3015	CONTAINS
	3016	!------------------------------------------------------------------------------!
	3017	! Description:
	3018	! ------------
	3019	!> Restores surfacle elements back on its respective type.
	3020	!------------------------------------------------------------------------------!
	3021	SUBROUTINE restore_surface_elements( surf_target, m_target, &
	3022	surf_file, m_file )
	3023
	3024	IMPLICIT NONE
	3025
	3026	INTEGER(iwp) :: m_file !< respective surface-element index of current surface array
	3027	INTEGER(iwp) :: m_target !< respecitve surface-element index of surface array on file
	3028
	3029	TYPE( surf_type ) :: surf_target !< target surface type
	3030	TYPE( surf_type ) :: surf_file !< surface type on file
	3031
	3032	IF ( SCAN( TRIM( field_chr ), '%us' ) /= 0 ) THEN
	3033	IF ( ALLOCATED( surf_target%us ) .AND. &
	3034	ALLOCATED( surf_file%us ) ) &
	3035	surf_target%us(m_target) = surf_file%us(m_file)
	3036	ENDIF
	3037
	3038	IF ( SCAN( TRIM( field_chr ), '%ol' ) /= 0 ) THEN
	3039	IF ( ALLOCATED( surf_target%ol ) .AND. &
	3040	ALLOCATED( surf_file%ol ) ) &
	3041	surf_target%ol(m_target) = surf_file%ol(m_file)
	3042	ENDIF
	3043
	3044	IF ( SCAN( TRIM( field_chr ), '%usws' ) /= 0 ) THEN
	3045	IF ( ALLOCATED( surf_target%usws ) .AND. &
	3046	ALLOCATED( surf_file%usws ) ) &
	3047	surf_target%usws(m_target) = surf_file%usws(m_file)
	3048	ENDIF
	3049
	3050	IF ( SCAN( TRIM( field_chr ), '%vsws' ) /= 0 ) THEN
	3051	IF ( ALLOCATED( surf_target%vsws ) .AND. &
	3052	ALLOCATED( surf_file%vsws ) ) &
	3053	surf_target%vsws(m_target) = surf_file%vsws(m_file)
	3054	ENDIF
	3055
	3056	IF ( SCAN( TRIM( field_chr ), '%ts' ) /= 0 ) THEN
	3057	IF ( ALLOCATED( surf_target%ts ) .AND. &
	3058	ALLOCATED( surf_file%ts ) ) &
	3059	surf_target%ts(m_target) = surf_file%ts(m_file)
	3060	ENDIF
	3061
	3062	IF ( SCAN( TRIM( field_chr ), '%shf' ) /= 0 ) THEN
	3063	IF ( ALLOCATED( surf_target%shf ) .AND. &
	3064	ALLOCATED( surf_file%shf ) ) &
	3065	surf_target%shf(m_target) = surf_file%shf(m_file)
	3066	ENDIF
	3067
	3068	IF ( SCAN( TRIM( field_chr ), '%qs' ) /= 0 ) THEN
	3069	IF ( ALLOCATED( surf_target%qs ) .AND. &
	3070	ALLOCATED( surf_file%qs ) ) &
	3071	surf_target%qs(m_target) = surf_file%qs(m_file)
	3072	ENDIF
	3073
	3074	IF ( SCAN( TRIM( field_chr ), '%qsws' ) /= 0 ) THEN
	3075	IF ( ALLOCATED( surf_target%qsws ) .AND. &
	3076	ALLOCATED( surf_file%qsws ) ) &
	3077	surf_target%qsws(m_target) = surf_file%qsws(m_file)
	3078	ENDIF
	3079
	3080	IF ( SCAN( TRIM( field_chr ), '%ss' ) /= 0 ) THEN
	3081	IF ( ALLOCATED( surf_target%ss ) .AND. &
	3082	ALLOCATED( surf_file%ss ) ) &
	3083	surf_target%ss(m_target) = surf_file%ss(m_file)
	3084	ENDIF
	3085
	3086	IF ( SCAN( TRIM( field_chr ), '%ssws' ) /= 0 ) THEN
	3087	IF ( ALLOCATED( surf_target%ssws ) .AND. &
	3088	ALLOCATED( surf_file%ssws ) ) &
	3089	surf_target%ssws(m_target) = surf_file%ssws(m_file)
	3090	ENDIF
	3091
[2292]	3092	IF ( SCAN( TRIM( field_chr ), '%qcs' ) /= 0 ) THEN
	3093	IF ( ALLOCATED( surf_target%qcs ) .AND. &
	3094	ALLOCATED( surf_file%qcs ) ) &
	3095	surf_target%qcs(m_target) = surf_file%qcs(m_file)
	3096	ENDIF
	3097
	3098	IF ( SCAN( TRIM( field_chr ), '%qcsws' ) /= 0 ) THEN
	3099	IF ( ALLOCATED( surf_target%qcsws ) .AND. &
	3100	ALLOCATED( surf_file%qcsws ) ) &
	3101	surf_target%qcsws(m_target) = surf_file%qcsws(m_file)
	3102	ENDIF
	3103
	3104	IF ( SCAN( TRIM( field_chr ), '%ncs' ) /= 0 ) THEN
	3105	IF ( ALLOCATED( surf_target%ncs ) .AND. &
	3106	ALLOCATED( surf_file%ncs ) ) &
	3107	surf_target%ncs(m_target) = surf_file%ncs(m_file)
	3108	ENDIF
	3109
	3110	IF ( SCAN( TRIM( field_chr ), '%ncsws' ) /= 0 ) THEN
	3111	IF ( ALLOCATED( surf_target%ncsws ) .AND. &
	3112	ALLOCATED( surf_file%ncsws ) ) &
	3113	surf_target%ncsws(m_target) = surf_file%ncsws(m_file)
	3114	ENDIF
	3115
[2232]	3116	IF ( SCAN( TRIM( field_chr ), '%qrs' ) /= 0 ) THEN
	3117	IF ( ALLOCATED( surf_target%qrs ) .AND. &
	3118	ALLOCATED( surf_file%qrs ) ) &
	3119	surf_target%qrs(m_target) = surf_file%qrs(m_file)
	3120	ENDIF
	3121
	3122	IF ( SCAN( TRIM( field_chr ), '%qrsws' ) /= 0 ) THEN
	3123	IF ( ALLOCATED( surf_target%qrsws ) .AND. &
	3124	ALLOCATED( surf_file%qrsws ) ) &
	3125	surf_target%qrsws(m_target) = surf_file%qrsws(m_file)
	3126	ENDIF
	3127
	3128	IF ( SCAN( TRIM( field_chr ), '%nrs' ) /= 0 ) THEN
	3129	IF ( ALLOCATED( surf_target%nrs ) .AND. &
	3130	ALLOCATED( surf_file%nrs ) ) &
	3131	surf_target%nrs(m_target) = surf_file%nrs(m_file)
	3132	ENDIF
	3133
	3134	IF ( SCAN( TRIM( field_chr ), '%nrsws' ) /= 0 ) THEN
	3135	IF ( ALLOCATED( surf_target%nrsws ) .AND. &
	3136	ALLOCATED( surf_file%nrsws ) ) &
	3137	surf_target%nrsws(m_target) = surf_file%nrsws(m_file)
	3138	ENDIF
	3139
	3140	IF ( SCAN( TRIM( field_chr ), '%sasws' ) /= 0 ) THEN
	3141	IF ( ALLOCATED( surf_target%sasws ) .AND. &
	3142	ALLOCATED( surf_file%sasws ) ) &
	3143	surf_target%sasws(m_target) = surf_file%sasws(m_file)
	3144	ENDIF
	3145
	3146	IF ( SCAN( TRIM( field_chr ), '%mom_uv' ) /= 0 ) THEN
	3147	IF ( ALLOCATED( surf_target%mom_flux_uv ) .AND. &
	3148	ALLOCATED( surf_file%mom_flux_uv ) ) &
	3149	surf_target%mom_flux_uv(m_target) = &
	3150	surf_file%mom_flux_uv(m_file)
	3151	ENDIF
	3152
	3153	IF ( SCAN( TRIM( field_chr ), '%mom_w' ) /= 0 ) THEN
	3154	IF ( ALLOCATED( surf_target%mom_flux_w ) .AND. &
	3155	ALLOCATED( surf_file%mom_flux_w ) ) &
	3156	surf_target%mom_flux_w(m_target) = &
	3157	surf_file%mom_flux_w(m_file)
	3158	ENDIF
	3159
	3160	IF ( SCAN( TRIM( field_chr ), '%mom_tke' ) /= 0 ) THEN
	3161	IF ( ALLOCATED( surf_target%mom_flux_tke ) .AND. &
	3162	ALLOCATED( surf_file%mom_flux_tke ) ) &
	3163	surf_target%mom_flux_tke(0:1,m_target) = &
	3164	surf_file%mom_flux_tke(0:1,m_file)
	3165	ENDIF
	3166
	3167	END SUBROUTINE restore_surface_elements
	3168
	3169	END SUBROUTINE surface_read_restart_data
	3170
	3171
	3172	!------------------------------------------------------------------------------!
	3173	! Description:
	3174	! ------------
	3175	!> Counts the number of surface elements with the same facing, required for
	3176	!> reading and writing restart data.
	3177	!------------------------------------------------------------------------------!
	3178	SUBROUTINE surface_last_actions
	3179
	3180	IMPLICIT NONE
	3181	!
	3182	!-- Horizontal surfaces
	3183	ns_h_on_file(0) = surf_def_h(0)%ns + surf_lsm_h%ns + surf_usm_h%ns
	3184	ns_h_on_file(1) = surf_def_h(1)%ns
	3185	ns_h_on_file(2) = surf_def_h(2)%ns
	3186	!
	3187	!-- Vertical surfaces
	3188	ns_v_on_file(0) = surf_def_v(0)%ns + surf_lsm_v(0)%ns + surf_usm_v(0)%ns
	3189	ns_v_on_file(1) = surf_def_v(1)%ns + surf_lsm_v(1)%ns + surf_usm_v(1)%ns
	3190	ns_v_on_file(2) = surf_def_v(2)%ns + surf_lsm_v(2)%ns + surf_usm_v(2)%ns
	3191	ns_v_on_file(3) = surf_def_v(3)%ns + surf_lsm_v(3)%ns + surf_usm_v(3)%ns
	3192
	3193	END SUBROUTINE surface_last_actions
	3194
	3195
	3196	END MODULE surface_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |