Home

Context Navigation

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

Last change on this file since 3049 was 3049, checked in by Giersch, 6 years ago
Revision history corrected
Property svn:keywords set to `Id`
File size: 153.8 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |