Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |