Home

Context Navigation

source: palm/trunk/SOURCE/vertical_nesting_mod.f90 @ 4304

Last change on this file since 4304 was 4182, checked in by scharf, 5 years ago
corrected "Former revisions" section minor formatting in "Former revisions" section added "Author" section
Property svn:keywords set to `Id`
File size: 149.7 KB

Rev	Line
[2365]	1	!> @file vertical_nesting_mod.f90
	2	!------------------------------------------------------------------------------!
[2696]	3	! This file is part of the PALM model system.
[2365]	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	!
[3655]	17	! Copyright 1997-2019 Leibniz Universitaet Hannover
[2718]	18	! Copyright 2017-2018 Karlsruhe Institute of Technology
[2365]	19	!------------------------------------------------------------------------------!
	20	!
	21	! Current revisions:
	22	! -----------------
	23	!
[3049]	24	!
[2365]	25	! Former revisions:
	26	! -----------------
	27	! $Id: vertical_nesting_mod.f90 4182 2019-08-22 15:20:23Z banzhafs $
[4182]	28	! Corrected "Former revisions" section
	29	!
	30	! 4102 2019-07-17 16:00:03Z suehring
[4102]	31	! - Slightly revise setting of boundary conditions at horizontal walls, use
	32	! data-structure offset index instead of pre-calculate it for each facing
	33	!
	34	! 4101 2019-07-17 15:14:26Z gronemeier
[4101]	35	! remove old_dt
	36	!
	37	! 3802 2019-03-17 13:33:42Z raasch
[3802]	38	! unused subroutines commented out
	39	!
	40	! 3655 2019-01-07 16:51:22Z knoop
[3241]	41	! unused variables removed
	42	!
[4182]	43	! 2365 2017-08-21 14:59:59Z kanani
	44	! Initial revision (SadiqHuq)
[3232]	45	!
[4182]	46	!
[2365]	47	! Description:
	48	! ------------
[2374]	49	!> Module for vertical nesting.
	50	!>
	51	!> Definition of parameters and variables for vertical nesting
[2712]	52	!> The horizontal extent of the parent (Coarse Grid) and the child (Fine Grid)
	53	!> have to be identical. The vertical extent of the FG should be smaller than CG.
	54	!> Only integer nesting ratio supported. Odd nesting ratio preferred
	55	!> The code follows MPI-1 standards. The available processors are split into
	56	!> two groups using MPI_COMM_SPLIT. Exchange of data from CG to FG is called
	57	!> interpolation. FG initialization by interpolation is done once at the start.
	58	!> FG boundary conditions are set by interpolated at every timestep.
	59	!> Exchange of data from CG to FG is called anterpolation, the two-way interaction
	60	!> occurs at every timestep.
	61	!> vnest_start_time set in PARIN allows delayed start of the coupling
	62	!> after spin-up of the CG
[2374]	63	!>
[3065]	64	!> @todo Replace dz(1) appropriatly to account for grid stretching
[2374]	65	!> @todo Ensure that code can be compiled for serial and parallel mode. Please
	66	!> check the placement of the directive "__parallel".
	67	!> @todo Add descriptions for all declared variables/parameters, one declaration
	68	!> statement per variable
	69	!> @todo Add a descriptive header above each subroutine (see land_surface_model)
	70	!> @todo FORTRAN language statements must not be used as names for variables
	71	!> (e.g. if). Please rename it.
	72	!> @todo Revise code according to PALM Coding Standard
[2365]	73	!------------------------------------------------------------------------------!
	74	MODULE vertical_nesting_mod
	75
	76	USE kinds
	77
	78	IMPLICIT NONE
	79
[3232]	80	LOGICAL :: vnested = .FALSE. !> set to true if palmrun
	81	!> provides specific information via stdin
[2712]	82	LOGICAL :: vnest_init = .FALSE. !> set to true when FG is initialized
	83	REAL(wp) :: vnest_start_time = 9999999.9 !> simulated time when FG should be
	84	!> initialized. Should be
	85	!> identical in PARIN & PARIN_N
[2365]	86
	87
	88
[2712]	89	INTEGER(iwp),DIMENSION(3,2) :: bdims = 0 !> sub-domain grid topology of current PE
	90	INTEGER(iwp),DIMENSION(3,2) :: bdims_rem = 0 !> sub-domain grid topology of partner PE
[3232]	91	INTEGER(iwp) :: cg_nprocs !> no. of PE in CG. Set by palmrun -Y
	92	INTEGER(iwp) :: fg_nprocs !> no. of PE in FG. Set by palmrun -Y
[2712]	93	INTEGER(iwp) :: TYPE_VNEST_BC !> derived contiguous data type for interpolation
	94	INTEGER(iwp) :: TYPE_VNEST_ANTER !> derived contiguous data type for anterpolation
	95	INTEGER(iwp),DIMENSION(:,:,:),ALLOCATABLE :: c2f_dims_cg !> One CG PE sends data to multiple FG PEs
	96	!> list of grid-topology of partners
	97	INTEGER(iwp),DIMENSION(:,:,:),ALLOCATABLE :: f2c_dims_cg !> One CG PE receives data from multiple FG PEs
	98	!> list of grid-topology of partners
	99	INTEGER(iwp),DIMENSION(:),ALLOCATABLE :: c2f_dims_fg !> One FG PE sends data to multiple CG PE
	100	!> list of grid-topology of partner
	101	INTEGER(iwp),DIMENSION(:),ALLOCATABLE :: f2c_dims_fg !> One FG PE sends data to only one CG PE
	102	!> list of grid-topology of partner
[2365]	103
[2712]	104	INTEGER(iwp),DIMENSION(:,:),ALLOCATABLE :: f_rnk_lst !> list storing rank of FG PE denoted by pdims
	105	INTEGER(iwp),DIMENSION(:,:),ALLOCATABLE :: c_rnk_lst !> list storing rank of CG PE denoted by pdims
	106	INTEGER(iwp),DIMENSION(3) :: cfratio !> Nesting ratio in x,y and z-directions
[2365]	107
[2712]	108	INTEGER(iwp) :: nxc !> no. of CG grid points in x-direction
	109	INTEGER(iwp) :: nxf !> no. of FG grid points in x-direction
	110	INTEGER(iwp) :: nyc !> no. of CG grid points in y-direction
	111	INTEGER(iwp) :: nyf !> no. of FG grid points in y-direction
	112	INTEGER(iwp) :: nzc !> no. of CG grid points in z-direction
	113	INTEGER(iwp) :: nzf !> no. of FG grid points in z-direction
	114	INTEGER(iwp) :: ngp_c !> no. of CG grid points in one vertical level
	115	INTEGER(iwp) :: ngp_f !> no. of FG grid points in one vertical level
[2365]	116
[2712]	117	INTEGER(iwp) :: n_cell_c !> total no. of CG grid points in a PE
	118	INTEGER(iwp),DIMENSION(2) :: pdims_partner !> processor topology of partner PE
	119	INTEGER(iwp) :: target_idex !> temporary variable
	120	INTEGER(iwp),DIMENSION(2) :: offset !> temporary variable
	121	INTEGER(iwp),DIMENSION(2) :: map_coord !> temporary variable
[2365]	122
[2712]	123	REAL(wp) :: dxc !> CG grid pacing in x-direction
	124	REAL(wp) :: dyc !> FG grid pacing in x-direction
	125	REAL(wp) :: dxf !> CG grid pacing in y-direction
	126	REAL(wp) :: dyf !> FG grid pacing in y-direction
	127	REAL(wp) :: dzc !> CG grid pacing in z-direction
	128	REAL(wp) :: dzf !> FG grid pacing in z-direction
	129	REAL(wp) :: dtc !> dt calculated for CG
	130	REAL(wp) :: dtf !> dt calculated for FG
[2365]	131
[2712]	132	REAL(wp), DIMENSION(:), ALLOCATABLE :: zuc !> CG vertical u-levels
	133	REAL(wp), DIMENSION(:), ALLOCATABLE :: zuf !> FG vertical u-levels
	134	REAL(wp), DIMENSION(:), ALLOCATABLE :: zwc !> CG vertical w-levels
	135	REAL(wp), DIMENSION(:), ALLOCATABLE :: zwf !> FG vertical w-levels
	136	REAL(wp), DIMENSION(:,:,:), POINTER :: interpol3d !> pointers to simplify function calls
	137	REAL(wp), DIMENSION(:,:,:), POINTER :: anterpol3d !> pointers to simplify function calls
[2365]	138
	139
[2712]	140	REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: work3d !> temporary array for exchange of 3D data
	141	REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dusws !> temporary array for exchange of 2D data
	142	REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dvsws !> temporary array for exchange of 2D data
	143	REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dts !> temporary array for exchange of 2D data
	144	REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dus !> temporary array for exchange of 2D data
	145
[2365]	146	SAVE
	147
	148	!-- Public functions
	149	PUBLIC vnest_init_fine, vnest_boundary_conds, vnest_anterpolate, &
[2514]	150	vnest_boundary_conds_khkm, vnest_anterpolate_e, &
	151	vnest_init_pegrid_rank, vnest_init_pegrid_domain, vnest_init_grid, &
	152	vnest_timestep_sync, vnest_deallocate
[2365]	153
	154	!-- Public constants and variables
[2712]	155	PUBLIC vnested, vnest_init, vnest_start_time
[2365]	156
	157	PRIVATE bdims, bdims_rem, &
[2712]	158	work3d, work2dusws, work2dvsws, work2dts, work2dus, &
[2365]	159	dxc, dyc, dxf, dyf, dzc, dzf, dtc, dtf, &
	160	zuc, zuf, zwc, zwf, interpol3d, anterpol3d, &
	161	cg_nprocs, fg_nprocs, &
	162	c2f_dims_cg, c2f_dims_fg, f2c_dims_cg, f2c_dims_fg, &
	163	f_rnk_lst, c_rnk_lst, cfratio, pdims_partner, &
	164	nxc, nxf, nyc, nyf, nzc, nzf, &
[3241]	165	ngp_c, ngp_f, target_idex, n_cell_c, &
[2514]	166	offset, map_coord, TYPE_VNEST_BC, TYPE_VNEST_ANTER
[2365]	167
	168	INTERFACE vnest_anterpolate
	169	MODULE PROCEDURE vnest_anterpolate
	170	END INTERFACE vnest_anterpolate
	171
	172	INTERFACE vnest_anterpolate_e
	173	MODULE PROCEDURE vnest_anterpolate_e
	174	END INTERFACE vnest_anterpolate_e
	175
	176	INTERFACE vnest_boundary_conds
	177	MODULE PROCEDURE vnest_boundary_conds
	178	END INTERFACE vnest_boundary_conds
	179
	180	INTERFACE vnest_boundary_conds_khkm
	181	MODULE PROCEDURE vnest_boundary_conds_khkm
	182	END INTERFACE vnest_boundary_conds_khkm
	183
	184	INTERFACE vnest_check_parameters
	185	MODULE PROCEDURE vnest_check_parameters
	186	END INTERFACE vnest_check_parameters
	187
	188	INTERFACE vnest_deallocate
	189	MODULE PROCEDURE vnest_deallocate
	190	END INTERFACE vnest_deallocate
	191
	192	INTERFACE vnest_init_fine
	193	MODULE PROCEDURE vnest_init_fine
	194	END INTERFACE vnest_init_fine
	195
	196	INTERFACE vnest_init_grid
	197	MODULE PROCEDURE vnest_init_grid
	198	END INTERFACE vnest_init_grid
	199
	200	INTERFACE vnest_init_pegrid_domain
	201	MODULE PROCEDURE vnest_init_pegrid_domain
	202	END INTERFACE vnest_init_pegrid_domain
	203
	204	INTERFACE vnest_init_pegrid_rank
	205	MODULE PROCEDURE vnest_init_pegrid_rank
	206	END INTERFACE vnest_init_pegrid_rank
	207
	208	INTERFACE vnest_timestep_sync
	209	MODULE PROCEDURE vnest_timestep_sync
	210	END INTERFACE vnest_timestep_sync
	211
	212	CONTAINS
	213
	214
	215
	216	SUBROUTINE vnest_init_fine
[2712]	217	#if defined( __parallel )
[2365]	218
	219	!--------------------------------------------------------------------------------!
	220	! Description:
	221	! ------------
	222	! At the specified vnest_start_time initialize the Fine Grid based on the coarse
	223	! grid values
	224	!------------------------------------------------------------------------------!
	225
	226
	227	USE arrays_3d
	228	USE control_parameters
	229	USE grid_variables
	230	USE indices
	231	USE interfaces
	232	USE pegrid
[2712]	233	USE turbulence_closure_mod, &
[2696]	234	ONLY : tcm_diffusivities
[2712]	235
[2365]	236
	237	IMPLICIT NONE
	238
[2712]	239	REAL(wp) :: time_since_reference_point_rem
	240	INTEGER(iwp) :: i
	241	INTEGER(iwp) :: j
	242	INTEGER(iwp) :: iif
	243	INTEGER(iwp) :: jjf
	244	INTEGER(iwp) :: kkf
[2365]	245
	246
[2712]	247	if (myid ==0 ) print *, ' TIME TO INIT FINE from COARSE', simulated_time
[2365]	248
	249	!
	250	!-- In case of model termination initiated by the remote model
	251	!-- (terminate_coupled_remote > 0), initiate termination of the local model.
	252	!-- The rest of the coupler must then be skipped because it would cause an MPI
	253	!-- intercomminucation hang.
	254	!-- If necessary, the coupler will be called at the beginning of the next
	255	!-- restart run.
	256
	257	IF ( myid == 0) THEN
[3045]	258	CALL MPI_SENDRECV( terminate_coupled, 1, MPI_INTEGER, &
	259	target_id, 0, &
	260	terminate_coupled_remote, 1, MPI_INTEGER, &
	261	target_id, 0, &
[2365]	262	comm_inter, status, ierr )
	263	ENDIF
	264	CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
	265	ierr )
	266
	267	IF ( terminate_coupled_remote > 0 ) THEN
[3045]	268	WRITE( message_string, * ) 'remote model "', &
	269	TRIM( coupling_mode_remote ), &
	270	'" terminated', &
[3046]	271	'&with terminate_coupled_remote = ', &
[3045]	272	terminate_coupled_remote, &
[3046]	273	'&local model "', TRIM( coupling_mode ), &
[3045]	274	'" has', &
[3046]	275	'&terminate_coupled = ', &
[2365]	276	terminate_coupled
	277	CALL message( 'vnest_init_fine', 'PA0310', 1, 2, 0, 6, 0 )
	278	RETURN
	279	ENDIF
	280
	281
	282	!
	283	!-- Exchange the current simulated time between the models,
	284	!-- currently just for total_2ding
	285	IF ( myid == 0 ) THEN
	286
	287	CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
	288	11, comm_inter, ierr )
	289	CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, &
	290	target_id, 11, comm_inter, status, ierr )
	291
	292	ENDIF
	293
	294	CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
	295	ierr )
	296
	297
	298	IF ( coupling_mode == 'vnested_crse' ) THEN
	299	!-- Send data to fine grid for initialization
	300
	301	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	302	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	303
	304	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	305	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	306	map_coord(1) = i+offset(1)
	307	map_coord(2) = j+offset(2)
	308
	309	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	310
	311	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 10, &
	312	comm_inter,status, ierr )
	313
	314	bdims (1,1) = bdims_rem (1,1) / cfratio(1)
	315	bdims (1,2) = bdims_rem (1,2) / cfratio(1)
	316	bdims (2,1) = bdims_rem (2,1) / cfratio(2)
	317	bdims (2,2) = bdims_rem (2,2) / cfratio(2)
	318	bdims (3,1) = bdims_rem (3,1)
	319	bdims (3,2) = bdims_rem (3,2) / cfratio(3)
	320
	321
	322	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 9, &
	323	comm_inter, ierr )
	324
	325
	326	n_cell_c = (bdims(1,2)-bdims(1,1)+3) * &
	327	(bdims(2,2)-bdims(2,1)+3) * &
	328	(bdims(3,2)-bdims(3,1)+3)
	329
	330	CALL MPI_SEND( u( bdims(3,1):bdims(3,2)+2, &
	331	bdims(2,1)-1:bdims(2,2)+1, &
	332	bdims(1,1)-1:bdims(1,2)+1),&
	333	n_cell_c, MPI_REAL, target_idex, &
[2514]	334	101, comm_inter, ierr)
[2365]	335
	336	CALL MPI_SEND( v( bdims(3,1):bdims(3,2)+2, &
	337	bdims(2,1)-1:bdims(2,2)+1, &
	338	bdims(1,1)-1:bdims(1,2)+1),&
	339	n_cell_c, MPI_REAL, target_idex, &
[2514]	340	102, comm_inter, ierr)
[2365]	341
	342	CALL MPI_SEND( w( bdims(3,1):bdims(3,2)+2, &
	343	bdims(2,1)-1:bdims(2,2)+1, &
	344	bdims(1,1)-1:bdims(1,2)+1),&
	345	n_cell_c, MPI_REAL, target_idex, &
[2514]	346	103, comm_inter, ierr)
[2365]	347
	348	CALL MPI_SEND( pt(bdims(3,1):bdims(3,2)+2, &
	349	bdims(2,1)-1:bdims(2,2)+1, &
	350	bdims(1,1)-1:bdims(1,2)+1),&
	351	n_cell_c, MPI_REAL, target_idex, &
[2514]	352	105, comm_inter, ierr)
[2365]	353
	354	IF ( humidity ) THEN
	355	CALL MPI_SEND( q(bdims(3,1):bdims(3,2)+2, &
	356	bdims(2,1)-1:bdims(2,2)+1, &
	357	bdims(1,1)-1:bdims(1,2)+1),&
	358	n_cell_c, MPI_REAL, target_idex, &
[2514]	359	116, comm_inter, ierr)
[2365]	360	ENDIF
	361
	362	CALL MPI_SEND( e( bdims(3,1):bdims(3,2)+2, &
	363	bdims(2,1)-1:bdims(2,2)+1, &
	364	bdims(1,1)-1:bdims(1,2)+1),&
	365	n_cell_c, MPI_REAL, target_idex, &
[2514]	366	104, comm_inter, ierr)
[2365]	367
	368	CALL MPI_SEND(kh( bdims(3,1):bdims(3,2)+2, &
	369	bdims(2,1)-1:bdims(2,2)+1, &
	370	bdims(1,1)-1:bdims(1,2)+1),&
	371	n_cell_c, MPI_REAL, target_idex, &
[2514]	372	106, comm_inter, ierr)
[2365]	373
	374	CALL MPI_SEND(km( bdims(3,1):bdims(3,2)+2, &
	375	bdims(2,1)-1:bdims(2,2)+1, &
	376	bdims(1,1)-1:bdims(1,2)+1),&
	377	n_cell_c, MPI_REAL, target_idex, &
[2514]	378	107, comm_inter, ierr)
[2365]	379
	380	!-- Send Surface fluxes
	381	IF ( use_surface_fluxes ) THEN
	382
	383	n_cell_c = (bdims(1,2)-bdims(1,1)+3) * &
	384	(bdims(2,2)-bdims(2,1)+3)
	385
[2712]	386	!
	387	!-- shf and z0 for CG / FG need to initialized in input file or user_code
	388	!-- TODO
	389	!-- initialization of usws, vsws, ts and us not vital to vnest FG
	390	!-- variables are not compatible with the new surface layer module
	391	!
	392	! CALL MPI_SEND(surf_def_h(0)%usws( bdims(2,1)-1:bdims(2,2)+1, &
	393	! bdims(1,1)-1:bdims(1,2)+1),&
	394	! n_cell_c, MPI_REAL, target_idex, &
	395	! 110, comm_inter, ierr )
	396	!
	397	! CALL MPI_SEND(surf_def_h(0)%vsws( bdims(2,1)-1:bdims(2,2)+1, &
	398	! bdims(1,1)-1:bdims(1,2)+1),&
	399	! n_cell_c, MPI_REAL, target_idex, &
	400	! 111, comm_inter, ierr )
	401	!
	402	! CALL MPI_SEND(ts ( bdims(2,1)-1:bdims(2,2)+1, &
	403	! bdims(1,1)-1:bdims(1,2)+1),&
	404	! n_cell_c, MPI_REAL, target_idex, &
	405	! 112, comm_inter, ierr )
	406	!
	407	! CALL MPI_SEND(us ( bdims(2,1)-1:bdims(2,2)+1, &
	408	! bdims(1,1)-1:bdims(1,2)+1),&
	409	! n_cell_c, MPI_REAL, target_idex, &
	410	! 113, comm_inter, ierr )
	411	!
[2365]	412	ENDIF
	413
	414
	415
	416
	417	end do
	418	end do
	419
	420	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	421	!-- Receive data from coarse grid for initialization
	422
	423	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	424	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	425	map_coord(1) = offset(1)
	426	map_coord(2) = offset(2)
	427	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	428
	429	bdims (1,1) = nxl
	430	bdims (1,2) = nxr
	431	bdims (2,1) = nys
	432	bdims (2,2) = nyn
	433	bdims (3,1) = nzb
	434	bdims (3,2) = nzt
	435
	436	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 10, &
	437	comm_inter, ierr )
	438
	439	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 9, &
	440	comm_inter,status, ierr )
	441
	442	n_cell_c = (bdims_rem(1,2)-bdims_rem(1,1)+3) * &
	443	(bdims_rem(2,2)-bdims_rem(2,1)+3) * &
	444	(bdims_rem(3,2)-bdims_rem(3,1)+3)
	445
	446	ALLOCATE( work3d ( bdims_rem(3,1) :bdims_rem(3,2)+2, &
	447	bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	448	bdims_rem(1,1)-1:bdims_rem(1,2)+1))
	449
	450
	451	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 101, &
	452	comm_inter,status, ierr )
	453	interpol3d => u
[3241]	454	call interpolate_to_fine_u
[2365]	455
	456	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 102, &
	457	comm_inter,status, ierr )
	458	interpol3d => v
[3241]	459	call interpolate_to_fine_v
[2365]	460
	461	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 103, &
	462	comm_inter,status, ierr )
	463	interpol3d => w
[3241]	464	call interpolate_to_fine_w
[2365]	465
	466	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 105, &
	467	comm_inter,status, ierr )
	468	interpol3d => pt
[3241]	469	call interpolate_to_fine_s
[2365]	470
	471	IF ( humidity ) THEN
	472	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 116, &
	473	comm_inter,status, ierr )
	474	interpol3d => q
[3241]	475	call interpolate_to_fine_s
[2365]	476	ENDIF
	477
	478	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 104, &
	479	comm_inter,status, ierr )
	480	interpol3d => e
[3241]	481	call interpolate_to_fine_s
[2365]	482
	483	!-- kh,km no target attribute, use of pointer not possible
	484	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 106, &
	485	comm_inter,status, ierr )
[3241]	486	call interpolate_to_fine_kh
[2365]	487
	488	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 107, &
	489	comm_inter,status, ierr )
[3241]	490	call interpolate_to_fine_km
[2365]	491
	492	DEALLOCATE( work3d )
	493	NULLIFY ( interpol3d )
	494
[2514]	495	!-- Recv Surface Fluxes
[2365]	496	IF ( use_surface_fluxes ) THEN
	497	n_cell_c = (bdims_rem(1,2)-bdims_rem(1,1)+3) * &
	498	(bdims_rem(2,2)-bdims_rem(2,1)+3)
	499
	500	ALLOCATE( work2dusws ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	501	bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	502	ALLOCATE( work2dvsws ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	503	bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	504	ALLOCATE( work2dts ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	505	bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	506	ALLOCATE( work2dus ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	507	bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	508
[2712]	509	!
	510	!-- shf and z0 for CG / FG need to initialized in input file or user_code
	511	!-- TODO
	512	!-- initialization of usws, vsws, ts and us not vital to vnest FG
	513	!-- variables are not compatible with the new surface layer module
	514	!
	515	! CALL MPI_RECV( work2dusws,n_cell_c, MPI_REAL, target_idex, 110, &
	516	! comm_inter,status, ierr )
	517	!
	518	! CALL MPI_RECV( work2dvsws,n_cell_c, MPI_REAL, target_idex, 111, &
	519	! comm_inter,status, ierr )
	520	!
	521	! CALL MPI_RECV( work2dts ,n_cell_c, MPI_REAL, target_idex, 112, &
	522	! comm_inter,status, ierr )
	523	!
	524	! CALL MPI_RECV( work2dus ,n_cell_c, MPI_REAL, target_idex, 113, &
	525	! comm_inter,status, ierr )
	526	!
	527	! CALL interpolate_to_fine_flux ( 108 )
[2365]	528
	529	DEALLOCATE( work2dusws )
	530	DEALLOCATE( work2dvsws )
	531	DEALLOCATE( work2dts )
	532	DEALLOCATE( work2dus )
	533	ENDIF
	534
	535	IF ( .NOT. constant_diffusion ) THEN
[2712]	536	DO kkf = bdims(3,1)+1,bdims(3,2)+1
	537	DO jjf = bdims(2,1),bdims(2,2)
	538	DO iif = bdims(1,1),bdims(1,2)
[2365]	539
[2712]	540	IF ( e(kkf,jjf,iif) < 0.0 ) THEN
	541	e(kkf,jjf,iif) = 1E-15_wp
[2365]	542	END IF
	543
	544	END DO
	545	END DO
	546	END DO
	547	ENDIF
	548
	549	w(nzt+1,:,:) = w(nzt,:,:)
	550
	551	CALL exchange_horiz( u, nbgp )
	552	CALL exchange_horiz( v, nbgp )
	553	CALL exchange_horiz( w, nbgp )
	554	CALL exchange_horiz( pt, nbgp )
	555	IF ( .NOT. constant_diffusion ) CALL exchange_horiz( e, nbgp )
	556	IF ( humidity ) CALL exchange_horiz( q, nbgp )
	557
	558	!
	559	!-- Velocity boundary conditions at the bottom boundary
	560	IF ( ibc_uv_b == 0 ) THEN
	561	u(nzb,:,:) = 0.0_wp
	562	v(nzb,:,:) = 0.0_wp
	563	ELSE
	564	u(nzb,:,:) = u(nzb+1,:,:)
	565	v(nzb,:,:) = v(nzb+1,:,:)
	566	END IF
	567
	568
	569	w(nzb,:,:) = 0.0_wp
	570
[2514]	571	!
	572	!-- Temperature boundary conditions at the bottom boundary
[2365]	573	IF ( ibc_pt_b /= 0 ) THEN
	574	pt(nzb,:,:) = pt(nzb+1,:,:)
	575	END IF
	576
	577	!
	578	!-- Bottom boundary condition for the turbulent kinetic energy
	579	!-- Generally a Neumann condition with de/dz=0 is assumed
	580	IF ( .NOT. constant_diffusion ) THEN
	581	e(nzb,:,:) = e(nzb+1,:,:)
	582	END IF
	583
	584	!
	585	!-- Bottom boundary condition for turbulent diffusion coefficients
	586	km(nzb,:,:) = km(nzb+1,:,:)
	587	kh(nzb,:,:) = kh(nzb+1,:,:)
	588
	589	!diffusivities required
	590	IF ( .NOT. humidity ) THEN
[2696]	591	CALL tcm_diffusivities( pt, pt_reference )
[2365]	592	ELSE
[2696]	593	CALL tcm_diffusivities( vpt, pt_reference )
[2365]	594	ENDIF
	595
	596
	597	!
	598	!-- Reset Fine Grid top Boundary Condition
	599	!-- At the top of the FG, the scalars always follow Dirichlet condition
	600
	601	ibc_pt_t = 0
	602
	603	!-- Initialize old time levels
	604	pt_p = pt; u_p = u; v_p = v; w_p = w
	605	IF ( .NOT. constant_diffusion ) e_p = e
	606	IF ( humidity ) THEN
	607	ibc_q_t = 0
	608	q_p = q
	609	ENDIF
	610
	611	ENDIF
	612
	613
	614	if (myid==0) print , '* Fine Initalized ** simulated_time:', simulated_time
[2712]	615
[2365]	616	CONTAINS
	617
[3241]	618	SUBROUTINE interpolate_to_fine_w
[2365]	619
	620	USE arrays_3d
	621	USE control_parameters
	622	USE grid_variables
	623	USE indices
	624	USE pegrid
	625
	626
	627	IMPLICIT NONE
	628
[2712]	629	INTEGER(iwp) :: i
	630	INTEGER(iwp) :: j
	631	INTEGER(iwp) :: k
	632	INTEGER(iwp) :: iif
	633	INTEGER(iwp) :: jjf
	634	INTEGER(iwp) :: kkf
	635	INTEGER(iwp) :: nzbottom
	636	INTEGER(iwp) :: nztop
	637	INTEGER(iwp) :: bottomx
	638	INTEGER(iwp) :: bottomy
	639	INTEGER(iwp) :: bottomz
	640	INTEGER(iwp) :: topx
	641	INTEGER(iwp) :: topy
	642	INTEGER(iwp) :: topz
	643	REAL(wp) :: eps
	644	REAL(wp) :: alpha
	645	REAL(wp) :: eminus
	646	REAL(wp) :: edot
	647	REAL(wp) :: eplus
	648	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs
	649	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprf
[2365]	650
	651
	652	nzbottom = bdims_rem (3,1)
	653	nztop = bdims_rem (3,2)
	654
	655	ALLOCATE( wprf(nzbottom:nztop, bdims_rem(2,1)-1: bdims_rem(2,2)+1,nxl:nxr) )
	656	ALLOCATE( wprs(nzbottom:nztop,nys:nyn,nxl:nxr) )
	657
	658
	659	!
	660	!-- Initialisation of the velocity component w
	661	!
	662	!-- Interpolation in x-direction
	663	DO k = nzbottom, nztop
	664	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	665	DO i = bdims_rem(1,1),bdims_rem(1,2)
	666
	667	bottomx = (nxf+1)/(nxc+1) * i
	668	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	669
[2712]	670	DO iif = bottomx, topx
[2365]	671
[2712]	672	eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
[2365]	673	alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	674	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	675	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	676	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	677
[2712]	678	wprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	679	+ edot * work3d(k,j,i) &
	680	+ eplus * work3d(k,j,i+1)
	681	END DO
	682
	683	END DO
	684	END DO
	685	END DO
	686
	687	!
	688	!-- Interpolation in y-direction (quadratic, Clark and Farley)
	689	DO k = nzbottom, nztop
	690	DO j = bdims_rem(2,1), bdims_rem(2,2)
	691
	692	bottomy = (nyf+1)/(nyc+1) * j
	693	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	694
[2712]	695	DO iif = nxl, nxr
	696	DO jjf = bottomy, topy
[2365]	697
[2712]	698	eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
[2365]	699	alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	700	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	701	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	702	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	703
[2712]	704	wprs(k,jjf,iif) = eminus * wprf(k,j-1,iif) &
	705	+ edot * wprf(k,j,iif) &
	706	+ eplus * wprf(k,j+1,iif)
[2365]	707
	708	END DO
	709	END DO
	710
	711	END DO
	712	END DO
	713
	714	!
	715	!-- Interpolation in z-direction (linear)
	716
	717	DO k = nzbottom, nztop-1
	718
	719	bottomz = (dzc/dzf) * k
	720	topz = (dzc/dzf) * (k+1) - 1
	721
[2712]	722	DO jjf = nys, nyn
	723	DO iif = nxl, nxr
	724	DO kkf = bottomz, topz
[2365]	725
[2712]	726	w(kkf,jjf,iif) = wprs(k,jjf,iif) + ( zwf(kkf) - zwc(k) ) &
	727	* ( wprs(k+1,jjf,iif) - wprs(k,jjf,iif) ) / dzc
[2365]	728
	729	END DO
	730	END DO
	731	END DO
	732
	733	END DO
	734
[2712]	735	DO jjf = nys, nyn
	736	DO iif = nxl, nxr
[2365]	737
[2712]	738	w(nzt,jjf,iif) = wprs(nztop,jjf,iif)
[2365]	739
	740	END DO
	741	END DO
	742	!
	743	! w(nzb:nzt+1,nys:nyn,nxl:nxr) = 0
	744
	745	DEALLOCATE( wprf, wprs )
	746
	747	END SUBROUTINE interpolate_to_fine_w
	748
[3241]	749	SUBROUTINE interpolate_to_fine_u
[2365]	750
	751
	752	USE arrays_3d
	753	USE control_parameters
	754	USE grid_variables
	755	USE indices
	756	USE pegrid
	757
	758
	759	IMPLICIT NONE
	760
[2712]	761	INTEGER(iwp) :: i
	762	INTEGER(iwp) :: j
	763	INTEGER(iwp) :: k
	764	INTEGER(iwp) :: iif
	765	INTEGER(iwp) :: jjf
	766	INTEGER(iwp) :: kkf
	767	INTEGER(iwp) :: nzbottom
	768	INTEGER(iwp) :: nztop
	769	INTEGER(iwp) :: bottomx
	770	INTEGER(iwp) :: bottomy
	771	INTEGER(iwp) :: bottomz
	772	INTEGER(iwp) :: topx
	773	INTEGER(iwp) :: topy
	774	INTEGER(iwp) :: topz
	775	REAL(wp) :: eps
	776	REAL(wp) :: alpha
	777	REAL(wp) :: eminus
	778	REAL(wp) :: edot
	779	REAL(wp) :: eplus
	780	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf
	781	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs
[2365]	782
	783
	784
	785	nzbottom = bdims_rem (3,1)
	786	nztop = bdims_rem (3,2)
	787
	788	ALLOCATE( uprf(nzbottom:nztop+2,nys:nyn,bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	789	ALLOCATE( uprs(nzb+1:nzt+1,nys:nyn,bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	790
	791	!
	792	!-- Initialisation of the velocity component uf
	793
	794	!
	795	!-- Interpolation in y-direction (quadratic, Clark and Farley)
	796
	797	DO k = nzbottom, nztop+2
	798	DO j = bdims_rem(2,1), bdims_rem(2,2)
	799
	800	bottomy = (nyf+1)/(nyc+1) * j
	801	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	802
	803	DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1
[2712]	804	DO jjf = bottomy, topy
[2365]	805
[2712]	806	eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
[2365]	807	alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	808	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	809	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	810	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	811
[2712]	812	uprf(k,jjf,i) = eminus * work3d(k,j-1,i) &
[2365]	813	+ edot * work3d(k,j,i) &
	814	+ eplus * work3d(k,j+1,i)
	815
	816	END DO
	817	END DO
	818
	819	END DO
	820	END DO
	821
	822	!
	823	!-- Interpolation in z-direction (quadratic, Clark and Farley)
	824
	825	DO k = nzbottom+1, nztop
	826
	827	bottomz = (dzc/dzf) * (k-1) + 1
	828	topz = (dzc/dzf) * k
	829
[2712]	830	DO jjf = nys, nyn
[2365]	831	DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1
[2712]	832	DO kkf = bottomz, topz
[2365]	833
[2712]	834	eps = ( zuf(kkf) - zuc(k) ) / dzc
[2365]	835	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	836	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	837	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	838	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	839
[2712]	840	uprs(kkf,jjf,i) = eminus * uprf(k-1,jjf,i) &
	841	+ edot * uprf(k,jjf,i) &
	842	+ eplus * uprf(k+1,jjf,i)
[2365]	843
	844	END DO
	845	END DO
	846	END DO
	847
	848	END DO
	849
[2712]	850	DO jjf = nys, nyn
[2365]	851	DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1
	852
	853	eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc
	854	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	855	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	856	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	857	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	858
[2712]	859	uprs(nzt+1,jjf,i) = eminus * uprf(nztop,jjf,i) &
	860	+ edot * uprf(nztop+1,jjf,i) &
	861	+ eplus * uprf(nztop+2,jjf,i)
[2365]	862
	863	END DO
	864	END DO
	865
	866	!
	867	!-- Interpolation in x-direction (linear)
	868
[2712]	869	DO kkf = nzb+1, nzt+1
	870	DO jjf = nys, nyn
[2365]	871	DO i = bdims_rem(1,1), bdims_rem(1,2)
	872
	873	bottomx = (nxf+1)/(nxc+1) * i
	874	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	875
[2712]	876	DO iif = bottomx, topx
	877	u(kkf,jjf,iif) = uprs(kkf,jjf,i) + ( iif * dxf - i * dxc ) &
	878	* ( uprs(kkf,jjf,i+1) - uprs(kkf,jjf,i) ) / dxc
[2365]	879	END DO
	880
	881	END DO
	882	END DO
	883	END DO
	884	!
	885	!-- Determination of uf at the bottom boundary
	886
	887
	888
	889	DEALLOCATE( uprf, uprs )
	890
	891	END SUBROUTINE interpolate_to_fine_u
	892
	893
[3241]	894	SUBROUTINE interpolate_to_fine_v
[2365]	895
	896
	897	USE arrays_3d
	898	USE control_parameters
	899	USE grid_variables
	900	USE indices
	901	USE pegrid
	902
	903
	904	IMPLICIT NONE
[2712]	905
	906	INTEGER(iwp) :: i
	907	INTEGER(iwp) :: j
	908	INTEGER(iwp) :: k
	909	INTEGER(iwp) :: iif
	910	INTEGER(iwp) :: jjf
	911	INTEGER(iwp) :: kkf
	912	INTEGER(iwp) :: nzbottom
	913	INTEGER(iwp) :: nztop
	914	INTEGER(iwp) :: bottomx
	915	INTEGER(iwp) :: bottomy
	916	INTEGER(iwp) :: bottomz
	917	INTEGER(iwp) :: topx
	918	INTEGER(iwp) :: topy
	919	INTEGER(iwp) :: topz
	920	REAL(wp) :: eps
	921	REAL(wp) :: alpha
	922	REAL(wp) :: eminus
	923	REAL(wp) :: edot
	924	REAL(wp) :: eplus
	925	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs
	926	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf
[2365]	927
	928
	929	nzbottom = bdims_rem (3,1)
	930	nztop = bdims_rem (3,2)
	931
	932	ALLOCATE( vprf(nzbottom:nztop+2,bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	933	ALLOCATE( vprs(nzb+1:nzt+1, bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	934	!
	935	!-- Initialisation of the velocity component vf
	936
	937	!
	938	!-- Interpolation in x-direction (quadratic, Clark and Farley)
	939
	940	DO k = nzbottom, nztop+2
	941	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	942	DO i = bdims_rem(1,1), bdims_rem(1,2)
	943
	944	bottomx = (nxf+1)/(nxc+1) * i
	945	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	946
[2712]	947	DO iif = bottomx, topx
[2365]	948
[2712]	949	eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
[2365]	950	alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	951	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	952	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	953	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	954
[2712]	955	vprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	956	+ edot * work3d(k,j,i) &
	957	+ eplus * work3d(k,j,i+1)
	958
	959	END DO
	960
	961	END DO
	962	END DO
	963	END DO
	964
	965	!
	966	!-- Interpolation in z-direction (quadratic, Clark and Farley)
	967
	968	DO k = nzbottom+1, nztop
	969
	970	bottomz = (dzc/dzf) * (k-1) + 1
	971	topz = (dzc/dzf) * k
	972
	973	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
[2712]	974	DO iif = nxl, nxr
	975	DO kkf = bottomz, topz
[2365]	976
[2712]	977	eps = ( zuf(kkf) - zuc(k) ) / dzc
[2365]	978	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	979	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	980	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	981	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	982
[2712]	983	vprs(kkf,j,iif) = eminus * vprf(k-1,j,iif) &
	984	+ edot * vprf(k,j,iif) &
	985	+ eplus * vprf(k+1,j,iif)
[2365]	986
	987	END DO
	988	END DO
	989	END DO
	990
	991	END DO
	992
	993	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
[2712]	994	DO iif = nxl, nxr
[2365]	995
	996	eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc
	997	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	998	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	999	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1000	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1001
[2712]	1002	vprs(nzt+1,j,iif) = eminus * vprf(nztop,j,iif) &
	1003	+ edot * vprf(nztop+1,j,iif) &
	1004	+ eplus * vprf(nztop+2,j,iif)
[2365]	1005
	1006	END DO
	1007	END DO
	1008
	1009	!
	1010	!-- Interpolation in y-direction (linear)
	1011
[2712]	1012	DO kkf = nzb+1, nzt+1
[2365]	1013	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1014
	1015	bottomy = (nyf+1)/(nyc+1) * j
	1016	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1017
[2712]	1018	DO iif = nxl, nxr
	1019	DO jjf = bottomy, topy
	1020	v (kkf,jjf,iif) = vprs(kkf,j,iif) + ( jjf * dyf - j * dyc ) &
	1021	* ( vprs(kkf,j+1,iif) - vprs(kkf,j,iif) ) / dyc
[2365]	1022	END DO
	1023	END DO
	1024
	1025	END DO
	1026	END DO
	1027
	1028	!
	1029	!-- Determination of vf at the bottom boundary
	1030
	1031
	1032	DEALLOCATE( vprf, vprs )
	1033
	1034	END SUBROUTINE interpolate_to_fine_v
	1035
	1036
[3241]	1037	SUBROUTINE interpolate_to_fine_s
[2365]	1038
	1039
	1040	USE arrays_3d
	1041	USE control_parameters
	1042	USE grid_variables
	1043	USE indices
	1044	USE pegrid
	1045
	1046
	1047	IMPLICIT NONE
[2712]	1048
	1049	INTEGER(iwp) :: i
	1050	INTEGER(iwp) :: j
	1051	INTEGER(iwp) :: k
	1052	INTEGER(iwp) :: iif
	1053	INTEGER(iwp) :: jjf
	1054	INTEGER(iwp) :: kkf
	1055	INTEGER(iwp) :: nzbottom
	1056	INTEGER(iwp) :: nztop
	1057	INTEGER(iwp) :: bottomx
	1058	INTEGER(iwp) :: bottomy
	1059	INTEGER(iwp) :: bottomz
	1060	INTEGER(iwp) :: topx
	1061	INTEGER(iwp) :: topy
	1062	INTEGER(iwp) :: topz
	1063	REAL(wp) :: eps
	1064	REAL(wp) :: alpha
	1065	REAL(wp) :: eminus
	1066	REAL(wp) :: edot
	1067	REAL(wp) :: eplus
	1068	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
	1069	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
[2365]	1070
	1071
	1072	nzbottom = bdims_rem (3,1)
	1073	nztop = bdims_rem (3,2)
	1074
	1075	ALLOCATE( ptprf(nzbottom:nztop+2,bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1076	ALLOCATE( ptprs(nzbottom:nztop+2,nys:nyn,nxl:nxr) )
	1077
	1078	!
	1079	!-- Initialisation of scalar variables
	1080
	1081	!
	1082	!-- Interpolation in x-direction (quadratic, Clark and Farley)
	1083
	1084	DO k = nzbottom, nztop+2
	1085	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	1086	DO i = bdims_rem(1,1), bdims_rem(1,2)
	1087
	1088	bottomx = (nxf+1)/(nxc+1) * i
	1089	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1090
[2712]	1091	DO iif = bottomx, topx
[2365]	1092
[2712]	1093	eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
[2365]	1094	alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	1095	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1096	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1097	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1098
[2712]	1099	ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	1100	+ edot * work3d(k,j,i) &
	1101	+ eplus * work3d(k,j,i+1)
	1102	END DO
	1103
	1104	END DO
	1105	END DO
	1106	END DO
	1107
	1108	!
	1109	!-- Interpolation in y-direction (quadratic, Clark and Farley)
	1110
	1111	DO k = nzbottom, nztop+2
	1112	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1113
	1114	bottomy = (nyf+1)/(nyc+1) * j
	1115	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1116
[2712]	1117	DO iif = nxl, nxr
	1118	DO jjf = bottomy, topy
[2365]	1119
[2712]	1120	eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
[2365]	1121	alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	1122	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1123	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1124	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1125
[2712]	1126	ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	1127	+ edot * ptprf(k,j,iif) &
	1128	+ eplus * ptprf(k,j+1,iif)
[2365]	1129
	1130	END DO
	1131	END DO
	1132
	1133	END DO
	1134	END DO
	1135
	1136	!
	1137	!-- Interpolation in z-direction (quadratic, Clark and Farley)
	1138
	1139	DO k = nzbottom+1, nztop
	1140
	1141	bottomz = (dzc/dzf) * (k-1) + 1
	1142	topz = (dzc/dzf) * k
	1143
[2712]	1144	DO jjf = nys, nyn
	1145	DO iif = nxl, nxr
	1146	DO kkf = bottomz, topz
[2365]	1147
[2712]	1148	eps = ( zuf(kkf) - zuc(k) ) / dzc
[2365]	1149	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1150	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1151	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1152	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1153
[2712]	1154	interpol3d(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) &
	1155	+ edot * ptprs(k,jjf,iif) &
	1156	+ eplus * ptprs(k+1,jjf,iif)
[2365]	1157
	1158	END DO
	1159	END DO
	1160	END DO
	1161
	1162	END DO
	1163
[2712]	1164	DO jjf = nys, nyn
	1165	DO iif = nxl, nxr
[2365]	1166
	1167	eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc
	1168	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1169	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1170	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1171	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1172
[2712]	1173	interpol3d(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) &
	1174	+ edot * ptprs(nztop+1,jjf,iif) &
	1175	+ eplus * ptprs(nztop+2,jjf,iif)
[2365]	1176
	1177	END DO
	1178	END DO
	1179
	1180
	1181	DEALLOCATE( ptprf, ptprs )
	1182
	1183	END SUBROUTINE interpolate_to_fine_s
	1184
	1185
[3241]	1186	SUBROUTINE interpolate_to_fine_kh
[2365]	1187
	1188
	1189	USE arrays_3d
	1190	USE control_parameters
	1191	USE grid_variables
	1192	USE indices
	1193	USE pegrid
	1194
	1195
	1196	IMPLICIT NONE
[2712]	1197
	1198	INTEGER(iwp) :: i
	1199	INTEGER(iwp) :: j
	1200	INTEGER(iwp) :: k
	1201	INTEGER(iwp) :: iif
	1202	INTEGER(iwp) :: jjf
	1203	INTEGER(iwp) :: kkf
	1204	INTEGER(iwp) :: nzbottom
	1205	INTEGER(iwp) :: nztop
	1206	INTEGER(iwp) :: bottomx
	1207	INTEGER(iwp) :: bottomy
	1208	INTEGER(iwp) :: bottomz
	1209	INTEGER(iwp) :: topx
	1210	INTEGER(iwp) :: topy
	1211	INTEGER(iwp) :: topz
	1212	REAL(wp) :: eps
	1213	REAL(wp) :: alpha
	1214	REAL(wp) :: eminus
	1215	REAL(wp) :: edot
	1216	REAL(wp) :: eplus
	1217	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
	1218	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
[2365]	1219
	1220
	1221	nzbottom = bdims_rem (3,1)
	1222	nztop = bdims_rem (3,2)
	1223	! nztop = blk_dim_rem (3,2)+1
	1224
	1225
	1226	ALLOCATE( ptprf(nzbottom:nztop+2,bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1227	ALLOCATE( ptprs(nzbottom:nztop+2,nys:nyn,nxl:nxr) )
	1228
	1229
	1230	!
	1231	!-- Initialisation of scalar variables
	1232
	1233	!
	1234	!-- Interpolation in x-direction (quadratic, Clark and Farley)
	1235
	1236	DO k = nzbottom, nztop+2
	1237	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	1238	DO i = bdims_rem(1,1), bdims_rem(1,2)
	1239
	1240	bottomx = (nxf+1)/(nxc+1) * i
	1241	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1242
[2712]	1243	DO iif = bottomx, topx
[2365]	1244
[2712]	1245	eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
[2365]	1246	alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	1247	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1248	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1249	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1250
[2712]	1251	ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	1252	+ edot * work3d(k,j,i) &
	1253	+ eplus * work3d(k,j,i+1)
	1254	END DO
	1255
	1256	END DO
	1257	END DO
	1258	END DO
	1259
	1260	!
	1261	!-- Interpolation in y-direction (quadratic, Clark and Farley)
	1262
	1263	DO k = nzbottom, nztop+2
	1264	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1265
	1266	bottomy = (nyf+1)/(nyc+1) * j
	1267	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1268
[2712]	1269	DO iif = nxl, nxr
	1270	DO jjf = bottomy, topy
[2365]	1271
[2712]	1272	eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
[2365]	1273	alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	1274	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1275	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1276	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1277
[2712]	1278	ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	1279	+ edot * ptprf(k,j,iif) &
	1280	+ eplus * ptprf(k,j+1,iif)
[2365]	1281
	1282	END DO
	1283	END DO
	1284
	1285	END DO
	1286	END DO
	1287
	1288	!
	1289	!-- Interpolation in z-direction (quadratic, Clark and Farley)
	1290
	1291	DO k = nzbottom+1, nztop
	1292
	1293	bottomz = (dzc/dzf) * (k-1) + 1
	1294	topz = (dzc/dzf) * k
	1295
[2712]	1296	DO jjf = nys, nyn
	1297	DO iif = nxl, nxr
	1298	DO kkf = bottomz, topz
[2365]	1299
[2712]	1300	eps = ( zuf(kkf) - zuc(k) ) / dzc
[2365]	1301	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1302	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1303	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1304	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1305
[2712]	1306	kh(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) &
	1307	+ edot * ptprs(k,jjf,iif) &
	1308	+ eplus * ptprs(k+1,jjf,iif)
[2365]	1309
	1310	END DO
	1311	END DO
	1312	END DO
	1313
	1314	END DO
	1315
[2712]	1316	DO jjf = nys, nyn
	1317	DO iif = nxl, nxr
[2365]	1318
	1319	eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc
	1320	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1321	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1322	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1323	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1324
[2712]	1325	kh(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) &
	1326	+ edot * ptprs(nztop+1,jjf,iif) &
	1327	+ eplus * ptprs(nztop+2,jjf,iif)
[2365]	1328
	1329	END DO
	1330	END DO
	1331
	1332
	1333	DEALLOCATE( ptprf, ptprs )
	1334
	1335	END SUBROUTINE interpolate_to_fine_kh
	1336
[3241]	1337	SUBROUTINE interpolate_to_fine_km
[2365]	1338
	1339
	1340	USE arrays_3d
	1341	USE control_parameters
	1342	USE grid_variables
	1343	USE indices
	1344	USE pegrid
	1345
	1346
	1347	IMPLICIT NONE
[2712]	1348
	1349	INTEGER(iwp) :: i
	1350	INTEGER(iwp) :: j
	1351	INTEGER(iwp) :: k
	1352	INTEGER(iwp) :: iif
	1353	INTEGER(iwp) :: jjf
	1354	INTEGER(iwp) :: kkf
	1355	INTEGER(iwp) :: nzbottom
	1356	INTEGER(iwp) :: nztop
	1357	INTEGER(iwp) :: bottomx
	1358	INTEGER(iwp) :: bottomy
	1359	INTEGER(iwp) :: bottomz
	1360	INTEGER(iwp) :: topx
	1361	INTEGER(iwp) :: topy
	1362	INTEGER(iwp) :: topz
	1363	REAL(wp) :: eps
	1364	REAL(wp) :: alpha
	1365	REAL(wp) :: eminus
	1366	REAL(wp) :: edot
	1367	REAL(wp) :: eplus
	1368	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
	1369	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
[2365]	1370
	1371
	1372	nzbottom = bdims_rem (3,1)
	1373	nztop = bdims_rem (3,2)
	1374	! nztop = blk_dim_rem (3,2)+1
	1375
	1376
	1377	ALLOCATE( ptprf(nzbottom:nztop+2,bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1378	ALLOCATE( ptprs(nzbottom:nztop+2,nys:nyn,nxl:nxr) )
	1379
	1380
	1381	!
	1382	!-- Initialisation of scalar variables
	1383
	1384	!
	1385	!-- Interpolation in x-direction (quadratic, Clark and Farley)
	1386
	1387	DO k = nzbottom, nztop+2
	1388	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	1389	DO i = bdims_rem(1,1), bdims_rem(1,2)
	1390
	1391	bottomx = (nxf+1)/(nxc+1) * i
	1392	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1393
[2712]	1394	DO iif = bottomx, topx
[2365]	1395
[2712]	1396	eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
[2365]	1397	alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	1398	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1399	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1400	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1401
[2712]	1402	ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	1403	+ edot * work3d(k,j,i) &
	1404	+ eplus * work3d(k,j,i+1)
	1405	END DO
	1406
	1407	END DO
	1408	END DO
	1409	END DO
	1410
	1411	!
	1412	!-- Interpolation in y-direction (quadratic, Clark and Farley)
	1413
	1414	DO k = nzbottom, nztop+2
	1415	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1416
	1417	bottomy = (nyf+1)/(nyc+1) * j
	1418	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1419
[2712]	1420	DO iif = nxl, nxr
	1421	DO jjf = bottomy, topy
[2365]	1422
[2712]	1423	eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
[2365]	1424	alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	1425	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1426	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1427	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1428
[2712]	1429	ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	1430	+ edot * ptprf(k,j,iif) &
	1431	+ eplus * ptprf(k,j+1,iif)
[2365]	1432
	1433	END DO
	1434	END DO
	1435
	1436	END DO
	1437	END DO
	1438
	1439	!
	1440	!-- Interpolation in z-direction (quadratic, Clark and Farley)
	1441
	1442	DO k = nzbottom+1, nztop
	1443
	1444	bottomz = (dzc/dzf) * (k-1) + 1
	1445	topz = (dzc/dzf) * k
	1446
[2712]	1447	DO jjf = nys, nyn
	1448	DO iif = nxl, nxr
	1449	DO kkf = bottomz, topz
[2365]	1450
[2712]	1451	eps = ( zuf(kkf) - zuc(k) ) / dzc
[2365]	1452	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1453	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1454	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1455	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1456
[2712]	1457	km(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) &
	1458	+ edot * ptprs(k,jjf,iif) &
	1459	+ eplus * ptprs(k+1,jjf,iif)
[2365]	1460
	1461	END DO
	1462	END DO
	1463	END DO
	1464
	1465	END DO
	1466
[2712]	1467	DO jjf = nys, nyn
	1468	DO iif = nxl, nxr
[2365]	1469
	1470	eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc
	1471	alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0
	1472	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1473	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1474	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1475
[2712]	1476	km(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) &
	1477	+ edot * ptprs(nztop+1,jjf,iif) &
	1478	+ eplus * ptprs(nztop+2,jjf,iif)
[2365]	1479
	1480	END DO
	1481	END DO
	1482
	1483
	1484	DEALLOCATE( ptprf, ptprs )
	1485
	1486	END SUBROUTINE interpolate_to_fine_km
	1487
	1488
	1489
	1490
[3802]	1491	! SUBROUTINE interpolate_to_fine_flux
	1492	!
	1493	!
	1494	! USE arrays_3d
	1495	! USE control_parameters
	1496	! USE grid_variables
	1497	! USE indices
	1498	! USE pegrid
	1499	!
	1500	!
	1501	! IMPLICIT NONE
	1502	!
	1503	! INTEGER(iwp) :: i
	1504	! INTEGER(iwp) :: j
	1505	! INTEGER(iwp) :: iif
	1506	! INTEGER(iwp) :: jjf
	1507	! INTEGER(iwp) :: bottomx
	1508	! INTEGER(iwp) :: bottomy
	1509	! INTEGER(iwp) :: topx
	1510	! INTEGER(iwp) :: topy
	1511	! REAL(wp) :: eps
	1512	! REAL(wp) :: alpha
	1513	! REAL(wp) :: eminus
	1514	! REAL(wp) :: edot
	1515	! REAL(wp) :: eplus
	1516	! REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uswspr
	1517	! REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vswspr
	1518	! REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr
	1519	! REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uspr
	1520	!
	1521	! ALLOCATE( uswspr(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1522	! ALLOCATE( vswspr(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1523	! ALLOCATE( tspr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1524	! ALLOCATE( uspr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1525	!
	1526	! !
	1527	! !-- Initialisation of scalar variables (2D)
	1528	!
	1529	! !
	1530	! !-- Interpolation in x-direction (quadratic, Clark and Farley)
	1531	!
	1532	! DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	1533	! DO i = bdims_rem(1,1), bdims_rem(1,2)
	1534	!
	1535	! bottomx = (nxf+1)/(nxc+1) * i
	1536	! topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1537	!
	1538	! DO iif = bottomx, topx
	1539	!
	1540	! eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc
	1541	! alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0
	1542	! eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1543	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1544	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1545	!
	1546	! uswspr(j,iif) = eminus * work2dusws(j,i-1) &
	1547	! + edot * work2dusws(j,i) &
	1548	! + eplus * work2dusws(j,i+1)
	1549	!
	1550	! vswspr(j,iif) = eminus * work2dvsws(j,i-1) &
	1551	! + edot * work2dvsws(j,i) &
	1552	! + eplus * work2dvsws(j,i+1)
	1553	!
	1554	! tspr(j,iif) = eminus * work2dts(j,i-1) &
	1555	! + edot * work2dts(j,i) &
	1556	! + eplus * work2dts(j,i+1)
	1557	!
	1558	! uspr(j,iif) = eminus * work2dus(j,i-1) &
	1559	! + edot * work2dus(j,i) &
	1560	! + eplus * work2dus(j,i+1)
	1561	!
	1562	! END DO
	1563	!
	1564	! END DO
	1565	! END DO
	1566	!
	1567	! !
	1568	! !-- Interpolation in y-direction (quadratic, Clark and Farley)
	1569	!
	1570	! DO j = bdims_rem(2,1), bdims_rem(2,2)
	1571	!
	1572	! bottomy = (nyf+1)/(nyc+1) * j
	1573	! topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1574	!
	1575	! DO iif = nxl, nxr
	1576	! DO jjf = bottomy, topy
	1577	!
	1578	! eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc
	1579	! alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0
	1580	! eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1581	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1582	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1583	!
[2712]	1584	!!
	1585	!!-- TODO
	1586	!-- variables are not compatible with the new surface layer module
	1587	!
	1588	! surf_def_h(0)%usws(jjf,iif) = eminus * uswspr(j-1,if) &
	1589	! + edot * uswspr(j,iif) &
	1590	! + eplus * uswspr(j+1,iif)
	1591	!
	1592	! surf_def_h(0)%vsws(jjf,iif) = eminus * vswspr(j-1,if) &
	1593	! + edot * vswspr(j,iif) &
	1594	! + eplus * vswspr(j+1,iif)
	1595	!
	1596	! ts(jjf,iif) = eminus * tspr(j-1,if) &
	1597	! + edot * tspr(j,iif) &
	1598	! + eplus * tspr(j+1,iif)
	1599	!
	1600	! us(jjf,iif) = eminus * uspr(j-1,if) &
	1601	! + edot * uspr(j,iif) &
	1602	! + eplus * uspr(j+1,iif)
[3802]	1603	!
	1604	! END DO
	1605	! END DO
	1606	!
	1607	! END DO
	1608	!
	1609	!
	1610	! DEALLOCATE( uswspr, vswspr )
	1611	! DEALLOCATE( tspr, uspr )
	1612	!
	1613	!
	1614	! END SUBROUTINE interpolate_to_fine_flux
[2365]	1615
	1616
[2712]	1617	#endif
[2365]	1618	END SUBROUTINE vnest_init_fine
	1619
	1620	SUBROUTINE vnest_boundary_conds
[2712]	1621	#if defined( __parallel )
[2365]	1622	!------------------------------------------------------------------------------!
	1623	! Description:
	1624	! ------------
	1625	! Boundary conditions for the prognostic quantities.
	1626	! One additional bottom boundary condition is applied for the TKE (=(u)*2)
	1627	! in prandtl_fluxes. The cyclic lateral boundary conditions are implicitly
	1628	! handled in routine exchange_horiz. Pressure boundary conditions are
	1629	! explicitly set in routines pres, poisfft, poismg and sor.
	1630	!------------------------------------------------------------------------------!
	1631
	1632	USE arrays_3d
	1633	USE control_parameters
	1634	USE grid_variables
	1635	USE indices
	1636	USE pegrid
	1637
	1638
	1639	IMPLICIT NONE
	1640
[2712]	1641	INTEGER(iwp) :: i
	1642	INTEGER(iwp) :: j
	1643	INTEGER(iwp) :: iif
	1644	INTEGER(iwp) :: jjf
[2365]	1645
	1646
	1647	!
	1648	!-- vnest: top boundary conditions
	1649
	1650	IF ( coupling_mode == 'vnested_crse' ) THEN
	1651	!-- Send data to fine grid for TOP BC
	1652
	1653	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	1654	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	1655
	1656	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	1657	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	1658	map_coord(1) = i+offset(1)
	1659	map_coord(2) = j+offset(2)
	1660
	1661	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	1662
	1663	bdims (1,1) = c2f_dims_cg (0,map_coord(1),map_coord(2))
	1664	bdims (1,2) = c2f_dims_cg (1,map_coord(1),map_coord(2))
	1665	bdims (2,1) = c2f_dims_cg (2,map_coord(1),map_coord(2))
	1666	bdims (2,2) = c2f_dims_cg (3,map_coord(1),map_coord(2))
	1667	bdims (3,1) = c2f_dims_cg (4,map_coord(1),map_coord(2))
	1668	bdims (3,2) = c2f_dims_cg (5,map_coord(1),map_coord(2))
	1669
	1670	n_cell_c = ( (bdims(1,2)-bdims(1,1)) + 3 ) * &
	1671	( (bdims(2,2)-bdims(2,1)) + 3 ) * &
	1672	( (bdims(3,2)-bdims(3,1)) + 1 )
	1673
	1674	CALL MPI_SEND(u (bdims(3,1), bdims(2,1)-1, bdims(1,1)-1), &
	1675	1, TYPE_VNEST_BC, target_idex, &
	1676	201, comm_inter, ierr)
	1677
	1678	CALL MPI_SEND(v(bdims(3,1), bdims(2,1)-1, bdims(1,1)-1),&
	1679	1, TYPE_VNEST_BC, target_idex, &
	1680	202, comm_inter, ierr)
	1681
	1682	CALL MPI_SEND(w(bdims(3,1), bdims(2,1)-1, bdims(1,1)-1),&
	1683	1, TYPE_VNEST_BC, target_idex, &
	1684	203, comm_inter, ierr)
	1685
	1686	CALL MPI_SEND(pt(bdims(3,1), bdims(2,1)-1, bdims(1,1)-1),&
	1687	1, TYPE_VNEST_BC, target_idex, &
	1688	205, comm_inter, ierr)
	1689
	1690	IF ( humidity ) THEN
	1691	CALL MPI_SEND(q(bdims(3,1), bdims(2,1)-1, bdims(1,1)-1),&
	1692	1, TYPE_VNEST_BC, target_idex, &
	1693	209, comm_inter, ierr)
	1694	ENDIF
	1695
	1696	end do
	1697	end do
	1698
	1699	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	1700	!-- Receive data from coarse grid for TOP BC
	1701
	1702	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	1703	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	1704	map_coord(1) = offset(1)
	1705	map_coord(2) = offset(2)
	1706	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	1707
	1708	bdims_rem (1,1) = c2f_dims_fg(0)
	1709	bdims_rem (1,2) = c2f_dims_fg(1)
	1710	bdims_rem (2,1) = c2f_dims_fg(2)
	1711	bdims_rem (2,2) = c2f_dims_fg(3)
	1712	bdims_rem (3,1) = c2f_dims_fg(4)
	1713	bdims_rem (3,2) = c2f_dims_fg(5)
	1714
	1715	n_cell_c = &
	1716	( (bdims_rem(1,2)-bdims_rem(1,1)) + 3 ) * &
	1717	( (bdims_rem(2,2)-bdims_rem(2,1)) + 3 ) * &
	1718	( (bdims_rem(3,2)-bdims_rem(3,1)) + 1 )
	1719
	1720	ALLOCATE( work3d ( &
	1721	bdims_rem(3,1) :bdims_rem(3,2) , &
	1722	bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	1723	bdims_rem(1,1)-1:bdims_rem(1,2)+1))
	1724
	1725
	1726	CALL MPI_RECV( work3d ,n_cell_c, MPI_REAL, target_idex, 201, &
	1727	comm_inter,status, ierr )
	1728	interpol3d => u
	1729	call vnest_set_topbc_u
	1730
	1731	CALL MPI_RECV( work3d ,n_cell_c, MPI_REAL, target_idex, 202, &
	1732	comm_inter,status, ierr )
	1733	interpol3d => v
	1734	call vnest_set_topbc_v
	1735
	1736	CALL MPI_RECV( work3d ,n_cell_c, MPI_REAL, target_idex, 203, &
	1737	comm_inter,status, ierr )
	1738	interpol3d => w
	1739	call vnest_set_topbc_w
	1740
	1741	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 205, &
	1742	comm_inter,status, ierr )
	1743	interpol3d => pt
	1744	call vnest_set_topbc_s
	1745
	1746	IF ( humidity ) THEN
[2514]	1747	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 209, &
[2365]	1748	comm_inter,status, ierr )
	1749	interpol3d => q
	1750	call vnest_set_topbc_s
	1751
	1752	CALL exchange_horiz_2d(q (nzt+1,:,:) )
	1753	ENDIF
	1754
	1755	!-- TKE Neumann BC for FG top
[2712]	1756	DO jjf = nys, nyn
	1757	DO iif = nxl, nxr
	1758	e(nzt+1,jjf,iif) = e(nzt,jjf,iif)
[2365]	1759	END DO
	1760	END DO
	1761
	1762	!
	1763	!-- w level nzt+1 does not impact results. Only to avoid jumps while
	1764	!-- plotting profiles
	1765	w(nzt+1,:,:) = w(nzt,:,:)
	1766
	1767	CALL exchange_horiz_2d(u (nzt+1,:,:) )
	1768	CALL exchange_horiz_2d(v (nzt+1,:,:) )
	1769	CALL exchange_horiz_2d(pt(nzt+1,:,:) )
	1770	CALL exchange_horiz_2d(e (nzt+1,:,:) )
	1771	CALL exchange_horiz_2d(w (nzt+1,:,:) )
	1772	CALL exchange_horiz_2d(w (nzt ,:,:) )
	1773
	1774	NULLIFY ( interpol3d )
	1775	DEALLOCATE ( work3d )
	1776
	1777	ENDIF
	1778
	1779
	1780	CONTAINS
	1781
	1782	SUBROUTINE vnest_set_topbc_w
	1783
	1784
	1785	USE arrays_3d
	1786	USE control_parameters
	1787	USE grid_variables
	1788	USE indices
	1789	USE pegrid
	1790
	1791
	1792	IMPLICIT NONE
[2712]	1793
	1794	INTEGER(iwp) :: i
	1795	INTEGER(iwp) :: j
	1796	INTEGER(iwp) :: iif
	1797	INTEGER(iwp) :: jjf
	1798	INTEGER(iwp) :: bottomx
	1799	INTEGER(iwp) :: bottomy
	1800	INTEGER(iwp) :: topx
	1801	INTEGER(iwp) :: topy
	1802	REAL(wp) :: eps
	1803	REAL(wp) :: alpha
	1804	REAL(wp) :: eminus
	1805	REAL(wp) :: edot
	1806	REAL(wp) :: eplus
	1807	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wprf
[2365]	1808
	1809
	1810	ALLOCATE( wprf(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	1811
	1812	!
	1813	!-- Determination of a boundary condition for the vertical velocity component w:
	1814	!-- In this case only interpolation in x- and y- direction is necessary, as the
	1815	!-- boundary w-node of the fine grid coincides with a w-node in the coarse grid.
	1816	!-- For both interpolations the scheme of Clark and Farley is used.
	1817
	1818	!
	1819	!-- Interpolation in x-direction
	1820
	1821	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	1822
	1823	DO i = bdims_rem(1,1), bdims_rem(1,2)
	1824
	1825	bottomx = (nxf+1)/(nxc+1) * i
	1826	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1827
[2712]	1828	DO iif = bottomx, topx
[2365]	1829
[2712]	1830	eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc
[2365]	1831	alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0
	1832	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1833	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1834	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
[2712]	1835	wprf(j,iif) = eminus * work3d(bdims_rem(3,1),j,i-1) &
[2365]	1836	+ edot * work3d(bdims_rem(3,1),j,i) &
	1837	+ eplus * work3d(bdims_rem(3,1),j,i+1)
	1838
	1839	END DO
	1840
	1841	END DO
	1842
	1843	END DO
	1844
	1845	!
	1846	!-- Interpolation in y-direction
	1847
	1848	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1849
	1850	bottomy = (nyf+1)/(nyc+1) * j
	1851	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1852
[2712]	1853	DO iif = nxl, nxr
[2365]	1854
[2712]	1855	DO jjf = bottomy, topy
[2365]	1856
[2712]	1857	eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc
[2365]	1858
	1859	alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0
	1860
	1861	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1862
	1863	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1864
	1865	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1866
[2712]	1867	w(nzt,jjf,iif) = eminus * wprf(j-1,iif) &
	1868	+ edot * wprf(j,iif) &
	1869	+ eplus * wprf(j+1,iif)
[2365]	1870
	1871	END DO
	1872
	1873	END DO
	1874
	1875	END DO
	1876
	1877	DEALLOCATE( wprf )
	1878
	1879	END SUBROUTINE vnest_set_topbc_w
	1880
	1881
	1882	SUBROUTINE vnest_set_topbc_u
	1883
	1884
	1885	USE arrays_3d
	1886	USE control_parameters
	1887	USE grid_variables
	1888	USE indices
	1889	USE pegrid
	1890
	1891
	1892	IMPLICIT NONE
[2712]	1893
	1894	INTEGER(iwp) :: i
	1895	INTEGER(iwp) :: j
	1896	INTEGER(iwp) :: k
	1897	INTEGER(iwp) :: iif
	1898	INTEGER(iwp) :: jjf
	1899	INTEGER(iwp) :: bottomx
	1900	INTEGER(iwp) :: bottomy
	1901	INTEGER(iwp) :: topx
	1902	INTEGER(iwp) :: topy
	1903	REAL(wp) :: eps
	1904	REAL(wp) :: alpha
	1905	REAL(wp) :: eminus
	1906	REAL(wp) :: edot
	1907	REAL(wp) :: eplus
	1908	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf
	1909	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uprs
[2365]	1910
	1911	ALLOCATE( uprf(bdims_rem(3,1):bdims_rem(3,2),nys:nyn,bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	1912	ALLOCATE( uprs(nys:nyn,bdims_rem(1,1)-1:bdims_rem(1,2)+1) )
	1913
	1914
	1915	!
	1916	!-- Interpolation in y-direction
	1917
	1918	DO k = bdims_rem(3,1), bdims_rem(3,2)
	1919	DO j = bdims_rem(2,1), bdims_rem(2,2)
	1920
	1921	bottomy = (nyf+1)/(nyc+1) * j
	1922	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	1923
	1924	DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1
[2712]	1925	DO jjf = bottomy, topy
[2365]	1926
[2712]	1927	eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc
[2365]	1928	alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0
	1929	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1930	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1931	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	1932
[2712]	1933	uprf(k,jjf,i) = eminus * work3d(k,j-1,i) &
[2365]	1934	+ edot * work3d(k,j,i) &
	1935	+ eplus * work3d(k,j+1,i)
	1936	END DO
	1937	END DO
	1938
	1939	END DO
	1940	END DO
	1941
	1942	!
	1943	!-- Interpolation in z-direction
	1944
[2712]	1945	DO jjf = nys, nyn
[2365]	1946	DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1
	1947	eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc
	1948	alpha = ( (dzf/dzc)**2.0 - 1.0) / 24.0
	1949	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	1950	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	1951	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
[2712]	1952	uprs(jjf,i) = eminus * uprf(bdims_rem(3,1),jjf,i) &
	1953	+ edot * uprf(bdims_rem(3,1)+1,jjf,i) &
	1954	+ eplus * uprf(bdims_rem(3,1)+2,jjf,i)
[2365]	1955	END DO
	1956	END DO
	1957
	1958	!
	1959	!-- Interpolation in x-direction
	1960
[2712]	1961	DO jjf = nys, nyn
[2365]	1962	DO i = bdims_rem(1,1), bdims_rem(1,2)
	1963
	1964	bottomx = (nxf+1)/(nxc+1) * i
	1965	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	1966
[2712]	1967	DO iif = bottomx, topx
	1968	u(nzt+1,jjf,iif) = uprs(jjf,i) + ( iif * dxf - i * dxc ) * ( uprs(jjf,i+1) - uprs(jjf,i) ) / dxc
[2365]	1969	END DO
	1970
	1971	END DO
	1972	END DO
	1973
	1974
	1975
	1976	DEALLOCATE ( uprf, uprs )
	1977
	1978	END SUBROUTINE vnest_set_topbc_u
	1979
	1980
	1981	SUBROUTINE vnest_set_topbc_v
	1982
	1983
	1984	USE arrays_3d
	1985	USE control_parameters
	1986	USE grid_variables
	1987	USE indices
	1988	USE pegrid
	1989
	1990
	1991	IMPLICIT NONE
[2712]	1992
	1993	INTEGER(iwp) :: i
	1994	INTEGER(iwp) :: j
	1995	INTEGER(iwp) :: k
	1996	INTEGER(iwp) :: iif
	1997	INTEGER(iwp) :: jjf
	1998	INTEGER(iwp) :: bottomx
	1999	INTEGER(iwp) :: bottomy
	2000	INTEGER(iwp) :: topx
	2001	INTEGER(iwp) :: topy
	2002	REAL(wp) :: eps
	2003	REAL(wp) :: alpha
	2004	REAL(wp) :: eminus
	2005	REAL(wp) :: edot
	2006	REAL(wp) :: eplus
	2007	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf
	2008	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vprs
[2365]	2009
	2010
	2011
	2012	ALLOCATE( vprf(bdims_rem(3,1):bdims_rem(3,2),bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	2013	ALLOCATE( vprs(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	2014	!
	2015	!-- Determination of a boundary condition for the horizontal velocity component v:
	2016	!-- Interpolation in x- and z-direction is carried out by using the scheme,
	2017	!-- which was derived by Clark and Farley (1984). In y-direction a
	2018	!-- linear interpolation is carried out.
	2019
	2020	!
	2021	!-- Interpolation in x-direction
	2022
	2023	DO k = bdims_rem(3,1), bdims_rem(3,2)
	2024	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	2025	DO i = bdims_rem(1,1), bdims_rem(1,2)
	2026
	2027	bottomx = (nxf+1)/(nxc+1) * i
	2028	topx = (nxf+1)/(nxc+1) * (i+1) - 1
	2029
[2712]	2030	DO iif = bottomx, topx
[2365]	2031
[2712]	2032	eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc
[2365]	2033	alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0
	2034	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	2035	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2036	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
[2712]	2037	vprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	2038	+ edot * work3d(k,j,i) &
	2039	+ eplus * work3d(k,j,i+1)
	2040	END DO
	2041
	2042	END DO
	2043	END DO
	2044	END DO
	2045
	2046	!
	2047	!-- Interpolation in z-direction
	2048
	2049	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
[2712]	2050	DO iif = nxl, nxr
[2365]	2051
	2052	eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc
	2053	alpha = ( (dzf/dzc)**2.0 - 1.0) / 24.0
	2054	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	2055	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2056	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
[2712]	2057	vprs(j,iif) = eminus * vprf(bdims_rem(3,1),j,iif) &
	2058	+ edot * vprf(bdims_rem(3,1)+1,j,iif) &
	2059	+ eplus * vprf(bdims_rem(3,1)+2,j,iif)
[2365]	2060
	2061	END DO
	2062	END DO
	2063
	2064	!
	2065	!-- Interpolation in y-direction
	2066
	2067	DO j = bdims_rem(2,1), bdims_rem(2,2)
[2712]	2068	DO iif = nxl, nxr
[2365]	2069
	2070	bottomy = (nyf+1)/(nyc+1) * j
	2071	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	2072
[2712]	2073	DO jjf = bottomy, topy
[2365]	2074
[2712]	2075	v(nzt+1,jjf,iif) = vprs(j,iif) + ( jjf * dyf - j * dyc ) * ( vprs(j+1,iif) - vprs(j,iif) ) / dyc
[2365]	2076
	2077	END DO
	2078	END DO
	2079	END DO
	2080
	2081
	2082	DEALLOCATE ( vprf, vprs)
	2083
	2084
	2085
	2086	END SUBROUTINE vnest_set_topbc_v
	2087
	2088
	2089	SUBROUTINE vnest_set_topbc_s
	2090
	2091
	2092	USE arrays_3d
	2093	USE control_parameters
	2094	USE grid_variables
	2095	USE indices
	2096	USE pegrid
	2097
	2098
	2099	IMPLICIT NONE
[2712]	2100
	2101	INTEGER(iwp) :: i
	2102	INTEGER(iwp) :: j
	2103	INTEGER(iwp) :: k
	2104	INTEGER(iwp) :: iif
	2105	INTEGER(iwp) :: jjf
	2106	INTEGER(iwp) :: bottomx
	2107	INTEGER(iwp) :: bottomy
	2108	INTEGER(iwp) :: topx
	2109	INTEGER(iwp) :: topy
	2110	REAL(wp) :: eps
	2111	REAL(wp) :: alpha
	2112	REAL(wp) :: eminus
	2113	REAL(wp) :: edot
	2114	REAL(wp) :: eplus
	2115	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
	2116	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
[2365]	2117
	2118
	2119
	2120	ALLOCATE( ptprf(bdims_rem(3,1):bdims_rem(3,2),bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	2121	ALLOCATE( ptprs(bdims_rem(3,1):bdims_rem(3,2),nys:nyn,nxl:nxr) )
	2122
	2123	!
	2124	!-- Determination of a boundary condition for the potential temperature pt:
	2125	!-- The scheme derived by Clark and Farley can be used in all three dimensions.
	2126
	2127	!
	2128	!-- Interpolation in x-direction
	2129
	2130	DO k = bdims_rem(3,1), bdims_rem(3,2)
	2131
	2132	DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	2133
	2134	DO i = bdims_rem(1,1), bdims_rem(1,2)
	2135
	2136	bottomx = (nxf+1)/(nxc+1) * i
	2137	topx = (nxf+1)/(nxc+1) *(i+1) - 1
	2138
[2712]	2139	DO iif = bottomx, topx
[2365]	2140
[2712]	2141	eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc
[2365]	2142
	2143	alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0
	2144
	2145	eminus = eps * (eps - 1.0 ) / 2.0 + alpha
	2146
	2147	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2148
	2149	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2150
[2712]	2151	ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
[2365]	2152	+ edot * work3d(k,j,i) &
	2153	+ eplus * work3d(k,j,i+1)
	2154	END DO
	2155
	2156	END DO
	2157
	2158	END DO
	2159
	2160	END DO
	2161
	2162	!
	2163	!-- Interpolation in y-direction
	2164
	2165	DO k = bdims_rem(3,1), bdims_rem(3,2)
	2166
	2167	DO j = bdims_rem(2,1), bdims_rem(2,2)
	2168
	2169	bottomy = (nyf+1)/(nyc+1) * j
	2170	topy = (nyf+1)/(nyc+1) * (j+1) - 1
	2171
[2712]	2172	DO iif = nxl, nxr
[2365]	2173
[2712]	2174	DO jjf = bottomy, topy
[2365]	2175
[2712]	2176	eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc
[2365]	2177
	2178	alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0
	2179
	2180	eminus = eps * (eps - 1.0) / 2.0 + alpha
	2181
	2182	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2183
	2184	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2185
[2712]	2186	ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	2187	+ edot * ptprf(k,j,iif) &
	2188	+ eplus * ptprf(k,j+1,iif)
[2365]	2189	END DO
	2190
	2191	END DO
	2192
	2193	END DO
	2194
	2195	END DO
	2196
	2197	!
	2198	!-- Interpolation in z-direction
	2199
[2712]	2200	DO jjf = nys, nyn
	2201	DO iif = nxl, nxr
[2365]	2202
	2203	eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc
	2204
	2205	alpha = ( (dzf/dzc)**2.0 - 1.0) / 24.0
	2206
	2207	eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	2208
	2209	edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2210
	2211	eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2212
[2712]	2213	interpol3d (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) &
	2214	+ edot * ptprs(bdims_rem(3,1)+1,jjf,iif) &
	2215	+ eplus * ptprs(bdims_rem(3,1)+2,jjf,iif)
[2365]	2216
	2217	END DO
	2218	END DO
	2219
	2220	DEALLOCATE ( ptprf, ptprs )
	2221
	2222
	2223
	2224	END SUBROUTINE vnest_set_topbc_s
[2712]	2225	#endif
[2365]	2226	END SUBROUTINE vnest_boundary_conds
	2227
	2228
	2229	SUBROUTINE vnest_boundary_conds_khkm
[2712]	2230	#if defined( __parallel )
[2365]	2231
	2232	!--------------------------------------------------------------------------------!
	2233	! Description:
	2234	! ------------
	2235	! Boundary conditions for the prognostic quantities.
	2236	! One additional bottom boundary condition is applied for the TKE (=(u)*2)
	2237	! in prandtl_fluxes. The cyclic lateral boundary conditions are implicitly
	2238	! handled in routine exchange_horiz. Pressure boundary conditions are
	2239	! explicitly set in routines pres, poisfft, poismg and sor.
	2240	!------------------------------------------------------------------------------!
	2241
	2242	USE arrays_3d
	2243	USE control_parameters
	2244	USE grid_variables
	2245	USE indices
	2246	USE pegrid
	2247
	2248
	2249	IMPLICIT NONE
	2250
[2712]	2251	INTEGER(iwp) :: i
	2252	INTEGER(iwp) :: j
	2253	INTEGER(iwp) :: iif
	2254	INTEGER(iwp) :: jjf
[2365]	2255
	2256
	2257	IF ( coupling_mode == 'vnested_crse' ) THEN
	2258	! Send data to fine grid for TOP BC
	2259
	2260	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	2261	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	2262
	2263	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	2264	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	2265	map_coord(1) = i+offset(1)
	2266	map_coord(2) = j+offset(2)
	2267
	2268	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	2269
	2270	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 10, &
	2271	comm_inter,status, ierr )
	2272
	2273	bdims (1,1) = bdims_rem (1,1) / cfratio(1)
	2274	bdims (1,2) = bdims_rem (1,2) / cfratio(1)
	2275	bdims (2,1) = bdims_rem (2,1) / cfratio(2)
	2276	bdims (2,2) = bdims_rem (2,2) / cfratio(2)
	2277	bdims (3,1) = bdims_rem (3,2) / cfratio(3)
	2278	bdims (3,2) = bdims (3,1) + 2
	2279
	2280	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 9, &
	2281	comm_inter, ierr )
	2282
	2283
	2284	n_cell_c = ( (bdims(1,2)-bdims(1,1)) + 3 ) * &
	2285	( (bdims(2,2)-bdims(2,1)) + 3 ) * &
	2286	( (bdims(3,2)-bdims(3,1)) + 1 )
	2287
	2288	CALL MPI_SEND(kh(bdims(3,1) :bdims(3,2) , &
	2289	bdims(2,1)-1:bdims(2,2)+1, &
	2290	bdims(1,1)-1:bdims(1,2)+1),&
	2291	n_cell_c, MPI_REAL, target_idex, &
[2514]	2292	207, comm_inter, ierr)
[2365]	2293
	2294	CALL MPI_SEND(km(bdims(3,1) :bdims(3,2) , &
	2295	bdims(2,1)-1:bdims(2,2)+1, &
	2296	bdims(1,1)-1:bdims(1,2)+1),&
	2297	n_cell_c, MPI_REAL, target_idex, &
[2514]	2298	208, comm_inter, ierr)
[2365]	2299
	2300
	2301
	2302	end do
	2303	end do
	2304
	2305	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	2306	! Receive data from coarse grid for TOP BC
	2307
	2308	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	2309	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	2310	map_coord(1) = offset(1)
	2311	map_coord(2) = offset(2)
	2312	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	2313
	2314	bdims (1,1) = nxl
	2315	bdims (1,2) = nxr
	2316	bdims (2,1) = nys
	2317	bdims (2,2) = nyn
	2318	bdims (3,1) = nzb
	2319	bdims (3,2) = nzt
	2320
	2321	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 10, &
	2322	comm_inter, ierr )
	2323
	2324	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 9, &
	2325	comm_inter,status, ierr )
	2326
	2327	n_cell_c = ( (bdims_rem(1,2)-bdims_rem(1,1)) + 3 ) * &
	2328	( (bdims_rem(2,2)-bdims_rem(2,1)) + 3 ) * &
	2329	( (bdims_rem(3,2)-bdims_rem(3,1)) + 1 )
	2330
	2331	ALLOCATE( work3d ( bdims_rem(3,1) :bdims_rem(3,2) , &
	2332	bdims_rem(2,1)-1:bdims_rem(2,2)+1, &
	2333	bdims_rem(1,1)-1:bdims_rem(1,2)+1))
	2334
	2335
	2336	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 207, &
	2337	comm_inter,status, ierr )
	2338
	2339	! Neumann BC for FG kh
[2712]	2340	DO jjf = nys, nyn
	2341	DO iif = nxl, nxr
	2342	kh(nzt+1,jjf,iif) = kh(nzt,jjf,iif)
[2365]	2343	END DO
	2344	END DO
	2345
	2346	CALL MPI_RECV( work3d,n_cell_c, MPI_REAL, target_idex, 208, &
	2347	comm_inter,status, ierr )
	2348
	2349	! Neumann BC for FG kh
[2712]	2350	DO jjf = nys, nyn
	2351	DO iif = nxl, nxr
	2352	km(nzt+1,jjf,iif) = km(nzt,jjf,iif)
[2365]	2353	END DO
	2354	END DO
	2355
	2356
	2357	!
	2358	!-- The following evaluation can only be performed, if the fine grid is situated below the inversion
[2712]	2359	!! DO jjf = nys-1, nyn+1
	2360	!! DO iif = nxl-1, nxr+1
[2365]	2361	!!
[2712]	2362	!! km(nzt+1,jjf,iif) = 0.1 * l_grid(nzt+1) * SQRT( e(nzt+1,jjf,iif) )
	2363	!! kh(nzt+1,jjf,iif) = 3.0 * km(nzt+1,jjf,iif)
[2365]	2364	!!
	2365	!! END DO
	2366	!! END DO
	2367
	2368	CALL exchange_horiz_2d(km(nzt+1,:,:) )
	2369	CALL exchange_horiz_2d(kh(nzt+1,:,:) )
	2370
	2371	DEALLOCATE ( work3d )
	2372
	2373	ENDIF
	2374
	2375
[3802]	2376	! CONTAINS
	2377	!
	2378	! SUBROUTINE vnest_set_topbc_kh
	2379	!
	2380	!
	2381	! USE arrays_3d
	2382	! USE control_parameters
	2383	! USE grid_variables
	2384	! USE indices
	2385	! USE pegrid
	2386	!
	2387	!
	2388	! IMPLICIT NONE
	2389	!
	2390	! INTEGER(iwp) :: i
	2391	! INTEGER(iwp) :: j
	2392	! INTEGER(iwp) :: k
	2393	! INTEGER(iwp) :: iif
	2394	! INTEGER(iwp) :: jjf
	2395	! INTEGER(iwp) :: bottomx
	2396	! INTEGER(iwp) :: bottomy
	2397	! INTEGER(iwp) :: topx
	2398	! INTEGER(iwp) :: topy
	2399	! REAL(wp) :: eps
	2400	! REAL(wp) :: alpha
	2401	! REAL(wp) :: eminus
	2402	! REAL(wp) :: edot
	2403	! REAL(wp) :: eplus
	2404	! REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
	2405	! REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
	2406	!
	2407	!
	2408	!
	2409	! ALLOCATE( ptprf(bdims_rem(3,1):bdims_rem(3,2),bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	2410	! ALLOCATE( ptprs(bdims_rem(3,1):bdims_rem(3,2),nys:nyn,nxl:nxr) )
	2411	!
	2412	! !
	2413	! !-- Determination of a boundary condition for the potential temperature pt:
	2414	! !-- The scheme derived by Clark and Farley can be used in all three dimensions.
	2415	!
	2416	! !
	2417	! !-- Interpolation in x-direction
	2418	!
	2419	! DO k = bdims_rem(3,1), bdims_rem(3,2)
	2420	!
	2421	! DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	2422	!
	2423	! DO i = bdims_rem(1,1), bdims_rem(1,2)
	2424	!
	2425	! bottomx = (nxf+1)/(nxc+1) * i
	2426	! topx = (nxf+1)/(nxc+1) *(i+1) - 1
	2427	!
	2428	! DO iif = bottomx, topx
	2429	!
	2430	! eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc
	2431	!
	2432	! alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0
	2433	!
	2434	! eminus = eps * (eps - 1.0 ) / 2.0 + alpha
	2435	!
	2436	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2437	!
	2438	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2439	!
	2440	! ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
	2441	! + edot * work3d(k,j,i) &
	2442	! + eplus * work3d(k,j,i+1)
	2443	! END DO
	2444	!
	2445	! END DO
	2446	!
	2447	! END DO
	2448	!
	2449	! END DO
	2450	!
	2451	! !
	2452	! !-- Interpolation in y-direction
	2453	!
	2454	! DO k = bdims_rem(3,1), bdims_rem(3,2)
	2455	!
	2456	! DO j = bdims_rem(2,1), bdims_rem(2,2)
	2457	!
	2458	! bottomy = (nyf+1)/(nyc+1) * j
	2459	! topy = (nyf+1)/(nyc+1) * (j+1) - 1
	2460	!
	2461	! DO iif = nxl, nxr
	2462	!
	2463	! DO jjf = bottomy, topy
	2464	!
	2465	! eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc
	2466	!
	2467	! alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0
	2468	!
	2469	! eminus = eps * (eps - 1.0) / 2.0 + alpha
	2470	!
	2471	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2472	!
	2473	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2474	!
	2475	! ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	2476	! + edot * ptprf(k,j,iif) &
	2477	! + eplus * ptprf(k,j+1,iif)
	2478	! END DO
	2479	!
	2480	! END DO
	2481	!
	2482	! END DO
	2483	!
	2484	! END DO
	2485	!
	2486	! !
	2487	! !-- Interpolation in z-direction
	2488	!
	2489	! DO jjf = nys, nyn
	2490	! DO iif = nxl, nxr
	2491	!
	2492	! eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc
	2493	!
	2494	! alpha = ( (dzf/dzc)**2.0 - 1.0) / 24.0
	2495	!
	2496	! eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	2497	!
	2498	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2499	!
	2500	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2501	!
	2502	! kh (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) &
	2503	! + edot * ptprs(bdims_rem(3,1)+1,jjf,iif) &
	2504	! + eplus * ptprs(bdims_rem(3,1)+2,jjf,iif)
	2505	!
	2506	! END DO
	2507	! END DO
	2508	!
	2509	! DEALLOCATE ( ptprf, ptprs )
	2510	!
	2511	!
	2512	!
	2513	! END SUBROUTINE vnest_set_topbc_kh
[2365]	2514
[3802]	2515	! SUBROUTINE vnest_set_topbc_km
	2516	!
	2517	!
	2518	! USE arrays_3d
	2519	! USE control_parameters
	2520	! USE grid_variables
	2521	! USE indices
	2522	! USE pegrid
	2523	!
	2524	!
	2525	! IMPLICIT NONE
	2526	!
	2527	! INTEGER(iwp) :: i
	2528	! INTEGER(iwp) :: j
	2529	! INTEGER(iwp) :: k
	2530	! INTEGER(iwp) :: iif
	2531	! INTEGER(iwp) :: jjf
	2532	! INTEGER(iwp) :: bottomx
	2533	! INTEGER(iwp) :: bottomy
	2534	! INTEGER(iwp) :: topx
	2535	! INTEGER(iwp) :: topy
	2536	! REAL(wp) :: eps
	2537	! REAL(wp) :: alpha
	2538	! REAL(wp) :: eminus
	2539	! REAL(wp) :: edot
	2540	! REAL(wp) :: eplus
	2541	! REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf
	2542	! REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs
	2543	!
	2544	!
	2545	!
	2546	! ALLOCATE( ptprf(bdims_rem(3,1):bdims_rem(3,2),bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )
	2547	! ALLOCATE( ptprs(bdims_rem(3,1):bdims_rem(3,2),nys:nyn,nxl:nxr) )
	2548	!
	2549	! !
	2550	! !-- Determination of a boundary condition for the potential temperature pt:
	2551	! !-- The scheme derived by Clark and Farley can be used in all three dimensions.
	2552	!
	2553	! !
	2554	! !-- Interpolation in x-direction
	2555	!
	2556	! DO k = bdims_rem(3,1), bdims_rem(3,2)
	2557	!
	2558	! DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1
	2559	!
	2560	! DO i = bdims_rem(1,1), bdims_rem(1,2)
	2561	!
	2562	! bottomx = (nxf+1)/(nxc+1) * i
	2563	! topx = (nxf+1)/(nxc+1) *(i+1) - 1
	2564	!
	2565	! DO iif = bottomx, topx
	2566	!
	2567	! eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc
	2568	!
	2569	! alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0
	2570	!
	2571	! eminus = eps * (eps - 1.0 ) / 2.0 + alpha
	2572	!
	2573	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2574	!
	2575	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2576	!
	2577	! ptprf(k,j,iif) = eminus * work3d(k,j,i-1) &
	2578	! + edot * work3d(k,j,i) &
	2579	! + eplus * work3d(k,j,i+1)
	2580	! END DO
	2581	!
	2582	! END DO
	2583	!
	2584	! END DO
	2585	!
	2586	! END DO
	2587	!
	2588	! !
	2589	! !-- Interpolation in y-direction
	2590	!
	2591	! DO k = bdims_rem(3,1), bdims_rem(3,2)
	2592	!
	2593	! DO j = bdims_rem(2,1), bdims_rem(2,2)
	2594	!
	2595	! bottomy = (nyf+1)/(nyc+1) * j
	2596	! topy = (nyf+1)/(nyc+1) * (j+1) - 1
	2597	!
	2598	! DO iif = nxl, nxr
	2599	!
	2600	! DO jjf = bottomy, topy
	2601	!
	2602	! eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc
	2603	!
	2604	! alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0
	2605	!
	2606	! eminus = eps * (eps - 1.0) / 2.0 + alpha
	2607	!
	2608	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2609	!
	2610	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2611	!
	2612	! ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) &
	2613	! + edot * ptprf(k,j,iif) &
	2614	! + eplus * ptprf(k,j+1,iif)
	2615	! END DO
	2616	!
	2617	! END DO
	2618	!
	2619	! END DO
	2620	!
	2621	! END DO
	2622	!
	2623	! !
	2624	! !-- Interpolation in z-direction
	2625	!
	2626	! DO jjf = nys, nyn
	2627	! DO iif = nxl, nxr
	2628	!
	2629	! eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc
	2630	!
	2631	! alpha = ( (dzf/dzc)**2.0 - 1.0) / 24.0
	2632	!
	2633	! eminus = eps * ( eps - 1.0 ) / 2.0 + alpha
	2634	!
	2635	! edot = ( 1.0 - eps*2.0 ) - 2.0 alpha
	2636	!
	2637	! eplus = eps * ( eps + 1.0 ) / 2.0 + alpha
	2638	!
	2639	! km (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) &
	2640	! + edot * ptprs(bdims_rem(3,1)+1,jjf,iif) &
	2641	! + eplus * ptprs(bdims_rem(3,1)+2,jjf,iif)
	2642	!
	2643	! END DO
	2644	! END DO
	2645	!
	2646	! DEALLOCATE ( ptprf, ptprs )
	2647	!
	2648	!
	2649	!
	2650	! END SUBROUTINE vnest_set_topbc_km
[2365]	2651
	2652
[2712]	2653	#endif
[2365]	2654	END SUBROUTINE vnest_boundary_conds_khkm
	2655
	2656
	2657
	2658	SUBROUTINE vnest_anterpolate
[2712]	2659
	2660	#if defined( __parallel )
[2365]	2661
	2662	!--------------------------------------------------------------------------------!
	2663	! Description:
	2664	! ------------
	2665	! Anterpolate data from fine grid to coarse grid.
	2666	!------------------------------------------------------------------------------!
	2667
	2668	USE arrays_3d
	2669	USE control_parameters
	2670	USE grid_variables
	2671	USE indices
	2672	USE interfaces
	2673	USE pegrid
	2674	USE surface_mod, &
	2675	ONLY : bc_h
	2676
	2677
	2678	IMPLICIT NONE
	2679
[2712]	2680	REAL(wp) :: time_since_reference_point_rem
	2681	INTEGER(iwp) :: i
	2682	INTEGER(iwp) :: j
[4102]	2683	INTEGER(iwp) :: k
[2712]	2684	INTEGER(iwp) :: l !< running index boundary type, for up- and downward-facing walls
	2685	INTEGER(iwp) :: m !< running index surface elements
[2365]	2686
	2687
	2688
	2689	!
	2690	!-- In case of model termination initiated by the remote model
	2691	!-- (terminate_coupled_remote > 0), initiate termination of the local model.
	2692	!-- The rest of the coupler must then be skipped because it would cause an MPI
	2693	!-- intercomminucation hang.
	2694	!-- If necessary, the coupler will be called at the beginning of the next
	2695	!-- restart run.
	2696
	2697	IF ( myid == 0) THEN
	2698	CALL MPI_SENDRECV( terminate_coupled, 1, MPI_INTEGER, &
	2699	target_id, 0, &
	2700	terminate_coupled_remote, 1, MPI_INTEGER, &
	2701	target_id, 0, &
	2702	comm_inter, status, ierr )
	2703	ENDIF
	2704	CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
	2705	ierr )
	2706
	2707	IF ( terminate_coupled_remote > 0 ) THEN
[3045]	2708	WRITE( message_string, * ) 'remote model "', &
	2709	TRIM( coupling_mode_remote ), &
	2710	'" terminated', &
[3046]	2711	'&with terminate_coupled_remote = ', &
[3045]	2712	terminate_coupled_remote, &
[3046]	2713	'&local model "', TRIM( coupling_mode ), &
[3045]	2714	'" has', &
[3046]	2715	'&terminate_coupled = ', &
[2365]	2716	terminate_coupled
	2717	CALL message( 'vnest_anterpolate', 'PA0310', 1, 2, 0, 6, 0 )
	2718	RETURN
	2719	ENDIF
	2720
	2721
	2722	!
	2723	!-- Exchange the current simulated time between the models
	2724
	2725	IF ( myid == 0 ) THEN
	2726
	2727	CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
	2728	11, comm_inter, ierr )
	2729	CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, &
	2730	target_id, 11, comm_inter, status, ierr )
	2731
	2732	ENDIF
	2733
	2734	CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
	2735	ierr )
	2736
	2737	IF ( coupling_mode == 'vnested_crse' ) THEN
	2738	! Receive data from fine grid for anterpolation
	2739
	2740	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	2741	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	2742
	2743	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	2744	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	2745	map_coord(1) = i+offset(1)
	2746	map_coord(2) = j+offset(2)
	2747
	2748	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	2749
	2750	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 10, &
	2751	comm_inter,status, ierr )
	2752
	2753	bdims (1,1) = bdims_rem (1,1) / cfratio(1)
	2754	bdims (1,2) = bdims_rem (1,2) / cfratio(1)
	2755	bdims (2,1) = bdims_rem (2,1) / cfratio(2)
	2756	bdims (2,2) = bdims_rem (2,2) / cfratio(2)
	2757	bdims (3,1) = bdims_rem (3,1)
	2758	bdims (3,2) = bdims_rem (3,2) / cfratio(3)
	2759
	2760	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 9, &
	2761	comm_inter, ierr )
	2762
	2763	n_cell_c = &
	2764	(bdims(1,2)-bdims(1,1)+1) * &
	2765	(bdims(2,2)-bdims(2,1)+1) * &
	2766	(bdims(3,2)-bdims(3,1)+0)
	2767
	2768	CALL MPI_RECV( u( &
	2769	bdims(3,1)+1:bdims(3,2), &
	2770	bdims(2,1) :bdims(2,2), &
	2771	bdims(1,1) :bdims(1,2)),&
	2772	n_cell_c, MPI_REAL, target_idex, 101, &
	2773	comm_inter,status, ierr )
	2774
	2775	CALL MPI_RECV( v( &
	2776	bdims(3,1)+1:bdims(3,2), &
	2777	bdims(2,1) :bdims(2,2), &
	2778	bdims(1,1) :bdims(1,2)),&
	2779	n_cell_c, MPI_REAL, target_idex, 102, &
	2780	comm_inter,status, ierr )
	2781
	2782	CALL MPI_RECV(pt( &
	2783	bdims(3,1)+1:bdims(3,2), &
	2784	bdims(2,1) :bdims(2,2), &
	2785	bdims(1,1) :bdims(1,2)),&
	2786	n_cell_c, MPI_REAL, target_idex, 105, &
	2787	comm_inter,status, ierr )
	2788
[2514]	2789	IF ( humidity ) THEN
[2365]	2790	CALL MPI_RECV(q( &
	2791	bdims(3,1)+1:bdims(3,2), &
	2792	bdims(2,1) :bdims(2,2), &
	2793	bdims(1,1) :bdims(1,2)),&
	2794	n_cell_c, MPI_REAL, target_idex, 106, &
	2795	comm_inter,status, ierr )
[2514]	2796	ENDIF
[2365]	2797
	2798	CALL MPI_RECV( w( &
	2799	bdims(3,1) :bdims(3,2)-1, &
	2800	bdims(2,1) :bdims(2,2), &
	2801	bdims(1,1) :bdims(1,2)), &
	2802	n_cell_c, MPI_REAL, target_idex, 103, &
	2803	comm_inter,status, ierr )
	2804
	2805	end do
	2806	end do
	2807
	2808
	2809
	2810	!
	2811	!-- Boundary conditions for the velocity components u and v
	2812
	2813
	2814	IF ( ibc_uv_b == 0 ) THEN
	2815	u(nzb,:,:) = 0.0_wp
	2816	v(nzb,:,:) = 0.0_wp
	2817	ELSE
	2818	u(nzb,:,:) = u(nzb+1,:,:)
	2819	v(nzb,:,:) = v(nzb+1,:,:)
	2820	END IF
	2821	!
	2822	!-- Boundary conditions for the velocity components w
	2823
	2824	w(nzb,:,:) = 0.0_wp
	2825
	2826	!
	2827	!-- Temperature at bottom boundary.
	2828	!-- Neumann, zero-gradient
	2829	IF ( ibc_pt_b == 1 ) THEN
	2830	DO l = 0, 1
[4102]	2831	DO m = 1, bc_h(l)%ns
	2832	i = bc_h(l)%i(m)
	2833	j = bc_h(l)%j(m)
	2834	k = bc_h(l)%k(m)
	2835	pt(k+bc_h(l)%koff,j,i) = pt(k,j,i)
	2836	ENDDO
	2837	ENDDO
[2365]	2838	ENDIF
	2839
	2840
	2841	CALL exchange_horiz( u, nbgp )
	2842	CALL exchange_horiz( v, nbgp )
	2843	CALL exchange_horiz( w, nbgp )
	2844	CALL exchange_horiz( pt, nbgp )
	2845
	2846
	2847	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	2848	! Send data to coarse grid for anterpolation
	2849
	2850	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	2851	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	2852	map_coord(1) = offset(1)
	2853	map_coord(2) = offset(2)
	2854	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	2855
	2856	!-- Limit anterpolation level to nzt - z nesting ratio (a pseudo-buffer layer)
	2857	bdims (1,1) = nxl
	2858	bdims (1,2) = nxr
	2859	bdims (2,1) = nys
	2860	bdims (2,2) = nyn
	2861	bdims (3,1) = nzb
	2862	bdims (3,2) = nzt-cfratio(3)
	2863
	2864	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 10, &
	2865	comm_inter, ierr )
	2866
	2867	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 9, &
	2868	comm_inter,status, ierr )
	2869
	2870
	2871	ALLOCATE( work3d ( &
	2872	bdims_rem(3,1)+1:bdims_rem(3,2), &
	2873	bdims_rem(2,1) :bdims_rem(2,2), &
	2874	bdims_rem(1,1) :bdims_rem(1,2)))
	2875
	2876
	2877	anterpol3d => u
[2514]	2878
[3241]	2879	CALL anterpolate_to_crse_u
[2365]	2880	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	2881	101, comm_inter, ierr)
[2365]	2882
	2883	anterpol3d => v
	2884
[3241]	2885	CALL anterpolate_to_crse_v
[2365]	2886	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	2887	102, comm_inter, ierr)
[2365]	2888
	2889	anterpol3d => pt
	2890
[3241]	2891	CALL anterpolate_to_crse_s
[2365]	2892	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	2893	105, comm_inter, ierr)
[2365]	2894
	2895
	2896	IF ( humidity ) THEN
	2897
	2898	anterpol3d => q
	2899
[3241]	2900	CALL anterpolate_to_crse_s
[2365]	2901	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	2902	106, comm_inter, ierr)
[2365]	2903	ENDIF
	2904
	2905
	2906	DEALLOCATE( work3d )
	2907	ALLOCATE( work3d ( bdims_rem(3,1) :bdims_rem(3,2)-1, &
	2908	bdims_rem(2,1) :bdims_rem(2,2), &
	2909	bdims_rem(1,1) :bdims_rem(1,2)))
	2910	anterpol3d => w
[3241]	2911	CALL anterpolate_to_crse_w
[2365]	2912	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	2913	103, comm_inter, ierr)
[2365]	2914
	2915	NULLIFY ( anterpol3d )
	2916	DEALLOCATE( work3d )
	2917
	2918	ENDIF
	2919
	2920
	2921
	2922	CONTAINS
[3241]	2923	SUBROUTINE anterpolate_to_crse_u
[2365]	2924
	2925
	2926	USE arrays_3d
	2927	USE control_parameters
	2928	USE grid_variables
	2929	USE indices
	2930	USE pegrid
	2931
	2932
	2933	IMPLICIT NONE
[2712]	2934
	2935	INTEGER(iwp) :: i
	2936	INTEGER(iwp) :: j
	2937	INTEGER(iwp) :: k
	2938	INTEGER(iwp) :: iif
	2939	INTEGER(iwp) :: jjf
	2940	INTEGER(iwp) :: kkf
	2941	INTEGER(iwp) :: bottomy
	2942	INTEGER(iwp) :: bottomz
	2943	INTEGER(iwp) :: topy
	2944	INTEGER(iwp) :: topz
	2945	REAL(wp) :: aweight
[2365]	2946
	2947	!
	2948	!-- Anterpolation of the velocity components u
	2949	!-- only values in yz-planes that coincide in the fine and
	2950	!-- the coarse grid are considered
	2951
	2952	DO k = bdims_rem(3,1)+1, bdims_rem(3,2)
	2953
	2954	bottomz = (dzc/dzf) * (k-1) + 1
	2955	topz = (dzc/dzf) * k
	2956
	2957	DO j = bdims_rem(2,1),bdims_rem(2,2)
	2958
	2959	bottomy = (nyf+1) / (nyc+1) * j
	2960	topy = (nyf+1) / (nyc+1) * (j+1) - 1
	2961
	2962	DO i = bdims_rem(1,1),bdims_rem(1,2)
	2963
[2712]	2964	iif = (nxf+1) / (nxc+1) * i
[2365]	2965
	2966	aweight = 0.0
	2967
[2712]	2968	DO kkf = bottomz, topz
	2969	DO jjf = bottomy, topy
[2365]	2970
[2712]	2971	aweight = aweight + anterpol3d(kkf,jjf,iif) * &
[2365]	2972	(dzf/dzc) * (dyf/dyc)
	2973
	2974	END DO
	2975	END DO
	2976
	2977	work3d(k,j,i) = aweight
	2978
	2979	END DO
	2980
	2981	END DO
	2982
	2983	END DO
	2984
	2985
	2986
	2987	END SUBROUTINE anterpolate_to_crse_u
	2988
	2989
[3241]	2990	SUBROUTINE anterpolate_to_crse_v
[2365]	2991
	2992
	2993	USE arrays_3d
	2994	USE control_parameters
	2995	USE grid_variables
	2996	USE indices
	2997	USE pegrid
	2998
	2999
	3000	IMPLICIT NONE
[2712]	3001
	3002	INTEGER(iwp) :: i
	3003	INTEGER(iwp) :: j
	3004	INTEGER(iwp) :: k
	3005	INTEGER(iwp) :: iif
	3006	INTEGER(iwp) :: jjf
	3007	INTEGER(iwp) :: kkf
	3008	INTEGER(iwp) :: bottomx
	3009	INTEGER(iwp) :: bottomz
	3010	INTEGER(iwp) :: topx
	3011	INTEGER(iwp) :: topz
	3012	REAL(wp) :: aweight
	3013
[2365]	3014	!
	3015	!-- Anterpolation of the velocity components v
	3016	!-- only values in xz-planes that coincide in the fine and
	3017	!-- the coarse grid are considered
	3018
	3019	DO k = bdims_rem(3,1)+1, bdims_rem(3,2)
	3020
	3021	bottomz = (dzc/dzf) * (k-1) + 1
	3022	topz = (dzc/dzf) * k
	3023
	3024	DO j = bdims_rem(2,1), bdims_rem(2,2)
	3025
[2712]	3026	jjf = (nyf+1) / (nyc+1) * j
[2365]	3027
	3028	DO i = bdims_rem(1,1), bdims_rem(1,2)
	3029
	3030	bottomx = (nxf+1) / (nxc+1) * i
	3031	topx = (nxf+1) / (nxc+1) * (i+1) - 1
	3032
	3033	aweight = 0.0
	3034
[2712]	3035	DO kkf = bottomz, topz
	3036	DO iif = bottomx, topx
[2365]	3037
[2712]	3038	aweight = aweight + anterpol3d(kkf,jjf,iif) * &
[2365]	3039	(dzf/dzc) * (dxf/dxc)
	3040
	3041
	3042	END DO
	3043	END DO
	3044
	3045	work3d(k,j,i) = aweight
	3046
	3047	END DO
	3048	END DO
	3049	END DO
	3050
	3051
	3052
	3053	END SUBROUTINE anterpolate_to_crse_v
	3054
	3055
[3241]	3056	SUBROUTINE anterpolate_to_crse_w
[2365]	3057
	3058
	3059	USE arrays_3d
	3060	USE control_parameters
	3061	USE grid_variables
	3062	USE indices
	3063	USE pegrid
	3064
	3065
	3066	IMPLICIT NONE
[2712]	3067
	3068	INTEGER(iwp) :: i
	3069	INTEGER(iwp) :: j
	3070	INTEGER(iwp) :: k
	3071	INTEGER(iwp) :: iif
	3072	INTEGER(iwp) :: jjf
	3073	INTEGER(iwp) :: kkf
	3074	INTEGER(iwp) :: bottomx
	3075	INTEGER(iwp) :: bottomy
	3076	INTEGER(iwp) :: topx
	3077	INTEGER(iwp) :: topy
	3078	REAL(wp) :: aweight
	3079
[2365]	3080	!
	3081	!-- Anterpolation of the velocity components w
	3082	!-- only values in xy-planes that coincide in the fine and
	3083	!-- the coarse grid are considered
	3084
	3085	DO k = bdims_rem(3,1), bdims_rem(3,2)-1
	3086
[2712]	3087	kkf = cfratio(3) * k
[2365]	3088
	3089	DO j = bdims_rem(2,1), bdims_rem(2,2)
	3090
	3091	bottomy = (nyf+1) / (nyc+1) * j
	3092	topy = (nyf+1) / (nyc+1) * (j+1) - 1
	3093
	3094	DO i = bdims_rem(1,1), bdims_rem(1,2)
	3095
	3096	bottomx = (nxf+1) / (nxc+1) * i
	3097	topx = (nxf+1) / (nxc+1) * (i+1) - 1
	3098
	3099	aweight = 0.0
	3100
[2712]	3101	DO jjf = bottomy, topy
	3102	DO iif = bottomx, topx
[2365]	3103
[2712]	3104	aweight = aweight + anterpol3d (kkf,jjf,iif) * &
[2365]	3105	(dxf/dxc) * (dyf/dyc)
	3106
	3107	END DO
	3108	END DO
	3109
	3110	work3d(k,j,i) = aweight
	3111
	3112	END DO
	3113
	3114	END DO
	3115
	3116	END DO
	3117
	3118
	3119	END SUBROUTINE anterpolate_to_crse_w
	3120
	3121
[3241]	3122	SUBROUTINE anterpolate_to_crse_s
[2365]	3123
	3124
	3125	USE arrays_3d
	3126	USE control_parameters
	3127	USE grid_variables
	3128	USE indices
	3129	USE pegrid
	3130
	3131
	3132	IMPLICIT NONE
[2712]	3133
	3134	INTEGER(iwp) :: i
	3135	INTEGER(iwp) :: j
	3136	INTEGER(iwp) :: k
	3137	INTEGER(iwp) :: iif
	3138	INTEGER(iwp) :: jjf
	3139	INTEGER(iwp) :: kkf
	3140	INTEGER(iwp) :: bottomx
	3141	INTEGER(iwp) :: bottomy
	3142	INTEGER(iwp) :: bottomz
	3143	INTEGER(iwp) :: topx
	3144	INTEGER(iwp) :: topy
	3145	INTEGER(iwp) :: topz
	3146	REAL(wp) :: aweight
[2365]	3147
	3148	!
	3149	!-- Anterpolation of the potential temperature pt
	3150	!-- all fine grid values are considered
	3151
	3152	DO k = bdims_rem(3,1)+1, bdims_rem(3,2)
	3153
	3154	bottomz = (dzc/dzf) * (k-1) + 1
	3155	topz = (dzc/dzf) * k
	3156
	3157	DO j = bdims_rem(2,1), bdims_rem(2,2)
	3158
	3159	bottomy = (nyf+1) / (nyc+1) * j
	3160	topy = (nyf+1) / (nyc+1) * (j+1) - 1
	3161
	3162	DO i = bdims_rem(1,1), bdims_rem(1,2)
	3163
	3164	bottomx = (nxf+1) / (nxc+1) * i
	3165	topx = (nxf+1) / (nxc+1) * (i+1) - 1
	3166
	3167	aweight = 0.0
	3168
[2712]	3169	DO kkf = bottomz, topz
	3170	DO jjf = bottomy, topy
	3171	DO iif = bottomx, topx
[2365]	3172
[2712]	3173	aweight = aweight + anterpol3d(kkf,jjf,iif) * &
[2365]	3174	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3175
	3176	END DO
	3177	END DO
	3178	END DO
	3179
	3180	work3d(k,j,i) = aweight
	3181
	3182	END DO
	3183
	3184	END DO
	3185
	3186	END DO
	3187
	3188
	3189	END SUBROUTINE anterpolate_to_crse_s
[2712]	3190	#endif
[2365]	3191	END SUBROUTINE vnest_anterpolate
	3192
	3193
	3194
	3195	SUBROUTINE vnest_anterpolate_e
[2712]	3196	#if defined( __parallel )
[2365]	3197
	3198	!--------------------------------------------------------------------------------!
	3199	! Description:
	3200	! ------------
	3201	! Anterpolate TKE from fine grid to coarse grid.
	3202	!------------------------------------------------------------------------------!
	3203
	3204	USE arrays_3d
	3205	USE control_parameters
	3206	USE grid_variables
	3207	USE indices
	3208	USE interfaces
	3209	USE pegrid
	3210
	3211
	3212	IMPLICIT NONE
	3213
[2712]	3214	REAL(wp) :: time_since_reference_point_rem
	3215	INTEGER(iwp) :: i
	3216	INTEGER(iwp) :: j
[2365]	3217
	3218	!
	3219	!-- In case of model termination initiated by the remote model
	3220	!-- (terminate_coupled_remote > 0), initiate termination of the local model.
	3221	!-- The rest of the coupler must then be skipped because it would cause an MPI
	3222	!-- intercomminucation hang.
	3223	!-- If necessary, the coupler will be called at the beginning of the next
	3224	!-- restart run.
	3225
	3226	IF ( myid == 0) THEN
[3045]	3227	CALL MPI_SENDRECV( terminate_coupled, 1, MPI_INTEGER, &
	3228	target_id, 0, &
	3229	terminate_coupled_remote, 1, MPI_INTEGER, &
	3230	target_id, 0, &
[2365]	3231	comm_inter, status, ierr )
	3232	ENDIF
[3045]	3233	CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
[2365]	3234	ierr )
	3235
	3236	IF ( terminate_coupled_remote > 0 ) THEN
[3045]	3237	WRITE( message_string, * ) 'remote model "', &
	3238	TRIM( coupling_mode_remote ), &
	3239	'" terminated', &
[3046]	3240	'&with terminate_coupled_remote = ', &
[3045]	3241	terminate_coupled_remote, &
[3046]	3242	'&local model "', TRIM( coupling_mode ), &
[3045]	3243	'" has', &
[3046]	3244	'&terminate_coupled = ', &
[2365]	3245	terminate_coupled
	3246	CALL message( 'vnest_anterpolate_e', 'PA0310', 1, 2, 0, 6, 0 )
	3247	RETURN
	3248	ENDIF
	3249
	3250
	3251	!
	3252	!-- Exchange the current simulated time between the models
	3253	IF ( myid == 0 ) THEN
	3254
	3255	CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
	3256	11, comm_inter, ierr )
	3257	CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, &
	3258	target_id, 11, comm_inter, status, ierr )
	3259
	3260	ENDIF
	3261
	3262	CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
	3263	ierr )
	3264
	3265	IF ( coupling_mode == 'vnested_crse' ) THEN
	3266	! Receive data from fine grid for anterpolation
	3267
	3268	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	3269	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	3270
	3271	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	3272	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	3273	map_coord(1) = i+offset(1)
	3274	map_coord(2) = j+offset(2)
	3275
	3276	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	3277
	3278	bdims (1,1) = f2c_dims_cg (0,map_coord(1),map_coord(2))
	3279	bdims (1,2) = f2c_dims_cg (1,map_coord(1),map_coord(2))
	3280	bdims (2,1) = f2c_dims_cg (2,map_coord(1),map_coord(2))
	3281	bdims (2,2) = f2c_dims_cg (3,map_coord(1),map_coord(2))
	3282	bdims (3,1) = f2c_dims_cg (4,map_coord(1),map_coord(2))
	3283	bdims (3,2) = f2c_dims_cg (5,map_coord(1),map_coord(2))
	3284
	3285
	3286	n_cell_c = (bdims(1,2)-bdims(1,1)+1) * &
	3287	(bdims(2,2)-bdims(2,1)+1) * &
	3288	(bdims(3,2)-bdims(3,1)+0)
	3289
	3290
	3291	CALL MPI_RECV( e( bdims(3,1)+1:bdims(3,2), &
	3292	bdims(2,1) :bdims(2,2), &
	3293	bdims(1,1) :bdims(1,2)),&
	3294	n_cell_c, MPI_REAL, target_idex, 104, &
	3295	comm_inter,status, ierr )
	3296	end do
	3297	end do
	3298
	3299
	3300	!
	3301	!-- Boundary conditions
	3302
	3303	IF ( .NOT. constant_diffusion ) THEN
	3304	e(nzb,:,:) = e(nzb+1,:,:)
	3305	END IF
	3306
	3307	IF ( .NOT. constant_diffusion ) CALL exchange_horiz( e, nbgp )
	3308
	3309	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	3310	! Send data to coarse grid for anterpolation
	3311
	3312	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	3313	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	3314	map_coord(1) = offset(1)
	3315	map_coord(2) = offset(2)
	3316	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	3317
	3318	bdims_rem (1,1) = f2c_dims_fg (0)
	3319	bdims_rem (1,2) = f2c_dims_fg (1)
	3320	bdims_rem (2,1) = f2c_dims_fg (2)
	3321	bdims_rem (2,2) = f2c_dims_fg (3)
	3322	bdims_rem (3,1) = f2c_dims_fg (4)
	3323	bdims_rem (3,2) = f2c_dims_fg (5)
	3324
	3325	ALLOCATE( work3d ( &
	3326	bdims_rem(3,1)+1:bdims_rem(3,2), &
	3327	bdims_rem(2,1) :bdims_rem(2,2), &
	3328	bdims_rem(1,1) :bdims_rem(1,2)))
	3329
	3330	anterpol3d => e
	3331
[3241]	3332	CALL anterpolate_to_crse_e
[2365]	3333
	3334	CALL MPI_SEND( work3d, 1, TYPE_VNEST_ANTER, target_idex, &
[2514]	3335	104, comm_inter, ierr)
[2365]	3336
	3337	NULLIFY ( anterpol3d )
	3338	DEALLOCATE( work3d )
	3339	ENDIF
	3340
	3341
	3342	CONTAINS
	3343
	3344
	3345
	3346
	3347
[3241]	3348	SUBROUTINE anterpolate_to_crse_e
[2365]	3349
	3350
	3351	USE arrays_3d
	3352	USE control_parameters
	3353	USE grid_variables
	3354	USE indices
	3355	USE pegrid
	3356
	3357
	3358	IMPLICIT NONE
[2712]	3359
	3360	INTEGER(iwp) :: i
	3361	INTEGER(iwp) :: j
	3362	INTEGER(iwp) :: k
	3363	INTEGER(iwp) :: iif
	3364	INTEGER(iwp) :: jjf
	3365	INTEGER(iwp) :: kkf
	3366	INTEGER(iwp) :: bottomx
	3367	INTEGER(iwp) :: bottomy
	3368	INTEGER(iwp) :: bottomz
	3369	INTEGER(iwp) :: topx
	3370	INTEGER(iwp) :: topy
	3371	INTEGER(iwp) :: topz
	3372	REAL(wp) :: aweight_a
	3373	REAL(wp) :: aweight_b
	3374	REAL(wp) :: aweight_c
	3375	REAL(wp) :: aweight_d
	3376	REAL(wp) :: aweight_e
	3377	REAL(wp) :: energ
[2365]	3378
	3379	DO k = bdims_rem(3,1)+1, bdims_rem(3,2)
	3380
	3381	bottomz = (dzc/dzf) * (k-1) + 1
	3382	topz = (dzc/dzf) * k
	3383
	3384	DO j = bdims_rem(2,1), bdims_rem(2,2)
	3385
	3386	bottomy = (nyf+1) / (nyc+1) * j
	3387	topy = (nyf+1) / (nyc+1) * (j+1) - 1
	3388
	3389	DO i = bdims_rem(1,1), bdims_rem(1,2)
	3390
	3391	bottomx = (nxf+1) / (nxc+1) * i
	3392	topx = (nxf+1) / (nxc+1) * (i+1) - 1
	3393
	3394	aweight_a = 0.0
	3395	aweight_b = 0.0
	3396	aweight_c = 0.0
	3397	aweight_d = 0.0
	3398	aweight_e = 0.0
	3399
[2712]	3400	DO kkf = bottomz, topz
	3401	DO jjf = bottomy, topy
	3402	DO iif = bottomx, topx
[2365]	3403
[2712]	3404	aweight_a = aweight_a + anterpol3d(kkf,jjf,iif) * &
[2365]	3405	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3406
	3407
[2712]	3408	energ = ( 0.5 * ( u(kkf,jjf,iif) + u(kkf,jjf,iif+1) ) )**2.0 + &
	3409	( 0.5 * ( v(kkf,jjf,iif) + v(kkf,jjf+1,iif) ) )**2.0 + &
	3410	( 0.5 * ( w(kkf-1,jjf,iif) + w(kkf,jjf,iif) ) )**2.0
[2365]	3411
	3412	aweight_b = aweight_b + energ * &
	3413	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3414
[2712]	3415	aweight_c = aweight_c + 0.5 * ( u(kkf,jjf,iif) + u(kkf,jjf,iif+1) ) * &
[2365]	3416	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3417
[2712]	3418	aweight_d = aweight_d + 0.5 * ( v(kkf,jjf,iif) + v(kkf,jjf+1,iif) ) * &
[2365]	3419	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3420
[2712]	3421	aweight_e = aweight_e + 0.5 * ( w(kkf-1,jjf,iif) + w(kkf,jjf,iif) ) * &
[2365]	3422	(dzf/dzc) * (dyf/dyc) * (dxf/dxc)
	3423
	3424
	3425	END DO
	3426	END DO
	3427	END DO
	3428
	3429	work3d(k,j,i) = aweight_a + 0.5 * ( aweight_b - &
	3430	aweight_c**2.0 - &
	3431	aweight_d**2.0 - &
	3432	aweight_e**2.0 )
	3433
	3434	END DO
	3435
	3436	END DO
	3437
	3438	END DO
	3439
	3440
	3441
	3442	END SUBROUTINE anterpolate_to_crse_e
[2712]	3443	#endif
[2365]	3444	END SUBROUTINE vnest_anterpolate_e
	3445
	3446	SUBROUTINE vnest_init_pegrid_rank
[2712]	3447	#if defined( __parallel )
[2365]	3448	! Domain decomposition and exchange of grid variables between coarse and fine
	3449	! Given processor coordinates as index f_rnk_lst(pcoord(1), pcoord(2))
	3450	! returns the rank. A single coarse block will have to send data to multiple
	3451	! fine blocks. In the coarse grid the pcoords of the remote block is first found and then using
	3452	! f_rnk_lst the target_idex is identified.
	3453	! blk_dim stores the index limits of a given block. blk_dim_remote is received
	3454	! from the asscoiated nest partner.
	3455	! cf_ratio(1:3) is the ratio between fine and coarse grid: nxc/nxf, nyc/nyf and
	3456	! ceiling(dxc/dxf)
	3457
	3458
[3241]	3459	USE control_parameters, &
	3460	ONLY: coupling_mode
[2365]	3461
	3462	USE kinds
	3463
	3464	USE pegrid
	3465
	3466
	3467	IMPLICIT NONE
	3468
[2712]	3469	INTEGER(iwp) :: dest_rnk
	3470	INTEGER(iwp) :: i !<
[2365]	3471
	3472	IF (myid == 0) THEN
	3473	IF ( coupling_mode == 'vnested_crse') THEN
	3474	CALL MPI_SEND( pdims, 2, MPI_INTEGER, numprocs, 33, comm_inter, &
	3475	ierr )
[3241]	3476	CALL MPI_RECV( pdims_partner, 2, MPI_INTEGER, numprocs, 66, &
[2365]	3477	comm_inter, status, ierr )
	3478	ELSEIF ( coupling_mode == 'vnested_fine') THEN
	3479	CALL MPI_RECV( pdims_partner, 2, MPI_INTEGER, 0, 33, &
	3480	comm_inter, status, ierr )
	3481	CALL MPI_SEND( pdims, 2, MPI_INTEGER, 0, 66, comm_inter, &
	3482	ierr )
	3483	ENDIF
	3484	ENDIF
	3485
	3486
	3487	IF ( coupling_mode == 'vnested_crse') THEN
	3488	CALL MPI_BCAST( pdims_partner, 2, MPI_INTEGER, 0, comm2d, ierr )
	3489	ALLOCATE( c_rnk_lst( 0:(pdims(1)-1) ,0:(pdims(2)-1) ) )
	3490	ALLOCATE( f_rnk_lst( 0:(pdims_partner(1)-1) ,0:(pdims_partner(2)-1) ) )
	3491	do i=0,numprocs-1
	3492	CALL MPI_CART_COORDS( comm2d, i, ndim, pcoord, ierr )
	3493	call MPI_Cart_rank(comm2d, pcoord, dest_rnk, ierr)
	3494	c_rnk_lst(pcoord(1),pcoord(2)) = dest_rnk
	3495	end do
	3496	ELSEIF ( coupling_mode == 'vnested_fine') THEN
	3497	CALL MPI_BCAST( pdims_partner, 2, MPI_INTEGER, 0, comm2d, ierr )
	3498	ALLOCATE( c_rnk_lst( 0:(pdims_partner(1)-1) ,0:(pdims_partner(2)-1) ) )
	3499	ALLOCATE( f_rnk_lst( 0:(pdims(1)-1) ,0:(pdims(2)-1) ) )
	3500
	3501	do i=0,numprocs-1
	3502	CALL MPI_CART_COORDS( comm2d, i, ndim, pcoord, ierr )
	3503	call MPI_Cart_rank(comm2d, pcoord, dest_rnk, ierr)
	3504	f_rnk_lst(pcoord(1),pcoord(2)) = dest_rnk
	3505	enddo
	3506	ENDIF
	3507
	3508
	3509	IF ( coupling_mode == 'vnested_crse') THEN
	3510	if (myid == 0) then
	3511	CALL MPI_SEND( c_rnk_lst, pdims(1)*pdims(2), MPI_INTEGER, numprocs, 0, comm_inter, ierr )
	3512	CALL MPI_RECV( f_rnk_lst, pdims_partner(1)*pdims_partner(2), MPI_INTEGER, numprocs, 4, comm_inter,status, ierr )
	3513	end if
	3514	CALL MPI_BCAST( f_rnk_lst, pdims_partner(1)*pdims_partner(2), MPI_INTEGER, 0, comm2d, ierr )
	3515	ELSEIF ( coupling_mode == 'vnested_fine') THEN
	3516	if (myid == 0) then
	3517	CALL MPI_RECV( c_rnk_lst, pdims_partner(1)*pdims_partner(2), MPI_INTEGER, 0, 0, comm_inter,status, ierr )
	3518	CALL MPI_SEND( f_rnk_lst, pdims(1)*pdims(2), MPI_INTEGER, 0, 4, comm_inter, ierr )
	3519	end if
	3520	CALL MPI_BCAST( c_rnk_lst, pdims_partner(1)*pdims_partner(2), MPI_INTEGER, 0, comm2d, ierr )
	3521	ENDIF
	3522
	3523	!-- Reason for MPI error unknown; solved if three lines duplicated
	3524	CALL MPI_CART_COORDS( comm2d, myid, ndim, pcoord, ierr )
	3525	CALL MPI_CART_SHIFT( comm2d, 0, 1, pleft, pright, ierr )
	3526	CALL MPI_CART_SHIFT( comm2d, 1, 1, psouth, pnorth, ierr )
	3527
	3528
[2712]	3529	#endif
[2365]	3530
	3531	END SUBROUTINE vnest_init_pegrid_rank
	3532
	3533
	3534	SUBROUTINE vnest_init_pegrid_domain
[2712]	3535	#if defined( __parallel )
[2365]	3536
[3241]	3537	USE control_parameters, &
	3538	ONLY: coupling_mode, coupling_topology, dz, &
[3065]	3539	dz_stretch_level_start, message_string
[2365]	3540
[3241]	3541	USE grid_variables, &
[2365]	3542	ONLY: dx, dy
	3543
[3241]	3544	USE indices, &
	3545	ONLY: nbgp, nx, ny, nz, nxl, nxr, nys, nyn, nzb, nzt
[2365]	3546
	3547	USE kinds
	3548
	3549	USE pegrid
	3550
	3551	IMPLICIT NONE
	3552
[2712]	3553	INTEGER(iwp) :: i !<
	3554	INTEGER(iwp) :: j !<
	3555	INTEGER(iwp) :: tempx
	3556	INTEGER(iwp) :: tempy
	3557	INTEGER(iwp) :: TYPE_INT_YZ
	3558	INTEGER(iwp) :: SIZEOFREAL
	3559	INTEGER(iwp) :: MTV_X
	3560	INTEGER(iwp) :: MTV_Y
	3561	INTEGER(iwp) :: MTV_Z
	3562	INTEGER(iwp) :: MTV_RX
	3563	INTEGER(iwp) :: MTV_RY
	3564	INTEGER(iwp) :: MTV_RZ
[2365]	3565
	3566	!
	3567	!-- Pass the number of grid points of the coarse model to
	3568	!-- the nested model and vice versa
	3569	IF ( coupling_mode == 'vnested_crse' ) THEN
	3570
	3571	nxc = nx
	3572	nyc = ny
	3573	nzc = nz
	3574	dxc = dx
	3575	dyc = dy
[3065]	3576	dzc = dz(1)
[2365]	3577	cg_nprocs = numprocs
	3578
	3579	IF ( myid == 0 ) THEN
	3580
[3065]	3581	CALL MPI_SEND( nxc, 1, MPI_INTEGER , numprocs, 1, comm_inter, &
[2365]	3582	ierr )
[3065]	3583	CALL MPI_SEND( nyc, 1, MPI_INTEGER , numprocs, 2, comm_inter, &
[2365]	3584	ierr )
[3065]	3585	CALL MPI_SEND( nzc, 1, MPI_INTEGER , numprocs, 3, comm_inter, &
[2365]	3586	ierr )
[3065]	3587	CALL MPI_SEND( dxc, 1, MPI_REAL , numprocs, 4, comm_inter, &
[2365]	3588	ierr )
[3065]	3589	CALL MPI_SEND( dyc, 1, MPI_REAL , numprocs, 5, comm_inter, &
[2365]	3590	ierr )
[3065]	3591	CALL MPI_SEND( dzc, 1, MPI_REAL , numprocs, 6, comm_inter, &
[2365]	3592	ierr )
[3065]	3593	CALL MPI_SEND( pdims, 2, MPI_INTEGER, numprocs, 7, comm_inter, &
[2365]	3594	ierr )
	3595	CALL MPI_SEND( cg_nprocs, 1, MPI_INTEGER, numprocs, 8, comm_inter, &
	3596	ierr )
[3065]	3597	CALL MPI_RECV( nxf, 1, MPI_INTEGER, numprocs, 21, comm_inter, &
[2365]	3598	status, ierr )
[3065]	3599	CALL MPI_RECV( nyf, 1, MPI_INTEGER, numprocs, 22, comm_inter, &
[2365]	3600	status, ierr )
[3065]	3601	CALL MPI_RECV( nzf, 1, MPI_INTEGER, numprocs, 23, comm_inter, &
[2365]	3602	status, ierr )
[3065]	3603	CALL MPI_RECV( dxf, 1, MPI_REAL, numprocs, 24, comm_inter, &
[2365]	3604	status, ierr )
[3065]	3605	CALL MPI_RECV( dyf, 1, MPI_REAL, numprocs, 25, comm_inter, &
[2365]	3606	status, ierr )
[3065]	3607	CALL MPI_RECV( dzf, 1, MPI_REAL, numprocs, 26, comm_inter, &
[2365]	3608	status, ierr )
[3065]	3609	CALL MPI_RECV( pdims_partner, 2, MPI_INTEGER, &
[2365]	3610	numprocs, 27, comm_inter, status, ierr )
[3065]	3611	CALL MPI_RECV( fg_nprocs, 1, MPI_INTEGER, &
[2365]	3612	numprocs, 28, comm_inter, status, ierr )
	3613	ENDIF
	3614
	3615	CALL MPI_BCAST( nxf, 1, MPI_INTEGER, 0, comm2d, ierr )
	3616	CALL MPI_BCAST( nyf, 1, MPI_INTEGER, 0, comm2d, ierr )
	3617	CALL MPI_BCAST( nzf, 1, MPI_INTEGER, 0, comm2d, ierr )
	3618	CALL MPI_BCAST( dxf, 1, MPI_REAL, 0, comm2d, ierr )
	3619	CALL MPI_BCAST( dyf, 1, MPI_REAL, 0, comm2d, ierr )
	3620	CALL MPI_BCAST( dzf, 1, MPI_REAL, 0, comm2d, ierr )
	3621	CALL MPI_BCAST( pdims_partner, 2, MPI_INTEGER, 0, comm2d, ierr )
	3622	CALL MPI_BCAST( fg_nprocs, 1, MPI_INTEGER, 0, comm2d, ierr )
[3065]	3623
	3624	!
	3625	!-- Check if stretching is used within the nested domain. ABS(...) is
	3626	!-- necessary because of the default value of -9999999.9_wp (negative)
	3627	IF ( ABS( dz_stretch_level_start(1) ) <= (nzf+1)*dzf ) THEN
	3628	message_string = 'Stretching in the parent domain is '// &
	3629	'only allowed above the nested domain'
[3066]	3630	CALL message( 'vnest_init_pegrid_domain', 'PA0497', 1, 2, 0, 6, 0 )
[3065]	3631	ENDIF
[2365]	3632
	3633	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	3634
	3635	nxf = nx
	3636	nyf = ny
	3637	nzf = nz
	3638	dxf = dx
	3639	dyf = dy
[3065]	3640	dzf = dz(1)
[2365]	3641	fg_nprocs = numprocs
	3642
	3643	IF ( myid == 0 ) THEN
	3644
	3645	CALL MPI_RECV( nxc, 1, MPI_INTEGER, 0, 1, comm_inter, status, &
	3646	ierr )
	3647	CALL MPI_RECV( nyc, 1, MPI_INTEGER, 0, 2, comm_inter, status, &
	3648	ierr )
	3649	CALL MPI_RECV( nzc, 1, MPI_INTEGER, 0, 3, comm_inter, status, &
	3650	ierr )
	3651	CALL MPI_RECV( dxc, 1, MPI_REAL, 0, 4, comm_inter, status, &
	3652	ierr )
	3653	CALL MPI_RECV( dyc, 1, MPI_REAL, 0, 5, comm_inter, status, &
	3654	ierr )
	3655	CALL MPI_RECV( dzc, 1, MPI_REAL, 0, 6, comm_inter, status, &
	3656	ierr )
	3657	CALL MPI_RECV( pdims_partner, 2, MPI_INTEGER, 0, 7, comm_inter, &
	3658	status, ierr )
	3659	CALL MPI_RECV( cg_nprocs, 1, MPI_INTEGER, 0, 8, comm_inter, &
	3660	status, ierr )
	3661	CALL MPI_SEND( nxf, 1, MPI_INTEGER, 0, 21, comm_inter, ierr )
	3662	CALL MPI_SEND( nyf, 1, MPI_INTEGER, 0, 22, comm_inter, ierr )
	3663	CALL MPI_SEND( nzf, 1, MPI_INTEGER, 0, 23, comm_inter, ierr )
	3664	CALL MPI_SEND( dxf, 1, MPI_REAL, 0, 24, comm_inter, ierr )
	3665	CALL MPI_SEND( dyf, 1, MPI_REAL, 0, 25, comm_inter, ierr )
	3666	CALL MPI_SEND( dzf, 1, MPI_REAL, 0, 26, comm_inter, ierr )
	3667	CALL MPI_SEND( pdims,2,MPI_INTEGER, 0, 27, comm_inter, ierr )
	3668	CALL MPI_SEND( fg_nprocs,1,MPI_INTEGER, 0, 28, comm_inter, ierr )
	3669	ENDIF
	3670
	3671	CALL MPI_BCAST( nxc, 1, MPI_INTEGER, 0, comm2d, ierr)
	3672	CALL MPI_BCAST( nyc, 1, MPI_INTEGER, 0, comm2d, ierr)
	3673	CALL MPI_BCAST( nzc, 1, MPI_INTEGER, 0, comm2d, ierr)
	3674	CALL MPI_BCAST( dxc, 1, MPI_REAL, 0, comm2d, ierr)
	3675	CALL MPI_BCAST( dyc, 1, MPI_REAL, 0, comm2d, ierr)
	3676	CALL MPI_BCAST( dzc, 1, MPI_REAL, 0, comm2d, ierr)
	3677	CALL MPI_BCAST( pdims_partner, 2, MPI_INTEGER, 0, comm2d, ierr)
	3678	CALL MPI_BCAST( cg_nprocs, 1, MPI_INTEGER, 0, comm2d, ierr)
	3679
	3680	ENDIF
[3065]	3681
[2365]	3682	ngp_c = ( nxc+1 + 2 * nbgp ) * ( nyc+1 + 2 * nbgp )
	3683	ngp_f = ( nxf+1 + 2 * nbgp ) * ( nyf+1 + 2 * nbgp )
	3684
	3685	IF ( coupling_mode(1:8) == 'vnested_') coupling_topology = 1
	3686
	3687
	3688	!-- Nesting Ratio: For each coarse grid cell how many fine grid cells exist
	3689	cfratio(1) = INT ( (nxf+1) / (nxc+1) )
	3690	cfratio(2) = INT ( (nyf+1) / (nyc+1) )
	3691	cfratio(3) = CEILING ( dzc / dzf )
	3692
	3693	!-- target_id is used only for exhange of information like simulated_time
	3694	!-- which are then MPI_BCAST to other processors in the group
	3695	IF ( myid == 0 ) THEN
	3696
	3697	IF ( TRIM( coupling_mode ) == 'vnested_crse' ) THEN
	3698	target_id = numprocs
	3699	ELSE IF ( TRIM( coupling_mode ) == 'vnested_fine' ) THEN
	3700	target_id = 0
	3701	ENDIF
	3702
	3703	ENDIF
	3704
	3705	!-- Store partner grid dimenstions and create MPI derived types
	3706
	3707	IF ( coupling_mode == 'vnested_crse' ) THEN
	3708
	3709	offset(1) = ( pdims_partner(1) / pdims(1) ) * pcoord(1)
	3710	offset(2) = ( pdims_partner(2) / pdims(2) ) * pcoord(2)
	3711
[2514]	3712	tempx = ( pdims_partner(1) / pdims(1) ) - 1
	3713	tempy = ( pdims_partner(2) / pdims(2) ) - 1
[2365]	3714	ALLOCATE( c2f_dims_cg (0:5,offset(1):tempx+offset(1),offset(2):tempy+offset(2) ) )
	3715	ALLOCATE( f2c_dims_cg (0:5,offset(1):tempx+offset(1),offset(2):tempy+offset(2) ) )
	3716
	3717	do j = 0, ( pdims_partner(2) / pdims(2) ) - 1
	3718	do i = 0, ( pdims_partner(1) / pdims(1) ) - 1
	3719	map_coord(1) = i+offset(1)
	3720	map_coord(2) = j+offset(2)
	3721
	3722	target_idex = f_rnk_lst(map_coord(1),map_coord(2)) + numprocs
	3723
	3724	CALL MPI_RECV( bdims_rem, 6, MPI_INTEGER, target_idex, 10, &
	3725	comm_inter,status, ierr )
	3726
	3727	!-- Store the CG dimensions that correspond to the FG partner; needed for FG top BC
	3728	!-- One CG can have multiple FG partners. The 3D array is mapped by partner proc co-ord
	3729	c2f_dims_cg (0,map_coord(1),map_coord(2)) = bdims_rem (1,1) / cfratio(1)
	3730	c2f_dims_cg (1,map_coord(1),map_coord(2)) = bdims_rem (1,2) / cfratio(1)
	3731	c2f_dims_cg (2,map_coord(1),map_coord(2)) = bdims_rem (2,1) / cfratio(2)
	3732	c2f_dims_cg (3,map_coord(1),map_coord(2)) = bdims_rem (2,2) / cfratio(2)
	3733	c2f_dims_cg (4,map_coord(1),map_coord(2)) = bdims_rem (3,2) / cfratio(3)
	3734	c2f_dims_cg (5,map_coord(1),map_coord(2)) =(bdims_rem (3,2) / cfratio(3)) + 2
	3735
	3736	!-- Store the CG dimensions that correspond to the FG partner; needed for anterpolation
	3737	f2c_dims_cg (0,map_coord(1),map_coord(2)) = bdims_rem (1,1) / cfratio(1)
	3738	f2c_dims_cg (1,map_coord(1),map_coord(2)) = bdims_rem (1,2) / cfratio(1)
	3739	f2c_dims_cg (2,map_coord(1),map_coord(2)) = bdims_rem (2,1) / cfratio(2)
	3740	f2c_dims_cg (3,map_coord(1),map_coord(2)) = bdims_rem (2,2) / cfratio(2)
	3741	f2c_dims_cg (4,map_coord(1),map_coord(2)) = bdims_rem (3,1)
	3742	f2c_dims_cg (5,map_coord(1),map_coord(2)) =(bdims_rem (3,2)-cfratio(3))/ cfratio(3)
	3743
	3744	CALL MPI_SEND( c2f_dims_cg (:,map_coord(1),map_coord(2)), 6, &
[2514]	3745	MPI_INTEGER, target_idex, 100, comm_inter, ierr )
[2365]	3746
	3747	CALL MPI_SEND( f2c_dims_cg (:,map_coord(1),map_coord(2)), 6, &
[2514]	3748	MPI_INTEGER, target_idex, 101, comm_inter, ierr )
[2365]	3749
	3750	end do
	3751	end do
	3752
	3753	!-- A derived data type to pack 3 Z-levels of CG to set FG top BC
	3754	MTV_X = ( nxr - nxl + 1 ) + 2*nbgp
	3755	MTV_Y = ( nyn - nys + 1 ) + 2*nbgp
	3756	MTV_Z = nzt+1 - nzb +1
	3757
	3758	MTV_RX = ( c2f_dims_cg (1,offset(1),offset(2)) - c2f_dims_cg (0,offset(1),offset(2)) ) +1+2
	3759	MTV_RY = ( c2f_dims_cg (3,offset(1),offset(2)) - c2f_dims_cg (2,offset(1),offset(2)) ) +1+2
	3760	MTV_RZ = ( c2f_dims_cg (5,offset(1),offset(2)) - c2f_dims_cg (4,offset(1),offset(2)) ) +1
	3761
	3762	CALL MPI_TYPE_EXTENT(MPI_REAL, SIZEOFREAL, IERR)
	3763
	3764	CALL MPI_TYPE_VECTOR ( MTV_RY, MTV_RZ, MTV_Z, MPI_REAL, TYPE_INT_YZ, IERR)
	3765	CALL MPI_TYPE_HVECTOR( MTV_RX, 1, MTV_ZMTV_YSIZEOFREAL, &
[2514]	3766	TYPE_INT_YZ, TYPE_VNEST_BC, IERR)
[2365]	3767	CALL MPI_TYPE_FREE(TYPE_INT_YZ, IERR)
[2514]	3768	CALL MPI_TYPE_COMMIT(TYPE_VNEST_BC, IERR)
[2365]	3769
	3770
	3771	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	3772
	3773	ALLOCATE( c2f_dims_fg (0:5) )
	3774	ALLOCATE( f2c_dims_fg (0:5) )
	3775
	3776	offset(1) = pcoord(1) / ( pdims(1)/pdims_partner(1) )
	3777	offset(2) = pcoord(2) / ( pdims(2)/pdims_partner(2) )
	3778	map_coord(1) = offset(1)
	3779	map_coord(2) = offset(2)
	3780	target_idex = c_rnk_lst(map_coord(1),map_coord(2))
	3781
	3782	bdims (1,1) = nxl
	3783	bdims (1,2) = nxr
	3784	bdims (2,1) = nys
	3785	bdims (2,2) = nyn
	3786	bdims (3,1) = nzb
	3787	bdims (3,2) = nzt
	3788
	3789	CALL MPI_SEND( bdims, 6, MPI_INTEGER, target_idex, 10, &
	3790	comm_inter, ierr )
	3791
	3792	!-- Store the CG dimensions that correspond to the FG partner; needed for FG top BC
	3793	!-- One FG can have only one CG partner
	3794	CALL MPI_RECV( c2f_dims_fg, 6, MPI_INTEGER, target_idex, 100, &
	3795	comm_inter,status, ierr )
	3796
	3797	CALL MPI_RECV( f2c_dims_fg, 6, MPI_INTEGER, target_idex, 101, &
	3798	comm_inter,status, ierr )
	3799
	3800	!-- Store the CG dimensions that correspond to the FG partner; needed for anterpolation
	3801
	3802	n_cell_c = (f2c_dims_fg(1)-f2c_dims_fg(0)+1) * &
	3803	(f2c_dims_fg(3)-f2c_dims_fg(2)+1) * &
	3804	(f2c_dims_fg(5)-f2c_dims_fg(4)+0)
	3805
	3806	CALL MPI_TYPE_CONTIGUOUS(n_cell_c, MPI_REAL, TYPE_VNEST_ANTER, IERR)
	3807	CALL MPI_TYPE_COMMIT(TYPE_VNEST_ANTER, ierr)
	3808
	3809	ENDIF
[2712]	3810	#endif
[2365]	3811	END SUBROUTINE vnest_init_pegrid_domain
	3812
	3813
	3814	SUBROUTINE vnest_init_grid
	3815
[2712]	3816	#if defined( __parallel )
[3065]	3817	USE arrays_3d, &
[2365]	3818	ONLY: zu, zw
	3819
[3065]	3820	USE control_parameters, &
	3821	ONLY: coupling_mode, message_string, number_stretch_level_start
[2365]	3822
[3065]	3823	USE indices, &
[2365]	3824	ONLY: nzt
	3825
	3826	USE kinds
	3827
	3828	USE pegrid
	3829
[2514]	3830	IMPLICIT NONE
[2365]	3831
[3065]	3832	!
	3833	!-- Allocate and Exchange zuc and zuf, zwc and zwf
[2365]	3834	IF ( coupling_mode(1:8) == 'vnested_' ) THEN
	3835
	3836	ALLOCATE( zuc(0:nzc+1), zuf(0:nzf+1) )
	3837	ALLOCATE( zwc(0:nzc+1), zwf(0:nzf+1) )
	3838
	3839	IF ( coupling_mode == 'vnested_crse' ) THEN
[3065]	3840
	3841	zuc = zu
	3842	zwc = zw
	3843
[2365]	3844	IF ( myid == 0 ) THEN
	3845
	3846	CALL MPI_SEND( zuc, nzt+2, MPI_REAL, numprocs, 41, comm_inter, &
	3847	ierr )
	3848	CALL MPI_RECV( zuf, nzf+2, MPI_REAL, numprocs, 42, comm_inter, &
	3849	status, ierr )
	3850
	3851	CALL MPI_SEND( zwc, nzt+2, MPI_REAL, numprocs, 43, comm_inter, &
	3852	ierr )
	3853	CALL MPI_RECV( zwf, nzf+2, MPI_REAL, numprocs, 44, comm_inter, &
	3854	status, ierr )
	3855
	3856	ENDIF
	3857
	3858	CALL MPI_BCAST( zuf,nzf+2,MPI_REAL, 0, comm2d, ierr )
	3859	CALL MPI_BCAST( zwf,nzf+2,MPI_REAL, 0, comm2d, ierr )
	3860
	3861	ELSEIF ( coupling_mode == 'vnested_fine' ) THEN
	3862
[3065]	3863	!
	3864	!-- Check if stretching is used within the nested domain
	3865	IF ( number_stretch_level_start > 0 ) THEN
	3866	message_string = 'Stretching in the nested domain is not '//&
	3867	'allowed'
[3066]	3868	CALL message( 'vnest_init_grid', 'PA0498', 1, 2, 0, 6, 0 )
[3065]	3869	ENDIF
	3870
	3871	zuf = zu
	3872	zwf = zw
	3873
[2365]	3874	IF ( myid == 0 ) THEN
	3875
	3876	CALL MPI_RECV( zuc,nzc+2, MPI_REAL, 0, 41, comm_inter, status, &
	3877	ierr )
	3878	CALL MPI_SEND( zuf,nzt+2, MPI_REAL, 0, 42, comm_inter, ierr )
	3879
	3880	CALL MPI_RECV( zwc,nzc+2, MPI_REAL, 0, 43, comm_inter, status, &
	3881	ierr )
	3882	CALL MPI_SEND( zwf,nzt+2, MPI_REAL, 0, 44, comm_inter, ierr )
	3883	ENDIF
	3884
	3885	CALL MPI_BCAST( zuc,nzc+2,MPI_REAL, 0, comm2d, ierr )
	3886	CALL MPI_BCAST( zwc,nzc+2,MPI_REAL, 0, comm2d, ierr )
	3887
	3888	ENDIF
	3889	ENDIF
	3890
[2712]	3891	#endif
[2365]	3892	END SUBROUTINE vnest_init_grid
	3893
	3894
	3895	SUBROUTINE vnest_check_parameters
[2712]	3896	#if defined( __parallel )
[2365]	3897
[2712]	3898	USE pegrid, &
	3899	ONLY: myid
[2365]	3900
[2514]	3901	IMPLICIT NONE
[2365]	3902
[2514]	3903	IF (myid==0) PRINT, ' vnest: check parameters not implemented yet *'
[2365]	3904
[2712]	3905	#endif
[2365]	3906	END SUBROUTINE vnest_check_parameters
	3907
	3908
	3909	SUBROUTINE vnest_timestep_sync
	3910
[2712]	3911	#if defined( __parallel )
[2365]	3912	USE control_parameters, &
[4101]	3913	ONLY: coupling_mode, dt_3d
[2365]	3914
	3915	USE interfaces
	3916
	3917	USE kinds
	3918
	3919	USE pegrid
	3920
	3921	IMPLICIT NONE
	3922
	3923	IF ( coupling_mode == 'vnested_crse') THEN
[2514]	3924	dtc = dt_3d
	3925	if (myid == 0) then
[2365]	3926	CALL MPI_SEND( dt_3d, 1, MPI_REAL, target_id, &
	3927	31, comm_inter, ierr )
	3928	CALL MPI_RECV( dtf, 1, MPI_REAL, &
	3929	target_id, 32, comm_inter, status, ierr )
	3930
[2514]	3931	endif
	3932	CALL MPI_BCAST( dtf, 1, MPI_REAL, 0, comm2d, ierr )
[2365]	3933	ELSE
[2514]	3934	dtf = dt_3d
	3935	if (myid == 0) then
[2365]	3936	CALL MPI_RECV( dtc, 1, MPI_REAL, &
	3937	target_id, 31, comm_inter, status, ierr )
	3938	CALL MPI_SEND( dt_3d, 1, MPI_REAL, target_id, &
	3939	32, comm_inter, ierr )
	3940
[2514]	3941	endif
	3942	CALL MPI_BCAST( dtc, 1, MPI_REAL, 0, comm2d, ierr )
[2365]	3943
	3944	ENDIF
	3945	!-- Identical timestep for coarse and fine grids
	3946	dt_3d = MIN( dtc, dtf )
[2712]	3947	#endif
[2365]	3948	END SUBROUTINE vnest_timestep_sync
	3949
	3950	SUBROUTINE vnest_deallocate
[2712]	3951	#if defined( __parallel )
[2365]	3952	USE control_parameters, &
	3953	ONLY: coupling_mode
	3954
	3955	IMPLICIT NONE
	3956
	3957	IF ( ALLOCATED(c_rnk_lst) ) DEALLOCATE (c_rnk_lst)
	3958	IF ( ALLOCATED(f_rnk_lst) ) DEALLOCATE (f_rnk_lst)
	3959
	3960	IF ( coupling_mode == 'vnested_crse') THEN
	3961	IF ( ALLOCATED (c2f_dims_cg) ) DEALLOCATE (c2f_dims_cg)
	3962	IF ( ALLOCATED (f2c_dims_cg) ) DEALLOCATE (f2c_dims_cg)
	3963	ELSEIF( coupling_mode == 'vnested_fine' ) THEN
	3964	IF ( ALLOCATED (c2f_dims_fg) ) DEALLOCATE (c2f_dims_fg)
	3965	IF ( ALLOCATED (f2c_dims_fg) ) DEALLOCATE (f2c_dims_fg)
	3966	ENDIF
[2712]	3967	#endif
[2365]	3968	END SUBROUTINE vnest_deallocate
	3969
[3802]	3970	END MODULE vertical_nesting_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |