Home

Context Navigation

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

Last change on this file since 2284 was 2270, checked in by maronga, 8 years ago
major revisions in land surface model
Property svn:keywords set to `Id`
File size: 146.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |