Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |