Home

Context Navigation

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

Last change on this file since 2550 was 2547, checked in by schwenkel, 7 years ago
extended by cloud_droplets option
Property svn:keywords set to `Id`
File size: 162.1 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |