Home

Context Navigation

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

Last change on this file since 2579 was 2575, checked in by maronga, 7 years ago
bugfix in radiation model and improvements in land surface scheme
Property svn:keywords set to `Id`
File size: 162.2 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |