Home

Context Navigation

pmc_interface_mod.f90 @ 3945

Last change on this file since 3945 was 3945, checked in by raasch, 6 years ago
messed document changes for r3932 cleaned up
Property svn:keywords set to `Id`
File size: 204.1 KB

Line
1	!> @file pmc_interface_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-2019 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! ------------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: pmc_interface_mod.f90 3945 2019-05-02 11:29:27Z raasch $
27	!
28	! 3932 2019-04-24 17:31:34Z suehring
29	! Add missing if statements for call of pmc_set_dataarray_name for TKE and
30	! dissipation.
31	!
32	! 3888 2019-04-12 09:18:10Z hellstea
33	! Variables renamed, commenting improved etc.
34	!
35	! 3885 2019-04-11 11:29:34Z kanani
36	! Changes related to global restructuring of location messages and introduction
37	! of additional debug messages
38	!
39	! 3883 2019-04-10 12:51:50Z hellstea
40	! Checks and error messages improved and extended. All the child index bounds in the
41	! parent-grid index space are made module variables. Function get_number_of_childs
42	! renamed get_number_of_children. A number of variables renamed
43	! and qite a lot of other code reshaping made all around the module.
44	!
45	! 3876 2019-04-08 18:41:49Z knoop
46	! Implemented nesting for salsa variables.
47	!
48	! 3833 2019-03-28 15:04:04Z forkel
49	! replaced USE chem_modules by USE chem_gasphase_mod
50	!
51	! 3822 2019-03-27 13:10:23Z hellstea
52	! Temporary increase of the vertical dimension of the parent-grid arrays and
53	! workarrc_t is cancelled as unnecessary.
54	!
55	! 3819 2019-03-27 11:01:36Z hellstea
56	! Adjustable anterpolation buffer introduced on all nest boundaries, it is controlled
57	! by the new nesting_parameters parameter anterpolation_buffer_width.
58	!
59	! 3804 2019-03-19 13:46:20Z hellstea
60	! Anterpolation domain is lowered from kct-1 to kct-3 to avoid exessive
61	! kinetic energy from building up in CBL flows.
62	!
63	! 3803 2019-03-19 13:44:40Z hellstea
64	! A bug fixed in lateral boundary interpolations. Dimension of val changed from
65	! 5 to 3 in pmci_setup_parent and pmci_setup_child.
66	!
67	! 3794 2019-03-15 09:36:33Z raasch
68	! two remaining unused variables removed
69	!
70	! 3792 2019-03-14 16:50:07Z hellstea
71	! Interpolations improved. Large number of obsolete subroutines removed.
72	! All unused variables removed.
73	!
74	! 3741 2019-02-13 16:24:49Z hellstea
75	! Interpolations and child initialization adjusted to handle set ups with child
76	! pe-subdomain dimension not integer divisible by the grid-spacing ratio in the
77	! respective direction. Set ups with pe-subdomain dimension smaller than the
78	! grid-spacing ratio in the respective direction are now forbidden.
79	!
80	! 3708 2019-01-30 12:58:13Z hellstea
81	! Checks for parent / child grid line matching introduced.
82	! Interpolation of nest-boundary-tangential velocity components revised.
83	!
84	! 3697 2019-01-24 17:16:13Z hellstea
85	! Bugfix: upper k-bound in the child initialization interpolation
86	! pmci_interp_1sto_all corrected.
87	! Copying of the nest boundary values into the redundant 2nd and 3rd ghost-node
88	! layers is added to the pmci_interp_1sto_*-routines.
89	!
90	! 3681 2019-01-18 15:06:05Z hellstea
91	! Linear interpolations are replaced by first order interpolations. The linear
92	! interpolation routines are still included but not called. In the child
93	! inititialization the interpolation is also changed to 1st order and the linear
94	! interpolation is not kept.
95	! Subroutine pmci_map_fine_to_coarse_grid is rewritten.
96	! Several changes in pmci_init_anterp_tophat.
97	! Child's parent-grid arrays (uc, vc,...) are made non-overlapping on the PE-
98	! subdomain boundaries in order to allow grid-spacing ratios higher than nbgp.
99	! Subroutine pmci_init_tkefactor is removed as unnecessary.
100	!
101	! 3655 2019-01-07 16:51:22Z knoop
102	! Remove unused variable simulated_time
103	!
104	! 3636 2018-12-19 13:48:34Z raasch
105	! nopointer option removed
106	!
107	! 3592 2018-12-03 12:38:40Z suehring
108	! Number of coupled arrays is determined dynamically (instead of a fixed value
109	! of 32)
110	!
111	! 3524 2018-11-14 13:36:44Z raasch
112	! declaration statements rearranged to avoid compile time errors
113	!
114	! 3484 2018-11-02 14:41:25Z hellstea
115	! Introduction of reversibility correction to the interpolation routines in order to
116	! guarantee mass and scalar conservation through the nest boundaries. Several errors
117	! are corrected and pmci_ensure_nest_mass_conservation is permanently removed.
118	!
119	! 3274 2018-09-24 15:42:55Z knoop
120	! Modularization of all bulk cloud physics code components
121	!
122	! 3241 2018-09-12 15:02:00Z raasch
123	! unused variables removed
124	!
125	! 3217 2018-08-29 12:53:59Z suehring
126	! Revise calculation of index bounds for array index_list, prevent compiler
127	! (cray) to delete the loop at high optimization level.
128	!
129	! 3215 2018-08-29 09:58:59Z suehring
130	! Apply an additional switch controlling the nesting of chemical species
131	!
132	! 3209 2018-08-27 16:58:37Z suehring
133	! Variable names for nest_bound_x replaced by bc_dirichlet_x.
134	! Remove commented prints into debug files.
135	!
136	! 3182 2018-07-27 13:36:03Z suehring
137	! dz was replaced by dz(1)
138	!
139	! 3049 2018-05-29 13:52:36Z Giersch
140	! Error messages revised
141	!
142	! 3045 2018-05-28 07:55:41Z Giersch
143	! Error messages revised
144	!
145	! 3022 2018-05-18 11:12:35Z suehring
146	! Minor fix - working precision added to real number
147	!
148	! 3021 2018-05-16 08:14:20Z maronga
149	! Bugfix: variable lcr was defined as INTENT(OUT) instead of INTENT(INOUT)
150	!
151	! 3020 2018-05-14 10:45:48Z hellstea
152	! Bugfix in pmci_define_loglaw_correction_parameters
153	!
154	! 3001 2018-04-20 12:27:13Z suehring
155	! Bugfix, replace MERGE function by an if-condition in the anterpolation (in
156	! order to avoid floating-point exceptions).
157	!
158	! 2997 2018-04-19 13:35:17Z suehring
159	! Mask topography grid points in anterpolation
160	!
161	! 2984 2018-04-18 11:51:30Z hellstea
162	! Bugfix in the log-law correction initialization. Pivot node cannot any more be
163	! selected from outside the subdomain in order to prevent array under/overflows.
164	!
165	! 2967 2018-04-13 11:22:08Z raasch
166	! bugfix: missing parallel cpp-directives added
167	!
168	! 2951 2018-04-06 09:05:08Z kanani
169	! Add log_point_s for pmci_model_configuration
170	!
171	! 2938 2018-03-27 15:52:42Z suehring
172	! - Nesting for RANS mode implemented
173	! - Interpolation of TKE onto child domain only if both parent and child are
174	! either in LES mode or in RANS mode
175	! - Interpolation of dissipation if both parent and child are in RANS mode
176	! and TKE-epsilon closure is applied
177	! - Enable anterpolation of TKE and dissipation rate in case parent and
178	! child operate in RANS mode
179	!
180	! - Some unused variables removed from ONLY list
181	! - Some formatting adjustments for particle nesting
182	!
183	! 2936 2018-03-27 14:49:27Z suehring
184	! Control logics improved to allow nesting also in cases with
185	! constant_flux_layer = .F. or constant_diffusion = .T.
186	!
187	! 2895 2018-03-15 10:26:12Z hellstea
188	! Change in the nest initialization (pmci_interp_tril_all). Bottom wall BC is no
189	! longer overwritten.
190	!
191	! 2868 2018-03-09 13:25:09Z hellstea
192	! Local conditional Neumann conditions for one-way coupling removed.
193	!
194	! 2853 2018-03-05 14:44:20Z suehring
195	! Bugfix in init log-law correction.
196	!
197	! 2841 2018-02-27 15:02:57Z knoop
198	! Bugfix: wrong placement of include 'mpif.h' corrected
199	!
200	! 2812 2018-02-16 13:40:25Z hellstea
201	! Bugfixes in computation of the interpolation loglaw-correction parameters
202	!
203	! 2809 2018-02-15 09:55:58Z schwenkel
204	! Bugfix for gfortran: Replace the function C_SIZEOF with STORAGE_SIZE
205	!
206	! 2806 2018-02-14 17:10:37Z thiele
207	! Bugfixing Write statements
208	!
209	! 2804 2018-02-14 16:57:03Z thiele
210	! preprocessor directive for c_sizeof in case of __gfortran added
211	!
212	! 2803 2018-02-14 16:56:32Z thiele
213	! Introduce particle transfer in nested models.
214	!
215	! 2795 2018-02-07 14:48:48Z hellstea
216	! Bugfix in computation of the anterpolation under-relaxation functions.
217	!
218	! 2773 2018-01-30 14:12:54Z suehring
219	! - Nesting for chemical species
220	! - Bugfix in setting boundary condition at downward-facing walls for passive
221	! scalar
222	! - Some formatting adjustments
223	!
224	! 2718 2018-01-02 08:49:38Z maronga
225	! Corrected "Former revisions" section
226	!
227	! 2701 2017-12-15 15:40:50Z suehring
228	! Changes from last commit documented
229	!
230	! 2698 2017-12-14 18:46:24Z suehring
231	! Bugfix in get_topography_top_index
232	!
233	! 2696 2017-12-14 17:12:51Z kanani
234	! Change in file header (GPL part)
235	! - Bugfix in init_tke_factor (MS)
236	!
237	! 2669 2017-12-06 16:03:27Z raasch
238	! file extension for nested domains changed to "_N##",
239	! created flag file in order to enable combine_plot_fields to process nest data
240	!
241	! 2663 2017-12-04 17:40:50Z suehring
242	! Bugfix, wrong coarse-grid index in init_tkefactor used.
243	!
244	! 2602 2017-11-03 11:06:41Z hellstea
245	! Index-limit related bug (occurred with nesting_mode='vertical') fixed in
246	! pmci_interp_tril_t. Check for too high nest domain added in pmci_setup_parent.
247	! Some cleaning up made.
248	!
249	! 2582 2017-10-26 13:19:46Z hellstea
250	! Resetting of e within buildings / topography in pmci_parent_datatrans removed
251	! as unnecessary since e is not anterpolated, and as incorrect since it overran
252	! the default Neumann condition (bc_e_b).
253	!
254	! 2359 2017-08-21 07:50:30Z hellstea
255	! A minor indexing error in pmci_init_loglaw_correction is corrected.
256	!
257	! 2351 2017-08-15 12:03:06Z kanani
258	! Removed check (PA0420) for nopointer version
259	!
260	! 2350 2017-08-15 11:48:26Z kanani
261	! Bugfix and error message for nopointer version.
262	!
263	! 2318 2017-07-20 17:27:44Z suehring
264	! Get topography top index via Function call
265	!
266	! 2317 2017-07-20 17:27:19Z suehring
267	! Set bottom boundary condition after anterpolation.
268	! Some variable description added.
269	!
270	! 2293 2017-06-22 12:59:12Z suehring
271	! In anterpolation, exclude grid points which are used for interpolation.
272	! This avoids the accumulation of numerical errors leading to increased
273	! variances for shallow child domains.
274	!
275	! 2292 2017-06-20 09:51:42Z schwenkel
276	! Implementation of new microphysic scheme: cloud_scheme = 'morrison'
277	! includes two more prognostic equations for cloud drop concentration (nc)
278	! and cloud water content (qc).
279	!
280	! 2285 2017-06-15 13:15:41Z suehring
281	! Consider multiple pure-vertical nest domains in overlap check
282	!
283	! 2283 2017-06-14 10:17:34Z suehring
284	! Bugfixes in inititalization of log-law correction concerning
285	! new topography concept
286	!
287	! 2281 2017-06-13 11:34:50Z suehring
288	! Bugfix, add pre-preprocessor directives to enable non-parrallel mode
289	!
290	! 2241 2017-06-01 13:46:13Z hellstea
291	! A minor indexing error in pmci_init_loglaw_correction is corrected.
292	!
293	! 2240 2017-06-01 13:45:34Z hellstea
294	!
295	! 2232 2017-05-30 17:47:52Z suehring
296	! Adjustments to new topography concept
297	!
298	! 2229 2017-05-30 14:52:52Z hellstea
299	! A minor indexing error in pmci_init_anterp_tophat is corrected.
300	!
301	! 2174 2017-03-13 08:18:57Z maronga
302	! Added support for cloud physics quantities, syntax layout improvements. Data
303	! transfer of qc and nc is prepared but currently deactivated until both
304	! quantities become prognostic variables.
305	! Some bugfixes.
306	!
307	! 2019 2016-09-30 13:40:09Z hellstea
308	! Bugfixes mainly in determining the anterpolation index bounds. These errors
309	! were detected when first time tested using 3:1 grid-spacing.
310	!
311	! 2003 2016-08-24 10:22:32Z suehring
312	! Humidity and passive scalar also separated in nesting mode
313	!
314	! 2000 2016-08-20 18:09:15Z knoop
315	! Forced header and separation lines into 80 columns
316	!
317	! 1938 2016-06-13 15:26:05Z hellstea
318	! Minor clean-up.
319	!
320	! 1901 2016-05-04 15:39:38Z raasch
321	! Initial version of purely vertical nesting introduced.
322	! Code clean up. The words server/client changed to parent/child.
323	!
324	! 1900 2016-05-04 15:27:53Z raasch
325	! unused variables removed
326	!
327	! 1894 2016-04-27 09:01:48Z raasch
328	! bugfix: pt interpolations are omitted in case that the temperature equation is
329	! switched off
330	!
331	! 1892 2016-04-26 13:49:47Z raasch
332	! bugfix: pt is not set as a data array in case that the temperature equation is
333	! switched off with neutral = .TRUE.
334	!
335	! 1882 2016-04-20 15:24:46Z hellstea
336	! The factor ijfc for nfc used in anterpolation is redefined as 2-D array
337	! and precomputed in pmci_init_anterp_tophat.
338	!
339	! 1878 2016-04-19 12:30:36Z hellstea
340	! Synchronization rewritten, logc-array index order changed for cache
341	! optimization
342	!
343	! 1850 2016-04-08 13:29:27Z maronga
344	! Module renamed
345	!
346	!
347	! 1808 2016-04-05 19:44:00Z raasch
348	! MPI module used by default on all machines
349	!
350	! 1801 2016-04-05 13:12:47Z raasch
351	! bugfix for r1797: zero setting of temperature within topography does not work
352	! and has been disabled
353	!
354	! 1797 2016-03-21 16:50:28Z raasch
355	! introduction of different datatransfer modes,
356	! further formatting cleanup, parameter checks added (including mismatches
357	! between root and nest model settings),
358	! +routine pmci_check_setting_mismatches
359	! comm_world_nesting introduced
360	!
361	! 1791 2016-03-11 10:41:25Z raasch
362	! routine pmci_update_new removed,
363	! pmc_get_local_model_info renamed pmc_get_model_info, some keywords also
364	! renamed,
365	! filling up redundant ghost points introduced,
366	! some index bound variables renamed,
367	! further formatting cleanup
368	!
369	! 1783 2016-03-06 18:36:17Z raasch
370	! calculation of nest top area simplified,
371	! interpolation and anterpolation moved to seperate wrapper subroutines
372	!
373	! 1781 2016-03-03 15:12:23Z raasch
374	! _p arrays are set zero within buildings too, t.._m arrays and respective
375	! settings within buildings completely removed
376	!
377	! 1779 2016-03-03 08:01:28Z raasch
378	! only the total number of PEs is given for the domains, npe_x and npe_y
379	! replaced by npe_total, two unused elements removed from array
380	! parent_grid_info_real,
381	! array management changed from linked list to sequential loop
382	!
383	! 1766 2016-02-29 08:37:15Z raasch
384	! modifications to allow for using PALM's pointer version,
385	! +new routine pmci_set_swaplevel
386	!
387	! 1764 2016-02-28 12:45:19Z raasch
388	! +cpl_parent_id,
389	! cpp-statements for nesting replaced by __parallel statements,
390	! errors output with message-subroutine,
391	! index bugfixes in pmci_interp_tril_all,
392	! some adjustments to PALM style
393	!
394	! 1762 2016-02-25 12:31:13Z hellstea
395	! Initial revision by A. Hellsten
396	!
397	! Description:
398	! ------------
399	! Domain nesting interface routines. The low-level inter-domain communication
400	! is conducted by the PMC-library routines.
401	!
402	! @todo Remove array_3d variables from USE statements thate not used in the
403	! routine
404	! @todo Data transfer of qc and nc is prepared but not activated
405	!------------------------------------------------------------------------------!
406	MODULE pmc_interface
407
408	USE ISO_C_BINDING
409
410
411	USE arrays_3d, &
412	ONLY: diss, diss_2, dzu, dzw, e, e_p, e_2, nc, nc_2, nc_p, nr, nr_2, &
413	pt, pt_2, q, q_2, qc, qc_2, qr, qr_2, s, s_2, &
414	u, u_p, u_2, v, v_p, v_2, w, w_p, w_2, zu, zw
415
416	USE control_parameters, &
417	ONLY: air_chemistry, bc_dirichlet_l, bc_dirichlet_n, bc_dirichlet_r, &
418	bc_dirichlet_s, child_domain, &
419	constant_diffusion, constant_flux_layer, &
420	coupling_char, &
421	debug_output, debug_string, &
422	dt_3d, dz, humidity, message_string, &
423	neutral, passive_scalar, rans_mode, rans_tke_e, &
424	roughness_length, salsa, topography, volume_flow
425
426	USE chem_gasphase_mod, &
427	ONLY: nspec
428
429	USE chem_modules, &
430	ONLY: chem_species
431
432	USE chemistry_model_mod, &
433	ONLY: nest_chemistry, spec_conc_2
434
435	USE cpulog, &
436	ONLY: cpu_log, log_point_s
437
438	USE grid_variables, &
439	ONLY: dx, dy
440
441	USE indices, &
442	ONLY: nbgp, nx, nxl, nxlg, nxlu, nxr, nxrg, ny, nyn, nyng, nys, nysg, &
443	nysv, nz, nzb, nzt, wall_flags_0
444
445	USE bulk_cloud_model_mod, &
446	ONLY: bulk_cloud_model, microphysics_morrison, microphysics_seifert
447
448	USE particle_attributes, &
449	ONLY: particle_advection
450
451	USE kinds
452
453	#if defined( __parallel )
454	#if !defined( __mpifh )
455	USE MPI
456	#endif
457
458	USE pegrid, &
459	ONLY: collective_wait, comm1dx, comm1dy, comm2d, myid, myidx, myidy, &
460	numprocs, pdims, pleft, pnorth, pright, psouth, status
461
462	USE pmc_child, &
463	ONLY: pmc_childinit, pmc_c_clear_next_array_list, &
464	pmc_c_getnextarray, pmc_c_get_2d_index_list, pmc_c_getbuffer, &
465	pmc_c_putbuffer, pmc_c_setind_and_allocmem, &
466	pmc_c_set_dataarray, pmc_set_dataarray_name
467
468	USE pmc_general, &
469	ONLY: da_namelen
470
471	USE pmc_handle_communicator, &
472	ONLY: pmc_get_model_info, pmc_init_model, pmc_is_rootmodel, &
473	pmc_no_namelist_found, pmc_parent_for_child, m_couplers
474
475	USE pmc_mpi_wrapper, &
476	ONLY: pmc_bcast, pmc_recv_from_child, pmc_recv_from_parent, &
477	pmc_send_to_child, pmc_send_to_parent
478
479	USE pmc_parent, &
480	ONLY: pmc_parentinit, pmc_s_clear_next_array_list, pmc_s_fillbuffer, &
481	pmc_s_getdata_from_buffer, pmc_s_getnextarray, &
482	pmc_s_setind_and_allocmem, pmc_s_set_active_data_array, &
483	pmc_s_set_dataarray, pmc_s_set_2d_index_list
484
485	#endif
486
487	USE salsa_mod, &
488	ONLY: aerosol_mass, aerosol_number, gconc_2, mconc_2, nbins_aerosol, &
489	ncomponents_mass, nconc_2, nest_salsa, ngases_salsa, salsa_gas, &
490	salsa_gases_from_chem
491
492	USE surface_mod, &
493	ONLY: get_topography_top_index_ji, surf_def_h, surf_lsm_h, surf_usm_h
494
495	IMPLICIT NONE
496
497	#if defined( __parallel )
498	#if defined( __mpifh )
499	INCLUDE "mpif.h"
500	#endif
501	#endif
502
503	PRIVATE
504	!
505	!-- Constants
506	INTEGER(iwp), PARAMETER :: child_to_parent = 2 !<
507	INTEGER(iwp), PARAMETER :: parent_to_child = 1 !<
508	INTEGER(iwp), PARAMETER :: interpolation_scheme_lrsn = 2 !< Interpolation scheme to be used on lateral boundaries
509	INTEGER(iwp), PARAMETER :: interpolation_scheme_t = 3 !< Interpolation scheme to be used on top boundary
510	!
511	!-- Coupler setup
512	INTEGER(iwp), SAVE :: comm_world_nesting !<
513	INTEGER(iwp), SAVE :: cpl_id = 1 !<
514	INTEGER(iwp), SAVE :: cpl_npe_total !<
515	INTEGER(iwp), SAVE :: cpl_parent_id !<
516
517	CHARACTER(LEN=32), SAVE :: cpl_name !<
518
519	!
520	!-- Control parameters
521	INTEGER(iwp), SAVE :: anterpolation_buffer_width = 2 !< Boundary buffer width for anterpolation
522	CHARACTER(LEN=7), SAVE :: nesting_datatransfer_mode = 'mixed' !< steering parameter for data-transfer mode
523	CHARACTER(LEN=8), SAVE :: nesting_mode = 'two-way' !< steering parameter for 1- or 2-way nesting
524
525	LOGICAL, SAVE :: nested_run = .FALSE. !< general switch
526	LOGICAL :: rans_mode_parent = .FALSE. !< mode of parent model (.F. - LES mode, .T. - RANS mode)
527	!
528	!-- Geometry
529	REAL(wp), SAVE, DIMENSION(:), ALLOCATABLE, PUBLIC :: coord_x !<
530	REAL(wp), SAVE, DIMENSION(:), ALLOCATABLE, PUBLIC :: coord_y !<
531	REAL(wp), SAVE, PUBLIC :: lower_left_coord_x !<
532	REAL(wp), SAVE, PUBLIC :: lower_left_coord_y !<
533	!
534	!-- Children's parent-grid arrays
535	INTEGER(iwp), SAVE, DIMENSION(5), PUBLIC :: coarse_bound !< subdomain index bounds for children's parent-grid arrays
536
537	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: dissc !< Parent-grid array on child domain - dissipation rate
538	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ec !< Parent-grid array on child domain - SGS TKE
539	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ptc !< Parent-grid array on child domain - potential temperature
540	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uc !< Parent-grid array on child domain - velocity component u
541	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vc !< Parent-grid array on child domain - velocity component v
542	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wc !< Parent-grid array on child domain - velocity component w
543	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: q_c !< Parent-grid array on child domain -
544	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: qcc !< Parent-grid array on child domain -
545	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: qrc !< Parent-grid array on child domain -
546	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: nrc !< Parent-grid array on child domain -
547	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ncc !< Parent-grid array on child domain -
548	REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: sc !< Parent-grid array on child domain -
549	INTEGER(idp), SAVE, DIMENSION(:,:), ALLOCATABLE, TARGET, PUBLIC :: nr_partc !<
550	INTEGER(idp), SAVE, DIMENSION(:,:), ALLOCATABLE, TARGET, PUBLIC :: part_adrc !<
551
552	REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: chem_spec_c !< Parent-grid array on child domain - chemical species
553
554	REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_mass_c !< Aerosol mass
555	REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_number_c !< Aerosol number
556	REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: salsa_gas_c !< SALSA gases
557	!
558	!-- Grid-spacing ratios.
559	INTEGER(iwp), SAVE :: igsr !< Integer grid-spacing ratio in i-direction
560	INTEGER(iwp), SAVE :: jgsr !< Integer grid-spacing ratio in j-direction
561	INTEGER(iwp), SAVE :: kgsr !< Integer grid-spacing ratio in k-direction
562	!
563	!-- Global parent-grid index bounds
564	INTEGER(iwp), SAVE :: iplg !< Leftmost parent-grid array ip index of the whole child domain
565	INTEGER(iwp), SAVE :: iprg !< Rightmost parent-grid array ip index of the whole child domain
566	INTEGER(iwp), SAVE :: jpsg !< Southmost parent-grid array jp index of the whole child domain
567	INTEGER(iwp), SAVE :: jpng !< Northmost parent-grid array jp index of the whole child domain
568	!
569	!-- Local parent-grid index bounds. Different sets of index bounds are needed for parent-grid arrays (uc, etc),
570	!-- for index mapping arrays (iflu, etc) and for work arrays (workarr_lr, etc). This is because these arrays
571	!-- have different dimensions depending on the location of the subdomain relative to boundaries and corners.
572	INTEGER(iwp), SAVE :: ipl !< Left index limit for children's parent-grid arrays
573	INTEGER(iwp), SAVE :: ipla !< Left index limit for allocation of index-mapping and other auxiliary arrays
574	INTEGER(iwp), SAVE :: iplw !< Left index limit for children's parent-grid work arrays
575	INTEGER(iwp), SAVE :: ipr !< Right index limit for children's parent-grid arrays
576	INTEGER(iwp), SAVE :: ipra !< Right index limit for allocation of index-mapping and other auxiliary arrays
577	INTEGER(iwp), SAVE :: iprw !< Right index limit for children's parent-grid work arrays
578	INTEGER(iwp), SAVE :: jpn !< North index limit for children's parent-grid arrays
579	INTEGER(iwp), SAVE :: jpna !< North index limit for allocation of index-mapping and other auxiliary arrays
580	INTEGER(iwp), SAVE :: jpnw !< North index limit for children's parent-grid work arrays
581	INTEGER(iwp), SAVE :: jps !< South index limit for children's parent-grid arrays
582	INTEGER(iwp), SAVE :: jpsa !< South index limit for allocation of index-mapping and other auxiliary arrays
583	INTEGER(iwp), SAVE :: jpsw !< South index limit for children's parent-grid work arrays
584	!
585	!-- Highest prognostic parent-grid k-indices.
586	INTEGER(iwp), SAVE :: kcto !< Upper bound for k in anterpolation of variables other than w.
587	INTEGER(iwp), SAVE :: kctw !< Upper bound for k in anterpolation of w.
588	!
589	!-- Child-array indices to be precomputed and stored for anterpolation.
590	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflu !< child index indicating left bound of parent grid box on u-grid
591	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuu !< child index indicating right bound of parent grid box on u-grid
592	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflo !< child index indicating left bound of parent grid box on scalar-grid
593	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuo !< child index indicating right bound of parent grid box on scalar-grid
594	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflv !< child index indicating south bound of parent grid box on v-grid
595	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuv !< child index indicating north bound of parent grid box on v-grid
596	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflo !< child index indicating south bound of parent grid box on scalar-grid
597	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuo !< child index indicating north bound of parent grid box on scalar-grid
598	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflw !< child index indicating lower bound of parent grid box on w-grid
599	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuw !< child index indicating upper bound of parent grid box on w-grid
600	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflo !< child index indicating lower bound of parent grid box on scalar-grid
601	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !< child index indicating upper bound of parent grid box on scalar-grid
602	!
603	!-- Number of child-grid nodes within anterpolation cells to be precomputed for anterpolation.
604	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_u !< number of child grid points contributing to a parent grid
605	!< node in anterpolation, u-grid
606	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_v !< number of child grid points contributing to a parent grid
607	!< node in anterpolation, v-grid
608	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_w !< number of child grid points contributing to a parent grid
609	!< node in anterpolation, w-grid
610	INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_s !< number of child grid points contributing to a parent grid
611	!< node in anterpolation, scalar-grid
612	!
613	!-- Work arrays for interpolation and user-defined type definitions for horizontal work-array exchange
614	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_lr
615	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_sn
616	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_t
617
618	INTEGER(iwp) :: workarr_lr_exchange_type
619	INTEGER(iwp) :: workarr_sn_exchange_type
620	INTEGER(iwp) :: workarr_t_exchange_type_x
621	INTEGER(iwp) :: workarr_t_exchange_type_y
622
623	INTEGER(iwp), DIMENSION(3) :: parent_grid_info_int !<
624
625	REAL(wp), DIMENSION(7) :: parent_grid_info_real !<
626	REAL(wp), DIMENSION(2) :: zmax_coarse !<
627
628	TYPE parentgrid_def
629	INTEGER(iwp) :: nx !<
630	INTEGER(iwp) :: ny !<
631	INTEGER(iwp) :: nz !<
632	REAL(wp) :: dx !<
633	REAL(wp) :: dy !<
634	REAL(wp) :: dz !<
635	REAL(wp) :: lower_left_coord_x !<
636	REAL(wp) :: lower_left_coord_y !<
637	REAL(wp) :: xend !<
638	REAL(wp) :: yend !<
639	REAL(wp), DIMENSION(:), ALLOCATABLE :: coord_x !<
640	REAL(wp), DIMENSION(:), ALLOCATABLE :: coord_y !<
641	REAL(wp), DIMENSION(:), ALLOCATABLE :: dzu !<
642	REAL(wp), DIMENSION(:), ALLOCATABLE :: dzw !<
643	REAL(wp), DIMENSION(:), ALLOCATABLE :: zu !<
644	REAL(wp), DIMENSION(:), ALLOCATABLE :: zw !<
645	END TYPE parentgrid_def
646
647	TYPE(parentgrid_def), SAVE, PUBLIC :: pg !< Parent-grid information package of type parentgrid_def
648	!
649	!-- Variables for particle coupling
650	TYPE, PUBLIC :: childgrid_def
651	INTEGER(iwp) :: nx !<
652	INTEGER(iwp) :: ny !<
653	INTEGER(iwp) :: nz !<
654	REAL(wp) :: dx !<
655	REAL(wp) :: dy !<
656	REAL(wp) :: dz !<
657	REAL(wp) :: lx_coord, lx_coord_b !< ! split onto separate lines
658	REAL(wp) :: rx_coord, rx_coord_b !<
659	REAL(wp) :: sy_coord, sy_coord_b !<
660	REAL(wp) :: ny_coord, ny_coord_b !<
661	REAL(wp) :: uz_coord, uz_coord_b !<
662	END TYPE childgrid_def
663
664	TYPE(childgrid_def), SAVE, ALLOCATABLE, DIMENSION(:), PUBLIC :: childgrid !<
665
666	INTEGER(idp), ALLOCATABLE,DIMENSION(:,:), PUBLIC,TARGET :: nr_part !<
667	INTEGER(idp), ALLOCATABLE,DIMENSION(:,:), PUBLIC,TARGET :: part_adr !<
668
669
670	INTERFACE pmci_boundary_conds
671	MODULE PROCEDURE pmci_boundary_conds
672	END INTERFACE pmci_boundary_conds
673
674	INTERFACE pmci_check_setting_mismatches
675	MODULE PROCEDURE pmci_check_setting_mismatches
676	END INTERFACE
677
678	INTERFACE pmci_child_initialize
679	MODULE PROCEDURE pmci_child_initialize
680	END INTERFACE
681
682	INTERFACE pmci_synchronize
683	MODULE PROCEDURE pmci_synchronize
684	END INTERFACE
685
686	INTERFACE pmci_datatrans
687	MODULE PROCEDURE pmci_datatrans
688	END INTERFACE pmci_datatrans
689
690	INTERFACE pmci_init
691	MODULE PROCEDURE pmci_init
692	END INTERFACE
693
694	INTERFACE pmci_modelconfiguration
695	MODULE PROCEDURE pmci_modelconfiguration
696	END INTERFACE
697
698	INTERFACE pmci_parent_initialize
699	MODULE PROCEDURE pmci_parent_initialize
700	END INTERFACE
701
702	INTERFACE get_number_of_children
703	MODULE PROCEDURE get_number_of_children
704	END INTERFACE get_number_of_children
705
706	INTERFACE get_childid
707	MODULE PROCEDURE get_childid
708	END INTERFACE get_childid
709
710	INTERFACE get_child_edges
711	MODULE PROCEDURE get_child_edges
712	END INTERFACE get_child_edges
713
714	INTERFACE get_child_gridspacing
715	MODULE PROCEDURE get_child_gridspacing
716	END INTERFACE get_child_gridspacing
717
718	INTERFACE pmci_set_swaplevel
719	MODULE PROCEDURE pmci_set_swaplevel
720	END INTERFACE pmci_set_swaplevel
721
722	PUBLIC child_to_parent, comm_world_nesting, cpl_id, nested_run, &
723	nesting_datatransfer_mode, nesting_mode, parent_to_child, rans_mode_parent
724
725	PUBLIC pmci_boundary_conds
726	PUBLIC pmci_child_initialize
727	PUBLIC pmci_datatrans
728	PUBLIC pmci_init
729	PUBLIC pmci_modelconfiguration
730	PUBLIC pmci_parent_initialize
731	PUBLIC pmci_synchronize
732	PUBLIC pmci_set_swaplevel
733	PUBLIC get_number_of_children, get_childid, get_child_edges, get_child_gridspacing
734
735
736	CONTAINS
737
738
739	SUBROUTINE pmci_init( world_comm )
740
741	USE control_parameters, &
742	ONLY: message_string
743
744	IMPLICIT NONE
745
746	INTEGER(iwp), INTENT(OUT) :: world_comm !<
747
748	#if defined( __parallel )
749
750	INTEGER(iwp) :: pmc_status !<
751
752
753	CALL pmc_init_model( world_comm, nesting_datatransfer_mode, nesting_mode, &
754	anterpolation_buffer_width, pmc_status )
755
756	IF ( pmc_status == pmc_no_namelist_found ) THEN
757	!
758	!-- This is not a nested run
759	world_comm = MPI_COMM_WORLD
760	cpl_id = 1
761	cpl_name = ""
762
763	RETURN
764
765	ENDIF
766	!
767	!-- Check steering parameter values
768	IF ( TRIM( nesting_mode ) /= 'one-way' .AND. &
769	TRIM( nesting_mode ) /= 'two-way' .AND. &
770	TRIM( nesting_mode ) /= 'vertical' ) &
771	THEN
772	message_string = 'illegal nesting mode: ' // TRIM( nesting_mode )
773	CALL message( 'pmci_init', 'PA0417', 3, 2, 0, 6, 0 )
774	ENDIF
775
776	IF ( TRIM( nesting_datatransfer_mode ) /= 'cascade' .AND. &
777	TRIM( nesting_datatransfer_mode ) /= 'mixed' .AND. &
778	TRIM( nesting_datatransfer_mode ) /= 'overlap' ) &
779	THEN
780	message_string = 'illegal nesting datatransfer mode: ' &
781	// TRIM( nesting_datatransfer_mode )
782	CALL message( 'pmci_init', 'PA0418', 3, 2, 0, 6, 0 )
783	ENDIF
784	!
785	!-- Set the general steering switch which tells PALM that its a nested run
786	nested_run = .TRUE.
787	!
788	!-- Get some variables required by the pmc-interface (and in some cases in the
789	!-- PALM code out of the pmci) out of the pmc-core
790	CALL pmc_get_model_info( comm_world_nesting = comm_world_nesting, &
791	cpl_id = cpl_id, cpl_parent_id = cpl_parent_id, &
792	cpl_name = cpl_name, npe_total = cpl_npe_total, &
793	lower_left_x = lower_left_coord_x, &
794	lower_left_y = lower_left_coord_y )
795	!
796	!-- Set the steering switch which tells the models that they are nested (of
797	!-- course the root domain (cpl_id = 1) is not nested)
798	IF ( cpl_id >= 2 ) THEN
799	child_domain = .TRUE.
800	WRITE( coupling_char, '(A2,I2.2)') '_N', cpl_id
801	ENDIF
802
803	!
804	!-- Message that communicators for nesting are initialized.
805	!-- Attention: myid has been set at the end of pmc_init_model in order to
806	!-- guarantee that only PE0 of the root domain does the output.
807	CALL location_message( 'initialize model nesting', 'finished' )
808	!
809	!-- Reset myid to its default value
810	myid = 0
811	#else
812	!
813	!-- Nesting cannot be used in serial mode. cpl_id is set to root domain (1)
814	!-- because no location messages would be generated otherwise.
815	!-- world_comm is given a dummy value to avoid compiler warnings (INTENT(OUT)
816	!-- should get an explicit value)
817	!-- todo: why don't we print an error message instead of these settings?
818	cpl_id = 1
819	nested_run = .FALSE.
820	world_comm = 1
821	#endif
822
823	END SUBROUTINE pmci_init
824
825
826
827	SUBROUTINE pmci_modelconfiguration
828
829	IMPLICIT NONE
830
831	INTEGER(iwp) :: ncpl !< number of nest domains
832
833
834	#if defined( __parallel )
835	CALL location_message( 'setup the nested model configuration', 'start' )
836	CALL cpu_log( log_point_s(79), 'pmci_model_config', 'start' )
837	!
838	!-- Compute absolute coordinates for all models
839	CALL pmci_setup_coordinates ! CONTAIN THIS
840	!
841	!-- Determine the number of coupled arrays
842	CALL pmci_num_arrays ! CONTAIN THIS
843	!
844	!-- Initialize the child (must be called before pmc_setup_parent)
845	! EXTEND THIS COMMENT EXPLAINEIN WHY IT MUST BE CALLED BEFORE
846	CALL pmci_setup_child ! CONTAIN THIS
847	!
848	!-- Initialize PMC parent
849	CALL pmci_setup_parent ! CONTAIN THIS
850	!
851	!-- Check for mismatches between settings of master and child variables
852	!-- (e.g., all children have to follow the end_time settings of the root master)
853	CALL pmci_check_setting_mismatches ! CONTAIN THIS
854	!
855	!-- Set flag file for combine_plot_fields for precessing the nest output data
856	OPEN( 90, FILE='3DNESTING', FORM='FORMATTED' )
857	CALL pmc_get_model_info( ncpl = ncpl )
858	WRITE( 90, '(I2)' ) ncpl
859	CLOSE( 90 )
860
861	CALL cpu_log( log_point_s(79), 'pmci_model_config', 'stop' )
862	CALL location_message( 'setup the nested model configuration', 'finished' )
863	#endif
864
865	END SUBROUTINE pmci_modelconfiguration
866
867
868
869	SUBROUTINE pmci_setup_parent
870
871	#if defined( __parallel )
872	IMPLICIT NONE
873
874	INTEGER(iwp) :: child_id !< Child id-number for the child m
875	INTEGER(iwp) :: ierr !< MPI-error code
876	INTEGER(iwp) :: kp !< Parent-grid index n the z-direction
877	INTEGER(iwp) :: lb = 1 !< Running index for aerosol size bins
878	INTEGER(iwp) :: lc = 1 !< Running index for aerosol mass bins
879	INTEGER(iwp) :: lg = 1 !< Running index for SALSA gases
880	INTEGER(iwp) :: m !< Loop index over all children of the current parent
881	INTEGER(iwp) :: msib !< Loop index over all other children than m in case of siblings (parallel children)
882	INTEGER(iwp) :: n = 1 !< Running index for chemical species
883	INTEGER(iwp) :: nest_overlap = 0 !< Tag for parallel child-domains' overlap situation (>0 if overlap found)
884	INTEGER(iwp) :: nomatch = 0 !< Tag for child-domain mismatch situation (>0 if mismatch found)
885	INTEGER(iwp) :: nx_child !< Number of child-grid points in the x-direction
886	INTEGER(iwp) :: ny_child !< Number of child-grid points in the y-direction
887	INTEGER(iwp) :: nz_child !< Number of child-grid points in the z-direction
888
889	INTEGER(iwp), DIMENSION(3) :: child_grid_dim !< Array for receiving the child-grid dimensions from the children
890
891	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_x_left !< Minimum x-coordinate of the child domain including the ghost
892	!< point layers
893	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_x_right !< Maximum x-coordinate of the child domain including the ghost
894	!< point layers
895	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_y_south !< Minimum y-coordinate of the child domain including the ghost
896	!< point layers
897	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_y_north !< Maximum y-coordinate of the child domain including the ghost
898	!< point layers
899	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_coord_x !< Child domain x-coordinate array
900	REAL(wp), DIMENSION(:), ALLOCATABLE :: child_coord_y !< Child domain y-coordinate array
901
902	REAL(wp), DIMENSION(5) :: child_grid_info !< Array for receiving the child-grid spacings etc from the children
903
904	REAL(wp) :: child_height !< Height of the child domain defined on the child side as zw(nzt+1)
905	REAL(wp) :: dx_child !< Child-grid spacing in the x-direction
906	REAL(wp) :: dy_child !< Child-grid spacing in the y-direction
907	REAL(wp) :: dz_child !< Child-grid spacing in the z-direction
908	REAL(wp) :: left_limit !< Left limit for the absolute x-coordinate of the child left boundary
909	REAL(wp) :: north_limit !< North limit for the absolute y-coordinate of the child north boundary
910	REAL(wp) :: right_limit !< Right limit for the absolute x-coordinate of the child right boundary
911	REAL(wp) :: south_limit !< South limit for the absolute y-coordinate of the child south boundary
912	REAL(wp) :: upper_right_coord_x !< Absolute x-coordinate of the upper right corner of the child domain
913	REAL(wp) :: upper_right_coord_y !< Absolute y-coordinate of the upper right corner of the child domain
914	REAL(wp) :: xez !< Minimum separation in the x-direction required between the child and
915	!< parent boundaries (left or right)
916	REAL(wp) :: yez !< Minimum separation in the y-direction required between the child and
917	!< parent boundaries (south or north)
918	CHARACTER(LEN=32) :: myname !< String for variable name such as 'u'
919
920	LOGICAL :: m_left_in_msib !< Logical auxiliary parameter for the overlap test: true if the left border
921	!< of the child m is within the x-range of the child msib
922	LOGICAL :: m_right_in_msib !< Logical auxiliary parameter for the overlap test: true if the right border
923	!< of the child m is within the x-range of the child msib
924	LOGICAL :: msib_left_in_m !< Logical auxiliary parameter for the overlap test: true if the left border
925	!< of the child msib is within the x-range of the child m
926	LOGICAL :: msib_right_in_m !< Logical auxiliary parameter for the overlap test: true if the right border
927	!< of the child msib is within the x-range of the child m
928	LOGICAL :: m_south_in_msib !< Logical auxiliary parameter for the overlap test: true if the south border
929	!< of the child m is within the y-range of the child msib
930	LOGICAL :: m_north_in_msib !< Logical auxiliary parameter for the overlap test: true if the north border
931	!< of the child m is within the y-range of the child msib
932	LOGICAL :: msib_south_in_m !< Logical auxiliary parameter for the overlap test: true if the south border
933	!< of the child msib is within the y-range of the child m
934	LOGICAL :: msib_north_in_m !< Logical auxiliary parameter for the overlap test: true if the north border
935	!< of the child msib is within the y-range of the child m
936	!
937	!-- Initialize the current pmc parent
938	CALL pmc_parentinit
939	!
940	!-- Corners of all children of the present parent. Note that
941	!-- SIZE( pmc_parent_for_child ) = 1 if we have no children.
942	IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN
943	ALLOCATE( child_x_left(1:SIZE( pmc_parent_for_child ) - 1) )
944	ALLOCATE( child_x_right(1:SIZE( pmc_parent_for_child ) - 1) )
945	ALLOCATE( child_y_south(1:SIZE( pmc_parent_for_child ) - 1) )
946	ALLOCATE( child_y_north(1:SIZE( pmc_parent_for_child ) - 1) )
947	ENDIF
948	IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) ) THEN
949	ALLOCATE( childgrid(1:SIZE( pmc_parent_for_child ) - 1) )
950	ENDIF
951	!
952	!-- Get coordinates from all children and check that the children match the parent
953	!-- domain and each others. Note that SIZE( pmc_parent_for_child ) = 1
954	!-- if we have no children, thence the loop is not executed at all.
955	DO m = 1, SIZE( pmc_parent_for_child ) - 1
956
957	child_id = pmc_parent_for_child(m)
958
959	IF ( myid == 0 ) THEN
960
961	CALL pmc_recv_from_child( child_id, child_grid_dim, SIZE(child_grid_dim), 0, 123, ierr )
962	CALL pmc_recv_from_child( child_id, child_grid_info, SIZE(child_grid_info), 0, 124, ierr )
963
964	nx_child = child_grid_dim(1)
965	ny_child = child_grid_dim(2)
966	dx_child = child_grid_info(3)
967	dy_child = child_grid_info(4)
968	dz_child = child_grid_info(5)
969	child_height = child_grid_info(1)
970	!
971	!-- Find the highest child-domain level in the parent grid for the reduced z
972	!-- transfer
973	DO kp = 1, nz
974	IF ( zw(kp) > child_height ) THEN
975	nz_child = kp
976	EXIT
977	ENDIF
978	ENDDO
979	zmax_coarse = child_grid_info(1:2)
980	!
981	!-- Get absolute coordinates from the child
982	ALLOCATE( child_coord_x(-nbgp:nx_child+nbgp) )
983	ALLOCATE( child_coord_y(-nbgp:ny_child+nbgp) )
984
985	CALL pmc_recv_from_child( child_id, child_coord_x, SIZE( child_coord_x ), 0, 11, ierr )
986	CALL pmc_recv_from_child( child_id, child_coord_y, SIZE( child_coord_y ), 0, 12, ierr )
987
988	parent_grid_info_real(1) = lower_left_coord_x
989	parent_grid_info_real(2) = lower_left_coord_y
990	parent_grid_info_real(3) = dx
991	parent_grid_info_real(4) = dy
992
993	upper_right_coord_x = lower_left_coord_x + ( nx + 1 ) * dx
994	upper_right_coord_y = lower_left_coord_y + ( ny + 1 ) * dy
995	parent_grid_info_real(5) = upper_right_coord_x
996	parent_grid_info_real(6) = upper_right_coord_y
997	parent_grid_info_real(7) = dz(1)
998
999	parent_grid_info_int(1) = nx
1000	parent_grid_info_int(2) = ny
1001	parent_grid_info_int(3) = nz_child
1002	!
1003	!-- Check that the child domain matches its parent domain.
1004	IF ( nesting_mode == 'vertical' ) THEN
1005	!
1006	!-- In case of vertical nesting, the lateral boundaries must match exactly.
1007	right_limit = upper_right_coord_x
1008	north_limit = upper_right_coord_y
1009	IF ( ( child_coord_x(nx_child+1) /= right_limit ) .OR. &
1010	( child_coord_y(ny_child+1) /= north_limit ) ) THEN
1011	nomatch = 1
1012	ENDIF
1013	ELSE
1014	!
1015	!-- In case of 3-D nesting, check that the child domain is completely
1016	!-- inside its parent domain.
1017	xez = ( nbgp + 1 ) * dx
1018	yez = ( nbgp + 1 ) * dy
1019	left_limit = lower_left_coord_x + xez
1020	right_limit = upper_right_coord_x - xez
1021	south_limit = lower_left_coord_y + yez
1022	north_limit = upper_right_coord_y - yez
1023	IF ( ( child_coord_x(0) < left_limit ) .OR. &
1024	( child_coord_x(nx_child+1) > right_limit ) .OR. &
1025	( child_coord_y(0) < south_limit ) .OR. &
1026	( child_coord_y(ny_child+1) > north_limit ) ) THEN
1027	nomatch = 1
1028	ENDIF
1029	ENDIF
1030	!
1031	!-- Child domain must be lower than the parent domain such
1032	!-- that the top ghost layer of the child grid does not exceed
1033	!-- the parent domain top boundary.
1034	IF ( child_height > zw(nz) ) THEN
1035	nomatch = 1
1036	ENDIF
1037	!
1038	!-- If parallel child domains (siblings) do exist ( m > 1 ),
1039	!-- check that they do not overlap.
1040	child_x_left(m) = child_coord_x(-nbgp)
1041	child_x_right(m) = child_coord_x(nx_child+nbgp)
1042	child_y_south(m) = child_coord_y(-nbgp)
1043	child_y_north(m) = child_coord_y(ny_child+nbgp)
1044
1045	IF ( nesting_mode /= 'vertical' ) THEN
1046	!
1047	!-- Note that the msib-loop is executed only if ( m > 1 ).
1048	!-- Also note that the tests have to be made both ways (m vs msib and msib vs m)
1049	!-- in order to detect all the possible overlap situations.
1050	DO msib = 1, m - 1
1051	!
1052	!-- Set some logical auxiliary parameters to simplify the IF-condition.
1053	m_left_in_msib = ( child_x_left(m) >= child_x_left(msib) ) .AND. &
1054	( child_x_left(m) <= child_x_right(msib) )
1055	m_right_in_msib = ( child_x_right(m) >= child_x_left(msib) ) .AND. &
1056	( child_x_right(m) <= child_x_right(msib) )
1057	msib_left_in_m = ( child_x_left(msib) >= child_x_left(m) ) .AND. &
1058	( child_x_left(msib) <= child_x_right(m) )
1059	msib_right_in_m = ( child_x_right(msib) >= child_x_left(m) ) .AND. &
1060	( child_x_right(msib) <= child_x_right(m) )
1061	m_south_in_msib = ( child_y_south(m) >= child_y_south(msib) ) .AND. &
1062	( child_y_south(m) <= child_y_north(msib) )
1063	m_north_in_msib = ( child_y_north(m) >= child_y_south(msib) ) .AND. &
1064	( child_y_north(m) <= child_y_north(msib) )
1065	msib_south_in_m = ( child_y_south(msib) >= child_y_south(m) ) .AND. &
1066	( child_y_south(msib) <= child_y_north(m) )
1067	msib_north_in_m = ( child_y_north(msib) >= child_y_south(m) ) .AND. &
1068	( child_y_north(msib) <= child_y_north(m) )
1069
1070	IF ( ( m_left_in_msib .OR. m_right_in_msib .OR. &
1071	msib_left_in_m .OR. msib_right_in_m ) &
1072	.AND. &
1073	( m_south_in_msib .OR. m_north_in_msib .OR. &
1074	msib_south_in_m .OR. msib_north_in_m ) ) THEN
1075	nest_overlap = 1
1076	ENDIF
1077
1078	ENDDO
1079	ENDIF
1080
1081	CALL set_child_edge_coords
1082
1083	DEALLOCATE( child_coord_x )
1084	DEALLOCATE( child_coord_y )
1085	!
1086	!-- Send information about operating mode (LES or RANS) to child. This will be
1087	!-- used to control TKE nesting and setting boundary conditions properly.
1088	CALL pmc_send_to_child( child_id, rans_mode, 1, 0, 19, ierr )
1089	!
1090	!-- Send coarse grid information to child
1091	CALL pmc_send_to_child( child_id, parent_grid_info_real, &
1092	SIZE( parent_grid_info_real ), 0, 21, &
1093	ierr )
1094	CALL pmc_send_to_child( child_id, parent_grid_info_int, 3, 0, &
1095	22, ierr )
1096	!
1097	!-- Send local grid to child
1098	CALL pmc_send_to_child( child_id, coord_x, nx+1+2*nbgp, 0, 24, &
1099	ierr )
1100	CALL pmc_send_to_child( child_id, coord_y, ny+1+2*nbgp, 0, 25, &
1101	ierr )
1102	!
1103	!-- Also send the dzu-, dzw-, zu- and zw-arrays here
1104	CALL pmc_send_to_child( child_id, dzu, nz_child + 1, 0, 26, ierr )
1105	CALL pmc_send_to_child( child_id, dzw, nz_child + 1, 0, 27, ierr )
1106	CALL pmc_send_to_child( child_id, zu, nz_child + 2, 0, 28, ierr )
1107	CALL pmc_send_to_child( child_id, zw, nz_child + 2, 0, 29, ierr )
1108
1109	IF ( nomatch /= 0 ) THEN
1110	WRITE ( message_string, * ) 'nested child domain does not fit into its parent domain'
1111	CALL message( 'pmci_setup_parent', 'PA0425', 3, 2, 0, 6, 0 )
1112	ENDIF
1113
1114	IF ( nest_overlap /= 0 .AND. nesting_mode /= 'vertical' ) THEN
1115	WRITE ( message_string, * ) 'nested parallel child domains overlap'
1116	CALL message( 'pmci_setup_parent', 'PA0426', 3, 2, 0, 6, 0 )
1117	ENDIF
1118
1119	ENDIF ! ( myid == 0 )
1120
1121	CALL MPI_BCAST( nz_child, 1, MPI_INTEGER, 0, comm2d, ierr )
1122
1123	CALL MPI_BCAST( childgrid(m), STORAGE_SIZE(childgrid(1))/8, MPI_BYTE, 0, comm2d, ierr )
1124	!
1125	!-- Set up the index-list which is an integer array that maps the child index space on
1126	!-- the parent index- and subdomain spaces.
1127	CALL pmci_create_index_list
1128	!
1129	!-- Include couple arrays into parent content
1130	!-- The adresses of the PALM 2D or 3D array (here server coarse grid) which are candidates
1131	!-- for coupling are stored once into the pmc context. While data transfer, the array do not
1132	!-- have to be specified again
1133	CALL pmc_s_clear_next_array_list
1134	DO WHILE ( pmc_s_getnextarray( child_id, myname ) )
1135	IF ( INDEX( TRIM( myname ), 'chem_' ) /= 0 ) THEN
1136	CALL pmci_set_array_pointer( myname, child_id = child_id, &
1137	nz_child = nz_child, n = n )
1138	n = n + 1
1139	ELSEIF ( INDEX( TRIM( myname ), 'an_' ) /= 0 ) THEN
1140	CALL pmci_set_array_pointer( myname, child_id = child_id, &
1141	nz_child = nz_child, n = lb )
1142	lb = lb + 1
1143	ELSEIF ( INDEX( TRIM( myname ), 'am_' ) /= 0 ) THEN
1144	CALL pmci_set_array_pointer( myname, child_id = child_id, &
1145	nz_child = nz_child, n = lc )
1146	lc = lc + 1
1147	ELSEIF ( INDEX( TRIM( myname ), 'sg_' ) /= 0 .AND. &
1148	.NOT. salsa_gases_from_chem ) &
1149	THEN
1150	CALL pmci_set_array_pointer( myname, child_id = child_id, &
1151	nz_child = nz_child, n = lg )
1152	lg = lg + 1
1153	ELSE
1154	CALL pmci_set_array_pointer( myname, child_id = child_id, nz_child = nz_child )
1155	ENDIF
1156	ENDDO
1157
1158	CALL pmc_s_setind_and_allocmem( child_id )
1159
1160	ENDDO ! m
1161
1162	IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN
1163	DEALLOCATE( child_x_left )
1164	DEALLOCATE( child_x_right )
1165	DEALLOCATE( child_y_south )
1166	DEALLOCATE( child_y_north )
1167	ENDIF
1168
1169
1170	CONTAINS
1171
1172
1173	SUBROUTINE pmci_create_index_list
1174
1175	IMPLICIT NONE
1176
1177	INTEGER(iwp) :: ilist !< Index-list index running over the child's parent-grid jc,ic-space
1178	INTEGER(iwp) :: index_list_size !< Dimension 2 of the array index_list
1179	INTEGER(iwp) :: ierr !< MPI error code
1180	INTEGER(iwp) :: ip !< Running parent-grid index on the child domain in the x-direction
1181	INTEGER(iwp) :: jp !< Running parent-grid index on the child domain in the y-direction
1182	INTEGER(iwp) :: n !< Running index over child subdomains
1183	INTEGER(iwp) :: nrx !< Parent subdomain dimension in the x-direction
1184	INTEGER(iwp) :: nry !< Parent subdomain dimension in the y-direction
1185	INTEGER(iwp) :: pex !< Two-dimensional subdomain (pe) index in the x-direction
1186	INTEGER(iwp) :: pey !< Two-dimensional subdomain (pe) index in the y-direction
1187	INTEGER(iwp) :: parent_pe !< Parent subdomain index (one-dimensional)
1188
1189	INTEGER(iwp), DIMENSION(2) :: pe_indices_2d !< Array for two-dimensional subdomain (pe)
1190	!< indices needed for MPI_CART_RANK
1191	INTEGER(iwp), DIMENSION(2) :: size_of_childs_parent_grid_bounds_all !< Dimensions of childs_parent_grid_bounds_all
1192	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: childs_parent_grid_bounds_all !< Array that contains the child's
1193	!< parent-grid index bounds for all its
1194	!< subdomains (pes)
1195	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: index_list !< Array that maps the child index space on
1196	!< the parent index- and subdomain spaces
1197
1198	IF ( myid == 0 ) THEN
1199
1200	CALL pmc_recv_from_child( child_id, size_of_childs_parent_grid_bounds_all, &
1201	2, 0, 40, ierr )
1202	ALLOCATE( childs_parent_grid_bounds_all(size_of_childs_parent_grid_bounds_all(1), &
1203	size_of_childs_parent_grid_bounds_all(2)) )
1204	CALL pmc_recv_from_child( child_id, childs_parent_grid_bounds_all, &
1205	SIZE( childs_parent_grid_bounds_all ), 0, 41, ierr )
1206	!
1207	!-- Compute size (dimension) of the index_list.
1208	index_list_size = 0
1209	DO n = 1, size_of_childs_parent_grid_bounds_all(2)
1210	index_list_size = index_list_size + &
1211	( childs_parent_grid_bounds_all(4,n) - childs_parent_grid_bounds_all(3,n) + 1 ) * &
1212	( childs_parent_grid_bounds_all(2,n) - childs_parent_grid_bounds_all(1,n) + 1 )
1213	ENDDO
1214
1215	ALLOCATE( index_list(6,index_list_size) )
1216
1217	nrx = nxr - nxl + 1
1218	nry = nyn - nys + 1
1219	ilist = 0
1220	!
1221	!-- Loop over all children PEs
1222	DO n = 1, size_of_childs_parent_grid_bounds_all(2) !
1223	!
1224	!-- Subspace along y required by actual child PE
1225	DO jp = childs_parent_grid_bounds_all(3,n), childs_parent_grid_bounds_all(4,n) ! jp = jps, jpn of child PE# n
1226	!
1227	!-- Subspace along x required by actual child PE
1228	DO ip = childs_parent_grid_bounds_all(1,n), childs_parent_grid_bounds_all(2,n) ! ip = ipl, ipr of child PE# n
1229
1230	pex = ip / nrx
1231	pey = jp / nry
1232	pe_indices_2d(1) = pex
1233	pe_indices_2d(2) = pey
1234	CALL MPI_CART_RANK( comm2d, pe_indices_2d, parent_pe, ierr )
1235
1236	ilist = ilist + 1
1237	!
1238	!-- First index in parent array ! TO_DO: IMPROVE THIS COMMENT
1239	index_list(1,ilist) = ip - ( pex * nrx ) + 1 + nbgp
1240	!
1241	!-- Second index in parent array ! TO_DO: IMPROVE THIS COMMENT
1242	index_list(2,ilist) = jp - ( pey * nry ) + 1 + nbgp
1243	!
1244	!-- x index of child's parent grid
1245	index_list(3,ilist) = ip - childs_parent_grid_bounds_all(1,n) + 1
1246	!
1247	!-- y index of child's parent grid
1248	index_list(4,ilist) = jp - childs_parent_grid_bounds_all(3,n) + 1
1249	!
1250	!-- PE number of child
1251	index_list(5,ilist) = n - 1
1252	!
1253	!-- PE number of parent
1254	index_list(6,ilist) = parent_pe
1255
1256	ENDDO
1257	ENDDO
1258	ENDDO
1259	!
1260	!-- TO_DO: Klaus: comment what is done here
1261	CALL pmc_s_set_2d_index_list( child_id, index_list(:,1:ilist) )
1262
1263	ELSE
1264	!
1265	!-- TO_DO: Klaus: comment why this dummy allocation is required
1266	ALLOCATE( index_list(6,1) )
1267	CALL pmc_s_set_2d_index_list( child_id, index_list )
1268	ENDIF
1269
1270	DEALLOCATE(index_list)
1271
1272	END SUBROUTINE pmci_create_index_list
1273
1274
1275
1276	SUBROUTINE set_child_edge_coords
1277	IMPLICIT NONE
1278
1279	INTEGER(iwp) :: nbgp_lpm = 1 !<
1280
1281
1282	nbgp_lpm = MIN( nbgp_lpm, nbgp )
1283
1284	childgrid(m)%nx = nx_child
1285	childgrid(m)%ny = ny_child
1286	childgrid(m)%nz = nz_child
1287	childgrid(m)%dx = dx_child
1288	childgrid(m)%dy = dy_child
1289	childgrid(m)%dz = dz_child
1290
1291	childgrid(m)%lx_coord = child_coord_x(0)
1292	childgrid(m)%lx_coord_b = child_coord_x(-nbgp_lpm)
1293	childgrid(m)%rx_coord = child_coord_x(nx_child) + dx_child
1294	childgrid(m)%rx_coord_b = child_coord_x(nx_child+nbgp_lpm) + dx_child
1295	childgrid(m)%sy_coord = child_coord_y(0)
1296	childgrid(m)%sy_coord_b = child_coord_y(-nbgp_lpm)
1297	childgrid(m)%ny_coord = child_coord_y(ny_child) + dy_child
1298	childgrid(m)%ny_coord_b = child_coord_y(ny_child+nbgp_lpm) + dy_child
1299	childgrid(m)%uz_coord = zmax_coarse(2)
1300	childgrid(m)%uz_coord_b = zmax_coarse(1)
1301
1302	END SUBROUTINE set_child_edge_coords
1303
1304	#endif
1305	END SUBROUTINE pmci_setup_parent
1306
1307
1308
1309	SUBROUTINE pmci_setup_child
1310
1311	#if defined( __parallel )
1312	IMPLICIT NONE
1313
1314	INTEGER(iwp) :: ierr !< MPI error code
1315	INTEGER(iwp) :: lb !< Running index for aerosol size bins
1316	INTEGER(iwp) :: lc !< Running index for aerosol mass bins
1317	INTEGER(iwp) :: lg !< Running index for salsa gases
1318	INTEGER(iwp) :: n !< Running index for number of chemical species
1319	INTEGER(iwp), DIMENSION(3) :: child_grid_dim !< Array for sending the child-grid dimensions to parent
1320
1321	REAL(wp), DIMENSION(5) :: child_grid_info !< Array for sending the child-grid spacings etc to parent
1322
1323	CHARACTER( LEN=da_namelen ) :: myname !<
1324	CHARACTER(LEN=5) :: salsa_char !<
1325
1326	!
1327	!-- Child setup
1328	!-- Root model does not have a parent and is not a child, therefore no child setup on root model
1329	IF ( .NOT. pmc_is_rootmodel() ) THEN
1330	!
1331	!-- ADD A DESCRIPTION HERE WHAT PMC_CHILDINIT DOES
1332	CALL pmc_childinit
1333	!
1334	!-- The arrays, which actually will be exchanged between child and parent
1335	!-- are defined Here AND ONLY HERE.
1336	!-- If a variable is removed, it only has to be removed from here.
1337	!-- Please check, if the arrays are in the list of POSSIBLE exchange arrays
1338	!-- in subroutines:
1339	!-- pmci_set_array_pointer (for parent arrays)
1340	!-- pmci_create_childs_parent_grid_arrays (for child's parent-grid arrays)
1341	CALL pmc_set_dataarray_name( 'coarse', 'u' ,'fine', 'u', ierr )
1342	CALL pmc_set_dataarray_name( 'coarse', 'v' ,'fine', 'v', ierr )
1343	CALL pmc_set_dataarray_name( 'coarse', 'w' ,'fine', 'w', ierr )
1344	!
1345	!-- Set data array name for TKE. Please note, nesting of TKE is actually
1346	!-- only done if both parent and child are in LES or in RANS mode. Due to
1347	!-- design of model coupler, however, data array names must be already
1348	!-- available at this point.
1349	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
1350	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
1351	.NOT. constant_diffusion ) ) THEN
1352	CALL pmc_set_dataarray_name( 'coarse', 'e' ,'fine', 'e', ierr )
1353	ENDIF
1354	!
1355	!-- Nesting of dissipation rate only if both parent and child are in RANS
1356	!-- mode and TKE-epsilo closure is applied. Please see also comment for TKE
1357	!-- above.
1358	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
1359	CALL pmc_set_dataarray_name( 'coarse', 'diss' ,'fine', 'diss', ierr )
1360	ENDIF
1361
1362	IF ( .NOT. neutral ) THEN
1363	CALL pmc_set_dataarray_name( 'coarse', 'pt' ,'fine', 'pt', ierr )
1364	ENDIF
1365
1366	IF ( humidity ) THEN
1367
1368	CALL pmc_set_dataarray_name( 'coarse', 'q' ,'fine', 'q', ierr )
1369
1370	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
1371	CALL pmc_set_dataarray_name( 'coarse', 'qc' ,'fine', 'qc', ierr )
1372	CALL pmc_set_dataarray_name( 'coarse', 'nc' ,'fine', 'nc', ierr )
1373	ENDIF
1374
1375	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
1376	CALL pmc_set_dataarray_name( 'coarse', 'qr' ,'fine', 'qr', ierr )
1377	CALL pmc_set_dataarray_name( 'coarse', 'nr' ,'fine', 'nr', ierr )
1378	ENDIF
1379
1380	ENDIF
1381
1382	IF ( passive_scalar ) THEN
1383	CALL pmc_set_dataarray_name( 'coarse', 's' ,'fine', 's', ierr )
1384	ENDIF
1385
1386	IF ( particle_advection ) THEN
1387	CALL pmc_set_dataarray_name( 'coarse', 'nr_part' ,'fine', &
1388	'nr_part', ierr )
1389	CALL pmc_set_dataarray_name( 'coarse', 'part_adr' ,'fine', &
1390	'part_adr', ierr )
1391	ENDIF
1392
1393	IF ( air_chemistry .AND. nest_chemistry ) THEN
1394	DO n = 1, nspec
1395	CALL pmc_set_dataarray_name( 'coarse', &
1396	'chem_' // &
1397	TRIM( chem_species(n)%name ), &
1398	'fine', &
1399	'chem_' // &
1400	TRIM( chem_species(n)%name ), &
1401	ierr )
1402	ENDDO
1403	ENDIF
1404
1405	IF ( salsa .AND. nest_salsa ) THEN
1406	DO lb = 1, nbins_aerosol
1407	WRITE(salsa_char,'(i0)') lb
1408	CALL pmc_set_dataarray_name( 'coarse', &
1409	'an_' // &
1410	TRIM( salsa_char ), &
1411	'fine', &
1412	'an_' // &
1413	TRIM( salsa_char ), &
1414	ierr )
1415	ENDDO
1416	DO lc = 1, nbins_aerosol * ncomponents_mass
1417	WRITE(salsa_char,'(i0)') lc
1418	CALL pmc_set_dataarray_name( 'coarse', &
1419	'am_' // &
1420	TRIM( salsa_char ), &
1421	'fine', &
1422	'am_' // &
1423	TRIM( salsa_char ), &
1424	ierr )
1425	ENDDO
1426	IF ( .NOT. salsa_gases_from_chem ) THEN
1427	DO lg = 1, ngases_salsa
1428	WRITE(salsa_char,'(i0)') lg
1429	CALL pmc_set_dataarray_name( 'coarse', &
1430	'sg_' // &
1431	TRIM( salsa_char ), &
1432	'fine', &
1433	'sg_' // &
1434	TRIM( salsa_char ), &
1435	ierr )
1436	ENDDO
1437	ENDIF
1438	ENDIF
1439
1440	CALL pmc_set_dataarray_name( lastentry = .TRUE. )
1441	!
1442	!-- Send grid to parent
1443	child_grid_dim(1) = nx
1444	child_grid_dim(2) = ny
1445	child_grid_dim(3) = nz
1446	child_grid_info(1) = zw(nzt+1)
1447	child_grid_info(2) = zw(nzt)
1448	child_grid_info(3) = dx
1449	child_grid_info(4) = dy
1450	child_grid_info(5) = dz(1)
1451
1452	IF ( myid == 0 ) THEN
1453
1454	CALL pmc_send_to_parent( child_grid_dim, SIZE( child_grid_dim ), 0, 123, ierr )
1455	CALL pmc_send_to_parent( child_grid_info, SIZE( child_grid_info ), 0, 124, ierr )
1456	CALL pmc_send_to_parent( coord_x, nx + 1 + 2 * nbgp, 0, 11, ierr )
1457	CALL pmc_send_to_parent( coord_y, ny + 1 + 2 * nbgp, 0, 12, ierr )
1458
1459	CALL pmc_recv_from_parent( rans_mode_parent, 1, 0, 19, ierr )
1460	!
1461	!-- Receive Coarse grid information.
1462	CALL pmc_recv_from_parent( parent_grid_info_real, &
1463	SIZE(parent_grid_info_real), 0, 21, ierr )
1464	CALL pmc_recv_from_parent( parent_grid_info_int, 3, 0, 22, ierr )
1465
1466	ENDIF
1467
1468	CALL MPI_BCAST( parent_grid_info_real, SIZE(parent_grid_info_real), &
1469	MPI_REAL, 0, comm2d, ierr )
1470	CALL MPI_BCAST( parent_grid_info_int, 3, MPI_INTEGER, 0, comm2d, ierr )
1471
1472	pg%dx = parent_grid_info_real(3)
1473	pg%dy = parent_grid_info_real(4)
1474	pg%dz = parent_grid_info_real(7)
1475	pg%nx = parent_grid_info_int(1)
1476	pg%ny = parent_grid_info_int(2)
1477	pg%nz = parent_grid_info_int(3)
1478	!
1479	!-- Get parent coordinates on coarse grid
1480	ALLOCATE( pg%coord_x(-nbgp:pg%nx+nbgp) )
1481	ALLOCATE( pg%coord_y(-nbgp:pg%ny+nbgp) )
1482	ALLOCATE( pg%dzu(1:pg%nz+1) )
1483	ALLOCATE( pg%dzw(1:pg%nz+1) )
1484	ALLOCATE( pg%zu(0:pg%nz+1) )
1485	ALLOCATE( pg%zw(0:pg%nz+1) )
1486	!
1487	!-- Get coarse grid coordinates and values of the z-direction from the parent
1488	IF ( myid == 0) THEN
1489	CALL pmc_recv_from_parent( pg%coord_x, pg%nx+1+2*nbgp, 0, 24, ierr )
1490	CALL pmc_recv_from_parent( pg%coord_y, pg%ny+1+2*nbgp, 0, 25, ierr )
1491	CALL pmc_recv_from_parent( pg%dzu, pg%nz+1, 0, 26, ierr )
1492	CALL pmc_recv_from_parent( pg%dzw, pg%nz+1, 0, 27, ierr )
1493	CALL pmc_recv_from_parent( pg%zu, pg%nz+2, 0, 28, ierr )
1494	CALL pmc_recv_from_parent( pg%zw, pg%nz+2, 0, 29, ierr )
1495	ENDIF
1496	!
1497	!-- Broadcast this information
1498	CALL MPI_BCAST( pg%coord_x, pg%nx+1+2*nbgp, MPI_REAL, 0, comm2d, ierr )
1499	CALL MPI_BCAST( pg%coord_y, pg%ny+1+2*nbgp, MPI_REAL, 0, comm2d, ierr )
1500	CALL MPI_BCAST( pg%dzu, pg%nz+1, MPI_REAL, 0, comm2d, ierr )
1501	CALL MPI_BCAST( pg%dzw, pg%nz+1, MPI_REAL, 0, comm2d, ierr )
1502	CALL MPI_BCAST( pg%zu, pg%nz+2, MPI_REAL, 0, comm2d, ierr )
1503	CALL MPI_BCAST( pg%zw, pg%nz+2, MPI_REAL, 0, comm2d, ierr )
1504	CALL MPI_BCAST( rans_mode_parent, 1, MPI_LOGICAL, 0, comm2d, ierr )
1505
1506	! CHECK IF pmci_check_grid_matching COULD BE MOVED HERE.
1507
1508	!
1509	!-- Find the index bounds for the nest domain in the coarse-grid index space
1510	CALL pmci_map_fine_to_coarse_grid
1511	!
1512	!-- TO_DO: Klaus give a comment what is happening here
1513	CALL pmc_c_get_2d_index_list
1514	!
1515	!-- Include couple arrays into child content
1516	!-- TO_DO: Klaus: better explain the above comment (what is child content?)
1517	CALL pmc_c_clear_next_array_list
1518
1519	n = 1
1520	lb = 1
1521	lc = 1
1522	lg = 1
1523
1524	DO WHILE ( pmc_c_getnextarray( myname ) )
1525	!
1526	!-- Note that cg%nz is not the original nz of parent, but the highest
1527	!-- parent-grid level needed for nesting.
1528	!-- Please note, in case of chemical species an additional parameter
1529	!-- need to be passed, which is required to set the pointer correctly
1530	!-- to the chemical-species data structure. Hence, first check if current
1531	!-- variable is a chemical species. If so, pass index id of respective
1532	!-- species and increment this subsequently.
1533	IF ( INDEX( TRIM( myname ), 'chem_' ) /= 0 ) THEN
1534	CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, n )
1535	n = n + 1
1536	ELSEIF ( INDEX( TRIM( myname ), 'an_' ) /= 0 ) THEN
1537	CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lb )
1538	lb = lb + 1
1539	ELSEIF ( INDEX( TRIM( myname ), 'am_' ) /= 0 ) THEN
1540	CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lc )
1541	lc = lc + 1
1542	ELSEIF ( INDEX( TRIM( myname ), 'sg_' ) /= 0 .AND. .NOT. salsa_gases_from_chem ) THEN
1543	CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lg )
1544	lg = lg + 1
1545	ELSE
1546	CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz )
1547	ENDIF
1548	ENDDO
1549	CALL pmc_c_setind_and_allocmem
1550	!
1551	!-- Precompute the index-mapping arrays
1552	CALL pmci_define_index_mapping
1553	!
1554	!-- Check that the child and parent grid lines do match
1555	CALL pmci_check_grid_matching
1556
1557	ENDIF
1558
1559	CONTAINS
1560
1561
1562	SUBROUTINE pmci_map_fine_to_coarse_grid
1563	!
1564	!-- Determine index bounds of interpolation/anterpolation area in the coarse
1565	!-- grid index space
1566	IMPLICIT NONE
1567
1568	INTEGER(iwp), DIMENSION(5,numprocs) :: coarse_bound_all !< Transfer array for parent-grid index bounds
1569
1570	INTEGER(iwp), DIMENSION(4) :: parent_bound_global !< Transfer array for global parent-grid index bounds
1571	INTEGER(iwp), DIMENSION(2) :: size_of_array !<
1572
1573	INTEGER(iwp) :: i !<
1574	INTEGER(iwp) :: iauxl !<
1575	INTEGER(iwp) :: iauxr !<
1576	INTEGER(iwp) :: ijaux !<
1577	INTEGER(iwp) :: j !<
1578	INTEGER(iwp) :: jauxs !<
1579	INTEGER(iwp) :: jauxn !<
1580
1581	REAL(wp) :: xexl !< Parent-grid array exceedance behind the left edge of the child PE subdomain
1582	REAL(wp) :: xexr !< Parent-grid array exceedance behind the right edge of the child PE subdomain
1583	REAL(wp) :: yexs !< Parent-grid array exceedance behind the south edge of the child PE subdomain
1584	REAL(wp) :: yexn !< Parent-grid array exceedance behind the north edge of the child PE subdomain
1585	REAL(wp) :: xcs !< RENAME
1586	REAL(wp) :: xce !< RENAME
1587	REAL(wp) :: ycs !< RENAME
1588	REAL(wp) :: yce !< RENAME
1589
1590	!
1591	!-- Determine the anterpolation index limits. If at least half of the
1592	!-- parent-grid cell is within the current child sub-domain, then it
1593	!-- is included in the current sub-domain's anterpolation domain.
1594	!-- Else the parent-grid cell is included in the neighbouring subdomain's
1595	!-- anterpolation domain, or not included at all if we are at the outer
1596	!-- edge of the child domain. This may occur especially when a large grid-spacing
1597	!-- ratio is used.
1598	!
1599	!-- Left
1600	!-- EXPLAIN THE EXTENSION HERE AND IN THE OTHER OCCASIONS (r, s, n)
1601	IF ( bc_dirichlet_l ) THEN
1602	xexl = 2 * pg%dx
1603	iauxl = 0
1604	ELSE
1605	xexl = 0.0_wp
1606	iauxl = 1
1607	ENDIF
1608	xcs = coord_x(nxl) - xexl
1609	DO i = 0, pg%nx
1610	IF ( pg%coord_x(i) + 0.5_wp * pg%dx >= xcs ) THEN ! Consider changing >= to ==
1611	ipl = MAX( 0, i )
1612	EXIT
1613	ENDIF
1614	ENDDO
1615	!
1616	!-- Right
1617	IF ( bc_dirichlet_r ) THEN
1618	xexr = 2 * pg%dx
1619	iauxr = 0
1620	ELSE
1621	xexr = 0.0_wp
1622	iauxr = 1
1623	ENDIF
1624	xce = coord_x(nxr+1) + xexr
1625	DO i = pg%nx, 0 , -1
1626	IF ( pg%coord_x(i) + 0.5_wp * pg%dx <= xce ) THEN
1627	ipr = MIN( pg%nx, MAX( ipl, i ) )
1628	EXIT
1629	ENDIF
1630	ENDDO
1631	!
1632	!-- South
1633	IF ( bc_dirichlet_s ) THEN
1634	yexs = 2 * pg%dy
1635	jauxs = 0
1636	ELSE
1637	yexs = 0.0_wp
1638	jauxs = 1
1639	ENDIF
1640	ycs = coord_y(nys) - yexs
1641	DO j = 0, pg%ny
1642	IF ( pg%coord_y(j) + 0.5_wp * pg%dy >= ycs ) THEN
1643	jps = MAX( 0, j )
1644	EXIT
1645	ENDIF
1646	ENDDO
1647	!
1648	!-- North
1649	IF ( bc_dirichlet_n ) THEN
1650	yexn = 2 * pg%dy
1651	jauxn = 0
1652	ELSE
1653	yexn = 0.0_wp
1654	jauxn = 1
1655	ENDIF
1656	yce = coord_y(nyn+1) + yexn
1657	DO j = pg%ny, 0 , -1
1658	IF ( pg%coord_y(j) + 0.5_wp * pg%dy <= yce ) THEN
1659	jpn = MIN( pg%ny, MAX( jps, j ) )
1660	EXIT
1661	ENDIF
1662	ENDDO
1663	!
1664	!-- Make sure that the indexing is contiguous (no gaps, no overlaps).
1665	!-- This is a safety measure mainly for cases with high grid-spacing
1666	!-- ratio and narrow child subdomains.
1667	IF ( nxl == 0 ) THEN
1668	CALL MPI_SEND( ipr, 1, MPI_INTEGER, pright, 717, comm2d, ierr )
1669	ELSE IF ( nxr == nx ) THEN
1670	CALL MPI_RECV( ijaux, 1, MPI_INTEGER, pleft, 717, comm2d, status, ierr )
1671	ipl = ijaux + 1
1672	ELSE
1673	CALL MPI_SEND( ipr, 1, MPI_INTEGER, pright, 717, comm2d, ierr )
1674	CALL MPI_RECV( ijaux, 1, MPI_INTEGER, pleft, 717, comm2d, status, ierr )
1675	ipl = ijaux + 1
1676	ENDIF
1677	IF ( nys == 0 ) THEN
1678	CALL MPI_SEND( jpn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr )
1679	ELSE IF ( nyn == ny ) THEN
1680	CALL MPI_RECV( ijaux, 1, MPI_INTEGER, psouth, 719, comm2d, status, ierr )
1681	jps = ijaux + 1
1682	ELSE
1683	CALL MPI_SEND( jpn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr )
1684	CALL MPI_RECV( ijaux, 1, MPI_INTEGER, psouth, 719, comm2d, status, ierr )
1685	jps = ijaux + 1
1686	ENDIF
1687
1688	WRITE(9,"('Pmci_map_fine_to_coarse_grid. Parent-grid array bounds: ',4(i4,2x))") &
1689	ipl, ipr, jps, jpn
1690	FLUSH(9)
1691
1692	coarse_bound(1) = ipl
1693	coarse_bound(2) = ipr
1694	coarse_bound(3) = jps
1695	coarse_bound(4) = jpn
1696	coarse_bound(5) = myid
1697	!
1698	!-- The following index bounds are used for allocating index mapping and some other auxiliary arrays
1699	ipla = ipl - iauxl
1700	ipra = ipr + iauxr
1701	jpsa = jps - jauxs
1702	jpna = jpn + jauxn
1703	!
1704	!-- Note that MPI_Gather receives data from all processes in the rank order
1705	!-- This fact is exploited in creating the index list in pmci_create_index_list
1706	! IMPROVE THIS COMMENT. EXPLAIN WHERE THIS INFORMATION IS NEEDED.
1707	CALL MPI_GATHER( coarse_bound, 5, MPI_INTEGER, coarse_bound_all, 5, &
1708	MPI_INTEGER, 0, comm2d, ierr )
1709
1710	IF ( myid == 0 ) THEN
1711	size_of_array(1) = SIZE( coarse_bound_all, 1 )
1712	size_of_array(2) = SIZE( coarse_bound_all, 2 )
1713	CALL pmc_send_to_parent( size_of_array, 2, 0, 40, ierr )
1714	CALL pmc_send_to_parent( coarse_bound_all, SIZE( coarse_bound_all ), 0, 41, ierr )
1715	!
1716	!-- Determine the global parent-grid index bounds
1717	parent_bound_global(1) = MINVAL( coarse_bound_all(1,:) )
1718	parent_bound_global(2) = MAXVAL( coarse_bound_all(2,:) )
1719	parent_bound_global(3) = MINVAL( coarse_bound_all(3,:) )
1720	parent_bound_global(4) = MAXVAL( coarse_bound_all(4,:) )
1721	ENDIF
1722	!
1723	!-- Broadcat the global parent-grid index bounds to all current child processes
1724	CALL MPI_BCAST( parent_bound_global, 4, MPI_INTEGER, 0, comm2d, ierr )
1725	iplg = parent_bound_global(1)
1726	iprg = parent_bound_global(2)
1727	jpsg = parent_bound_global(3)
1728	jpng = parent_bound_global(4)
1729	WRITE( 9, "('Pmci_map_fine_to_coarse_grid. Global parent-grid index bounds iplg, iprg, jpsg, jpng: ',4(i4,2x))" ) &
1730	iplg, iprg, jpsg, jpng
1731	FLUSH( 9 )
1732
1733	END SUBROUTINE pmci_map_fine_to_coarse_grid
1734
1735
1736
1737	SUBROUTINE pmci_define_index_mapping
1738	!
1739	!-- Precomputation of the mapping of the child- and parent-grid indices.
1740
1741	IMPLICIT NONE
1742
1743	INTEGER(iwp) :: i !< Child-grid index in the x-direction
1744	INTEGER(iwp) :: ii !< Parent-grid index in the x-direction
1745	INTEGER(iwp) :: istart !<
1746	INTEGER(iwp) :: ir !<
1747	INTEGER(iwp) :: iw !< Child-grid index limited to -1 <= iw <= nx+1 for wall_flags_0
1748	INTEGER(iwp) :: j !< Child-grid index in the y-direction
1749	INTEGER(iwp) :: jj !< Parent-grid index in the y-direction
1750	INTEGER(iwp) :: jstart !<
1751	INTEGER(iwp) :: jr !<
1752	INTEGER(iwp) :: jw !< Child-grid index limited to -1 <= jw <= ny+1 for wall_flags_0
1753	INTEGER(iwp) :: k !< Child-grid index in the z-direction
1754	INTEGER(iwp) :: kk !< Parent-grid index in the z-direction
1755	INTEGER(iwp) :: kstart !<
1756	INTEGER(iwp) :: kw !< Child-grid index limited to kw <= nzt+1 for wall_flags_0
1757
1758	!
1759	!-- Allocate child-grid work arrays for interpolation.
1760	igsr = NINT( pg%dx / dx, iwp )
1761	jgsr = NINT( pg%dy / dy, iwp )
1762	kgsr = NINT( pg%dzw(1) / dzw(1), iwp )
1763	WRITE(9,"('igsr, jgsr, kgsr: ',3(i3,2x))") igsr, jgsr, kgsr
1764	FLUSH(9)
1765	!
1766	!-- Determine index bounds for the parent-grid work arrays for
1767	!-- interpolation and allocate them.
1768	CALL pmci_allocate_workarrays
1769	!
1770	!-- Define the MPI-datatypes for parent-grid work array
1771	!-- exchange between the PE-subdomains.
1772	CALL pmci_create_workarray_exchange_datatypes
1773	!
1774	!-- First determine kcto and kctw which refer to the uppermost
1775	!-- coarse-grid levels below the child top-boundary level.
1776	kk = 0
1777	DO WHILE ( pg%zu(kk) <= zu(nzt) )
1778	kk = kk + 1
1779	ENDDO
1780	kcto = kk - 1
1781
1782	kk = 0
1783	DO WHILE ( pg%zw(kk) <= zw(nzt-1) )
1784	kk = kk + 1
1785	ENDDO
1786	kctw = kk - 1
1787
1788	WRITE( 9, "('kcto, kctw = ', 2(i3,2x))" ) kcto, kctw
1789	FLUSH( 9 )
1790
1791	ALLOCATE( iflu(ipla:ipra) )
1792	ALLOCATE( iflo(ipla:ipra) )
1793	ALLOCATE( ifuu(ipla:ipra) )
1794	ALLOCATE( ifuo(ipla:ipra) )
1795	ALLOCATE( jflv(jpsa:jpna) )
1796	ALLOCATE( jflo(jpsa:jpna) )
1797	ALLOCATE( jfuv(jpsa:jpna) )
1798	ALLOCATE( jfuo(jpsa:jpna) )
1799	ALLOCATE( kflw(0:pg%nz+1) )
1800	ALLOCATE( kflo(0:pg%nz+1) )
1801	ALLOCATE( kfuw(0:pg%nz+1) )
1802	ALLOCATE( kfuo(0:pg%nz+1) )
1803	ALLOCATE( ijkfc_u(0:pg%nz+1,jpsa:jpna,ipla:ipra) )
1804	ALLOCATE( ijkfc_v(0:pg%nz+1,jpsa:jpna,ipla:ipra) )
1805	ALLOCATE( ijkfc_w(0:pg%nz+1,jpsa:jpna,ipla:ipra) )
1806	ALLOCATE( ijkfc_s(0:pg%nz+1,jpsa:jpna,ipla:ipra) )
1807
1808	ijkfc_u = 0
1809	ijkfc_v = 0
1810	ijkfc_w = 0
1811	ijkfc_s = 0
1812	!
1813	!-- i-indices of u for each ii-index value
1814	istart = nxlg
1815	DO ii = ipla, ipra
1816	!
1817	!-- The parent and child grid lines do always match in x, hence we
1818	!-- use only the local k,j-child-grid plane for the anterpolation.
1819	!-- However, icru still has to be stored separately as these index bounds
1820	!-- are passed as arguments to the interpolation and anterpolation
1821	!-- subroutines.
1822	i = istart
1823	DO WHILE ( coord_x(i) < pg%coord_x(ii) .AND. i < nxrg )
1824	i = i + 1
1825	ENDDO
1826	iflu(ii) = MIN( MAX( i, nxlg ), nxrg )
1827	ifuu(ii) = iflu(ii)
1828	istart = iflu(ii)
1829	!
1830	!-- Print out the index bounds for checking and debugging purposes
1831	WRITE( 9, "('pmci_define_index_mapping, ii, iflu, ifuu: ', 3(i4,2x))" ) &
1832	ii, iflu(ii), ifuu(ii)
1833	FLUSH( 9 )
1834	ENDDO
1835	WRITE( 9, * )
1836	!
1837	!-- i-indices of others for each ii-index value
1838	istart = nxlg
1839	DO ii = ipla, ipra
1840	i = istart
1841	DO WHILE ( ( coord_x(i) + 0.5_wp * dx < pg%coord_x(ii) ) .AND. &
1842	( i < nxrg ) )
1843	i = i + 1
1844	ENDDO
1845	iflo(ii) = MIN( MAX( i, nxlg ), nxrg )
1846	ir = i
1847	DO WHILE ( ( coord_x(ir) + 0.5_wp * dx <= pg%coord_x(ii) + pg%dx ) &
1848	.AND. ( i < nxrg+1 ) )
1849	i = i + 1
1850	ir = MIN( i, nxrg )
1851	ENDDO
1852	ifuo(ii) = MIN( MAX( i-1, iflo(ii) ), nxrg )
1853	istart = iflo(ii)
1854	!
1855	!-- Print out the index bounds for checking and debugging purposes
1856	WRITE( 9, "('pmci_define_index_mapping, ii, iflo, ifuo: ', 3(i4,2x))" ) &
1857	ii, iflo(ii), ifuo(ii)
1858	FLUSH( 9 )
1859	ENDDO
1860	WRITE( 9, * )
1861	!
1862	!-- j-indices of v for each jj-index value
1863	jstart = nysg
1864	DO jj = jpsa, jpna
1865	!
1866	!-- The parent and child grid lines do always match in y, hence we
1867	!-- use only the local k,i-child-grid plane for the anterpolation.
1868	!-- However, jcnv still has to be stored separately as these index bounds
1869	!-- are passed as arguments to the interpolation and anterpolation
1870	!-- subroutines.
1871	j = jstart
1872	DO WHILE ( coord_y(j) < pg%coord_y(jj) .AND. j < nyng )
1873	j = j + 1
1874	ENDDO
1875	jflv(jj) = MIN( MAX( j, nysg ), nyng )
1876	jfuv(jj) = jflv(jj)
1877	jstart = jflv(jj)
1878	!
1879	!-- Print out the index bounds for checking and debugging purposes
1880	WRITE( 9, "('pmci_define_index_mapping, jj, jflv, jfuv: ', 3(i4,2x))" ) &
1881	jj, jflv(jj), jfuv(jj)
1882	FLUSH(9)
1883	ENDDO
1884	WRITE( 9, * )
1885	!
1886	!-- j-indices of others for each jj-index value
1887	jstart = nysg
1888	DO jj = jpsa, jpna
1889	j = jstart
1890	DO WHILE ( ( coord_y(j) + 0.5_wp * dy < pg%coord_y(jj) ) .AND. ( j < nyng ) )
1891	j = j + 1
1892	ENDDO
1893	jflo(jj) = MIN( MAX( j, nysg ), nyng )
1894	jr = j
1895	DO WHILE ( ( coord_y(jr) + 0.5_wp * dy <= pg%coord_y(jj) + pg%dy ) .AND. ( j < nyng+1 ) )
1896	j = j + 1
1897	jr = MIN( j, nyng )
1898	ENDDO
1899	jfuo(jj) = MIN( MAX( j-1, jflo(jj) ), nyng )
1900	jstart = jflo(jj)
1901	!
1902	!-- Print out the index bounds for checking and debugging purposes
1903	WRITE( 9, "('pmci_define_index_mapping, jj, jflo, jfuo: ', 3(i4,2x))" ) &
1904	jj, jflo(jj), jfuo(jj)
1905	FLUSH( 9 )
1906	ENDDO
1907	WRITE( 9, * )
1908	!
1909	!-- k-indices of w for each kk-index value
1910	!-- Note that anterpolation index limits are needed also for the top boundary
1911	!-- ghost cell level because they are used also in the interpolation.
1912	kstart = 0
1913	kflw(0) = 0
1914	kfuw(0) = 0
1915	DO kk = 1, pg%nz+1
1916	!
1917	!-- The parent and child grid lines do always match in z, hence we
1918	!-- use only the local j,i-child-grid plane for the anterpolation.
1919	!-- However, kctw still has to be stored separately as these index bounds
1920	!-- are passed as arguments to the interpolation and anterpolation
1921	!-- subroutines.
1922	k = kstart
1923	DO WHILE ( ( zw(k) < pg%zw(kk) ) .AND. ( k < nzt+1 ) )
1924	k = k + 1
1925	ENDDO
1926	kflw(kk) = MIN( MAX( k, 1 ), nzt + 1 )
1927	kfuw(kk) = kflw(kk)
1928	kstart = kflw(kk)
1929	!
1930	!-- Print out the index bounds for checking and debugging purposes
1931	WRITE( 9, "('pmci_define_index_mapping, kk, kflw, kfuw: ', 4(i4,2x), 2(e12.5,2x))" ) &
1932	kk, kflw(kk), kfuw(kk), nzt, pg%zu(kk), pg%zw(kk)
1933	FLUSH( 9 )
1934	ENDDO
1935	WRITE( 9, * )
1936	!
1937	!-- k-indices of others for each kk-index value
1938	kstart = 0
1939	kflo(0) = 0
1940	kfuo(0) = 0
1941	!
1942	!-- Note that anterpolation index limits are needed also for the top boundary
1943	!-- ghost cell level because they are used also in the interpolation.
1944	DO kk = 1, pg%nz+1
1945	k = kstart
1946	DO WHILE ( ( zu(k) < pg%zw(kk-1) ) .AND. ( k <= nzt ) )
1947	k = k + 1
1948	ENDDO
1949	kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 )
1950	DO WHILE ( ( zu(k) <= pg%zw(kk) ) .AND. ( k <= nzt+1 ) )
1951	k = k + 1
1952	IF ( k > nzt + 1 ) EXIT ! This EXIT is to prevent zu(k) from flowing over.
1953	ENDDO
1954	kfuo(kk) = MIN( MAX( k-1, kflo(kk) ), nzt + 1 )
1955	kstart = kflo(kk)
1956	ENDDO
1957	!
1958	!-- Set the k-index bounds separately for the parent-grid cells pg%nz and pg%nz+1
1959	!-- although they are not actually needed.
1960	! WHY IS THIS LIKE THIS? REVISE WITH CARE.
1961	kflo(pg%nz) = nzt+1
1962	kfuo(pg%nz) = nzt+kgsr
1963	kflo(pg%nz+1) = nzt+kgsr
1964	kfuo(pg%nz+1) = nzt+kgsr
1965	!
1966	!-- Print out the index bounds for checking and debugging purposes
1967	DO kk = 1, pg%nz+1
1968	WRITE( 9, "('pmci_define_index_mapping, kk, kflo, kfuo: ', 4(i4,2x), 2(e12.5,2x))" ) &
1969	kk, kflo(kk), kfuo(kk), nzt, pg%zu(kk), pg%zw(kk)
1970	FLUSH( 9 )
1971	ENDDO
1972	WRITE( 9, * )
1973	!
1974	!-- Precomputation of number of fine-grid nodes inside parent-grid cells.
1975	!-- Note that ii, jj, and kk are parent-grid indices.
1976	!-- This information is needed in the anterpolation.
1977	!-- The indices for wall_flags_0 (kw,jw,iw) must be limited to the range
1978	!-- [-1,...,nx/ny/nzt+1] in order to avoid zero values on the outer ghost nodes.
1979	DO ii = ipla, ipra
1980	DO jj = jpsa, jpna
1981	DO kk = 0, pg%nz+1
1982	!
1983	!-- u-component
1984	DO i = iflu(ii), ifuu(ii)
1985	iw = MAX( MIN( i, nx+1 ), -1 )
1986	DO j = jflo(jj), jfuo(jj)
1987	jw = MAX( MIN( j, ny+1 ), -1 )
1988	DO k = kflo(kk), kfuo(kk)
1989	kw = MIN( k, nzt+1 )
1990	ijkfc_u(kk,jj,ii) = ijkfc_u(kk,jj,ii) &
1991	+ MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 1 ) )
1992	ENDDO
1993	ENDDO
1994	ENDDO
1995	!
1996	!-- v-component
1997	DO i = iflo(ii), ifuo(ii)
1998	iw = MAX( MIN( i, nx+1 ), -1 )
1999	DO j = jflv(jj), jfuv(jj)
2000	jw = MAX( MIN( j, ny+1 ), -1 )
2001	DO k = kflo(kk), kfuo(kk)
2002	kw = MIN( k, nzt+1 )
2003	ijkfc_v(kk,jj,ii) = ijkfc_v(kk,jj,ii) &
2004	+ MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 2 ) )
2005	ENDDO
2006	ENDDO
2007	ENDDO
2008	!
2009	!-- scalars
2010	DO i = iflo(ii), ifuo(ii)
2011	iw = MAX( MIN( i, nx+1 ), -1 )
2012	DO j = jflo(jj), jfuo(jj)
2013	jw = MAX( MIN( j, ny+1 ), -1 )
2014	DO k = kflo(kk), kfuo(kk)
2015	kw = MIN( k, nzt+1 )
2016	ijkfc_s(kk,jj,ii) = ijkfc_s(kk,jj,ii) &
2017	+ MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 0 ) )
2018	ENDDO
2019	ENDDO
2020	ENDDO
2021	!
2022	!-- w-component
2023	DO i = iflo(ii), ifuo(ii)
2024	iw = MAX( MIN( i, nx+1 ), -1 )
2025	DO j = jflo(jj), jfuo(jj)
2026	jw = MAX( MIN( j, ny+1 ), -1 )
2027	DO k = kflw(kk), kfuw(kk)
2028	kw = MIN( k, nzt+1 )
2029	ijkfc_w(kk,jj,ii) = ijkfc_w(kk,jj,ii) &
2030	+ MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 3 ) )
2031	ENDDO
2032	ENDDO
2033	ENDDO
2034
2035	ENDDO ! kk
2036	ENDDO ! jj
2037	ENDDO ! ii
2038
2039	END SUBROUTINE pmci_define_index_mapping
2040
2041
2042
2043	SUBROUTINE pmci_allocate_workarrays
2044	!
2045	!-- Allocate parent-grid work-arrays for interpolation
2046	IMPLICIT NONE
2047
2048	!
2049	!-- Determine and store the PE-subdomain dependent index bounds
2050	IF ( bc_dirichlet_l ) THEN
2051	iplw = ipl + 1
2052	ELSE
2053	iplw = ipl - 1
2054	ENDIF
2055
2056	IF ( bc_dirichlet_r ) THEN
2057	iprw = ipr - 1
2058	ELSE
2059	iprw = ipr + 1
2060	ENDIF
2061
2062	IF ( bc_dirichlet_s ) THEN
2063	jpsw = jps + 1
2064	ELSE
2065	jpsw = jps - 1
2066	ENDIF
2067
2068	IF ( bc_dirichlet_n ) THEN
2069	jpnw = jpn - 1
2070	ELSE
2071	jpnw = jpn + 1
2072	ENDIF
2073	!
2074	!-- Left and right boundaries.
2075	ALLOCATE( workarr_lr(0:pg%nz+1,jpsw:jpnw,0:2) )
2076	!
2077	!-- South and north boundaries.
2078	ALLOCATE( workarr_sn(0:pg%nz+1,0:2,iplw:iprw) )
2079	!
2080	!-- Top boundary.
2081	ALLOCATE( workarr_t(0:2,jpsw:jpnw,iplw:iprw) )
2082
2083	END SUBROUTINE pmci_allocate_workarrays
2084
2085
2086
2087	SUBROUTINE pmci_create_workarray_exchange_datatypes
2088	!
2089	!-- Define specific MPI types for workarr-exchange.
2090	IMPLICIT NONE
2091
2092	!
2093	!-- For the left and right boundaries
2094	CALL MPI_TYPE_VECTOR( 3, pg%nz+2, (jpnw-jpsw+1)*(pg%nz+2), MPI_REAL, &
2095	workarr_lr_exchange_type, ierr )
2096	CALL MPI_TYPE_COMMIT( workarr_lr_exchange_type, ierr )
2097	!
2098	!-- For the south and north boundaries
2099	CALL MPI_TYPE_VECTOR( 1, 3(pg%nz+2), 3(pg%nz+2), MPI_REAL, &
2100	workarr_sn_exchange_type, ierr )
2101	CALL MPI_TYPE_COMMIT( workarr_sn_exchange_type, ierr )
2102	!
2103	!-- For the top-boundary x-slices
2104	CALL MPI_TYPE_VECTOR( iprw-iplw+1, 3, 3*(jpnw-jpsw+1), MPI_REAL, &
2105	workarr_t_exchange_type_x, ierr )
2106	CALL MPI_TYPE_COMMIT( workarr_t_exchange_type_x, ierr )
2107	!
2108	!-- For the top-boundary y-slices
2109	CALL MPI_TYPE_VECTOR( 1, 3(jpnw-jpsw+1), 3(jpnw-jpsw+1), MPI_REAL, &
2110	workarr_t_exchange_type_y, ierr )
2111	CALL MPI_TYPE_COMMIT( workarr_t_exchange_type_y, ierr )
2112
2113	END SUBROUTINE pmci_create_workarray_exchange_datatypes
2114
2115
2116
2117	SUBROUTINE pmci_check_grid_matching
2118	!
2119	!-- Check that the grid lines of child and parent do match.
2120	!-- Also check that the child subdomain width is not smaller than
2121	!-- the parent grid spacing in the respective direction.
2122	IMPLICIT NONE
2123
2124	INTEGER(iwp) :: non_matching_height = 0 !< Flag for non-matching child-domain height
2125	INTEGER(iwp) :: non_matching_lower_left_corner = 0 !< Flag for non-matching lower left corner
2126	INTEGER(iwp) :: non_matching_upper_right_corner = 0 !< Flag for non-matching upper right corner
2127	INTEGER(iwp) :: non_int_gsr_x = 0 !< Flag for non-integer grid-spacing ration in x-direction
2128	INTEGER(iwp) :: non_int_gsr_y = 0 !< Flag for non-integer grid-spacing ration in y-direction
2129	INTEGER(iwp) :: non_int_gsr_z = 0 !< Flag for non-integer grid-spacing ration in z-direction
2130	INTEGER(iwp) :: too_narrow_pesd_x = 0 !< Flag for too narrow pe-subdomain in x-direction
2131	INTEGER(iwp) :: too_narrow_pesd_y = 0 !< Flag for too narrow pe-subdomain in y-direction
2132
2133	REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp !< Relative tolerence for grid-line matching
2134
2135	REAL(wp) :: child_ngp_x_l !< Number of gridpoints in child subdomain in x-direction
2136	!< converted to REAL(wp)
2137	REAL(wp) :: child_ngp_y_l !< Number of gridpoints in child subdomain in y-direction
2138	!< converted to REAL(wp)
2139	REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction
2140	REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction
2141	REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction
2142	REAL(wp) :: upper_right_coord_x !< X-coordinate of the upper right corner of the child domain
2143	REAL(wp) :: upper_right_coord_y !< Y-coordinate of the upper right corner of the child domain
2144
2145
2146	IF ( myid == 0 ) THEN
2147
2148	tolex = tolefac * dx
2149	toley = tolefac * dy
2150	tolez = tolefac * dz(1)
2151	!
2152	!-- First check that the child domain lower left corner matches the parent grid lines.
2153	IF ( MOD( lower_left_coord_x, pg%dx ) > tolex ) non_matching_lower_left_corner = 1
2154	IF ( MOD( lower_left_coord_y, pg%dy ) > toley ) non_matching_lower_left_corner = 1
2155	!
2156	!-- Then check that the child doman upper right corner matches the parent grid lines.
2157	upper_right_coord_x = lower_left_coord_x + ( nx + 1 ) * dx
2158	upper_right_coord_y = lower_left_coord_y + ( ny + 1 ) * dy
2159	IF ( MOD( upper_right_coord_x, pg%dx ) > tolex ) non_matching_upper_right_corner = 1
2160	IF ( MOD( upper_right_coord_y, pg%dy ) > toley ) non_matching_upper_right_corner = 1
2161	!
2162	!-- Also check that the cild domain height matches the parent grid lines.
2163	IF ( MOD( zw(nzt), pg%dz ) > tolez ) non_matching_height = 1
2164	!
2165	!-- Check that the grid-spacing ratios in each direction are integer valued.
2166	IF ( MOD( pg%dx, dx ) > tolex ) non_int_gsr_x = 1
2167	IF ( MOD( pg%dy, dy ) > toley ) non_int_gsr_y = 1
2168	!
2169	!-- In the z-direction, all levels need to be checked separately against grid stretching
2170	!-- which is not allowed.
2171	DO n = 0, kctw+1
2172	IF ( ABS( pg%zw(n) - zw(kflw(n)) ) > tolez ) non_int_gsr_z = 1
2173	ENDDO
2174
2175	child_ngp_x_l = REAL( nxr - nxl + 1, KIND=wp )
2176	IF ( child_ngp_x_l / REAL( igsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_x = 1
2177	child_ngp_y_l = REAL( nyn - nys + 1, KIND=wp )
2178	IF ( child_ngp_y_l / REAL( jgsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_y = 1
2179
2180	IF ( non_matching_height > 0 ) THEN
2181	WRITE( message_string, * ) 'nested child domain height must match ', &
2182	'its parent grid lines'
2183	CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 )
2184	ENDIF
2185
2186	IF ( non_matching_lower_left_corner > 0 ) THEN
2187	WRITE( message_string, * ) 'nested child domain lower left ', &
2188	'corner must match its parent grid lines'
2189	CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 )
2190	ENDIF
2191
2192	IF ( non_matching_upper_right_corner > 0 ) THEN
2193	WRITE( message_string, * ) 'nested child domain upper right ', &
2194	'corner must match its parent grid lines'
2195	CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 )
2196	ENDIF
2197
2198	IF ( non_int_gsr_x > 0 ) THEN
2199	WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dx / child dx ) ', &
2200	'must have an integer value'
2201	CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 )
2202	ENDIF
2203
2204	IF ( non_int_gsr_y > 0 ) THEN
2205	WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dy / child dy ) ', &
2206	'must have an integer value'
2207	CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 )
2208	ENDIF
2209
2210	IF ( non_int_gsr_z > 0 ) THEN
2211	WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dz / child dz ) ', &
2212	'must have an integer value for each z-level'
2213	CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 )
2214	ENDIF
2215
2216	IF ( too_narrow_pesd_x > 0 ) THEN
2217	WRITE( message_string, * ) 'child subdomain width in x-direction must not be ', &
2218	'smaller than its parent grid dx. Change the PE-grid ', &
2219	'setting (npex, npey) to satisfy this requirement.'
2220	CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 )
2221	ENDIF
2222
2223	IF ( too_narrow_pesd_y > 0 ) THEN
2224	WRITE( message_string, * ) 'child subdomain width in y-direction must not be ', &
2225	'smaller than its parent grid dy. Change the PE-grid ', &
2226	'setting (npex, npey) to satisfy this requirement.'
2227	CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 )
2228	ENDIF
2229
2230	ENDIF ! ( myid == 0 )
2231
2232	END SUBROUTINE pmci_check_grid_matching
2233
2234	#endif
2235	END SUBROUTINE pmci_setup_child
2236
2237
2238
2239	SUBROUTINE pmci_setup_coordinates
2240
2241	#if defined( __parallel )
2242	IMPLICIT NONE
2243
2244	INTEGER(iwp) :: i !<
2245	INTEGER(iwp) :: j !<
2246
2247	!
2248	!-- Create coordinate arrays.
2249	ALLOCATE( coord_x(-nbgp:nx+nbgp) )
2250	ALLOCATE( coord_y(-nbgp:ny+nbgp) )
2251
2252	DO i = -nbgp, nx + nbgp
2253	coord_x(i) = lower_left_coord_x + i * dx
2254	ENDDO
2255
2256	DO j = -nbgp, ny + nbgp
2257	coord_y(j) = lower_left_coord_y + j * dy
2258	ENDDO
2259
2260	#endif
2261	END SUBROUTINE pmci_setup_coordinates
2262
2263	!------------------------------------------------------------------------------!
2264	! Description:
2265	! ------------
2266	!> In this subroutine the number of coupled arrays is determined.
2267	!------------------------------------------------------------------------------!
2268	SUBROUTINE pmci_num_arrays
2269
2270	#if defined( __parallel )
2271	USE pmc_general, &
2272	ONLY: pmc_max_array
2273
2274	IMPLICIT NONE
2275	!
2276	!-- The number of coupled arrays depends on the model settings. At least
2277	!-- 5 arrays need to be coupled (u, v, w, e, diss). Please note, actually
2278	!-- e and diss (TKE and dissipation rate) are only required if RANS-RANS
2279	!-- nesting is applied, but memory is allocated nevertheless. This is because
2280	!-- the information whether they are needed or not is retrieved at a later
2281	!-- point in time. In case e and diss are not needed, they are also not
2282	!-- exchanged between parent and child.
2283	pmc_max_array = 5
2284	!
2285	!-- pt
2286	IF ( .NOT. neutral ) pmc_max_array = pmc_max_array + 1
2287
2288	IF ( humidity ) THEN
2289	!
2290	!-- q
2291	pmc_max_array = pmc_max_array + 1
2292	!
2293	!-- qc, nc
2294	IF ( bulk_cloud_model .AND. microphysics_morrison ) &
2295	pmc_max_array = pmc_max_array + 2
2296	!
2297	!-- qr, nr
2298	IF ( bulk_cloud_model .AND. microphysics_seifert ) &
2299	pmc_max_array = pmc_max_array + 2
2300	ENDIF
2301	!
2302	!-- s
2303	IF ( passive_scalar ) pmc_max_array = pmc_max_array + 1
2304	!
2305	!-- nr_part, part_adr
2306	IF ( particle_advection ) pmc_max_array = pmc_max_array + 2
2307	!
2308	!-- Chemistry, depends on number of species
2309	IF ( air_chemistry .AND. nest_chemistry ) pmc_max_array = pmc_max_array + nspec
2310	!
2311	!-- SALSA, depens on the number aerosol size bins and chemical components +
2312	!-- the number of default gases
2313	IF ( salsa .AND. nest_salsa ) pmc_max_array = pmc_max_array + nbins_aerosol + &
2314	nbins_aerosol * ncomponents_mass
2315	IF ( .NOT. salsa_gases_from_chem ) pmc_max_array = pmc_max_array + ngases_salsa
2316
2317
2318	#endif
2319
2320	END SUBROUTINE pmci_num_arrays
2321
2322
2323	SUBROUTINE pmci_set_array_pointer( name, child_id, nz_child, n )
2324
2325	IMPLICIT NONE
2326
2327	INTEGER(iwp), INTENT(IN) :: child_id !<
2328	INTEGER(iwp), INTENT(IN) :: nz_child !<
2329
2330	INTEGER(iwp), INTENT(IN), OPTIONAL :: n !< index of chemical species
2331
2332	CHARACTER(LEN=*), INTENT(IN) :: name !<
2333	!
2334	!-- Local variables:
2335	INTEGER(iwp) :: ierr !< MPI error code
2336
2337	INTEGER(idp), POINTER, DIMENSION(:,:) :: i_2d !<
2338
2339	REAL(wp), POINTER, DIMENSION(:,:) :: p_2d !<
2340	REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d !<
2341	REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d_sec !<
2342
2343
2344	NULLIFY( p_3d )
2345	NULLIFY( p_2d )
2346	NULLIFY( i_2d )
2347	!
2348	!-- List of array names, which can be coupled.
2349	!-- In case of 3D please change also the second array for the pointer version
2350	IF ( TRIM(name) == "u" ) p_3d => u
2351	IF ( TRIM(name) == "v" ) p_3d => v
2352	IF ( TRIM(name) == "w" ) p_3d => w
2353	IF ( TRIM(name) == "e" ) p_3d => e
2354	IF ( TRIM(name) == "pt" ) p_3d => pt
2355	IF ( TRIM(name) == "q" ) p_3d => q
2356	IF ( TRIM(name) == "qc" ) p_3d => qc
2357	IF ( TRIM(name) == "qr" ) p_3d => qr
2358	IF ( TRIM(name) == "nr" ) p_3d => nr
2359	IF ( TRIM(name) == "nc" ) p_3d => nc
2360	IF ( TRIM(name) == "s" ) p_3d => s
2361	IF ( TRIM(name) == "diss" ) p_3d => diss
2362	IF ( TRIM(name) == "nr_part" ) i_2d => nr_part
2363	IF ( TRIM(name) == "part_adr" ) i_2d => part_adr
2364	IF ( INDEX( TRIM(name), "chem_" ) /= 0 ) p_3d => chem_species(n)%conc
2365	IF ( INDEX( TRIM(name), "an_" ) /= 0 ) p_3d => aerosol_number(n)%conc
2366	IF ( INDEX( TRIM(name), "am_" ) /= 0 ) p_3d => aerosol_mass(n)%conc
2367	IF ( INDEX( TRIM(name), "sg_" ) /= 0 .AND. .NOT. salsa_gases_from_chem ) &
2368	p_3d => salsa_gas(n)%conc
2369	!
2370	!-- Next line is just an example for a 2D array (not active for coupling!)
2371	!-- Please note, that z0 has to be declared as TARGET array in modules.f90.
2372	! IF ( TRIM(name) == "z0" ) p_2d => z0
2373	IF ( TRIM(name) == "u" ) p_3d_sec => u_2
2374	IF ( TRIM(name) == "v" ) p_3d_sec => v_2
2375	IF ( TRIM(name) == "w" ) p_3d_sec => w_2
2376	IF ( TRIM(name) == "e" ) p_3d_sec => e_2
2377	IF ( TRIM(name) == "pt" ) p_3d_sec => pt_2
2378	IF ( TRIM(name) == "q" ) p_3d_sec => q_2
2379	IF ( TRIM(name) == "qc" ) p_3d_sec => qc_2
2380	IF ( TRIM(name) == "qr" ) p_3d_sec => qr_2
2381	IF ( TRIM(name) == "nr" ) p_3d_sec => nr_2
2382	IF ( TRIM(name) == "nc" ) p_3d_sec => nc_2
2383	IF ( TRIM(name) == "s" ) p_3d_sec => s_2
2384	IF ( TRIM(name) == "diss" ) p_3d_sec => diss_2
2385	IF ( INDEX( TRIM(name), "chem_" ) /= 0 ) p_3d_sec => spec_conc_2(:,:,:,n)
2386	IF ( INDEX( TRIM(name), "an_" ) /= 0 ) p_3d_sec => nconc_2(:,:,:,n)
2387	IF ( INDEX( TRIM(name), "am_" ) /= 0 ) p_3d_sec => mconc_2(:,:,:,n)
2388	IF ( INDEX( TRIM(name), "sg_" ) /= 0 .AND. .NOT. salsa_gases_from_chem ) &
2389	p_3d_sec => gconc_2(:,:,:,n)
2390
2391	IF ( ASSOCIATED( p_3d ) ) THEN
2392	CALL pmc_s_set_dataarray( child_id, p_3d, nz_child, nz, array_2 = p_3d_sec )
2393	ELSEIF ( ASSOCIATED( p_2d ) ) THEN
2394	CALL pmc_s_set_dataarray( child_id, p_2d )
2395	ELSEIF ( ASSOCIATED( i_2d ) ) THEN
2396	CALL pmc_s_set_dataarray( child_id, i_2d )
2397	ELSE
2398	!
2399	!-- Give only one message for the root domain
2400	IF ( myid == 0 .AND. cpl_id == 1 ) THEN
2401	message_string = 'pointer for array "' // TRIM( name ) // '" can''t be associated'
2402	CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 )
2403	ELSE
2404	!
2405	!-- Avoid others to continue
2406	CALL MPI_BARRIER( comm2d, ierr )
2407	ENDIF
2408
2409	ENDIF
2410
2411	END SUBROUTINE pmci_set_array_pointer
2412
2413
2414
2415	INTEGER FUNCTION get_number_of_children()
2416
2417	IMPLICIT NONE
2418
2419
2420	#if defined( __parallel )
2421	get_number_of_children = SIZE( pmc_parent_for_child ) - 1
2422	#else
2423	get_number_of_children = 0
2424	#endif
2425
2426	RETURN
2427
2428	END FUNCTION get_number_of_children
2429
2430
2431
2432	INTEGER FUNCTION get_childid( id_index )
2433
2434	IMPLICIT NONE
2435
2436	INTEGER, INTENT(IN) :: id_index !<
2437
2438
2439	#if defined( __parallel )
2440	get_childid = pmc_parent_for_child(id_index)
2441	#else
2442	get_childid = 0
2443	#endif
2444
2445	RETURN
2446
2447	END FUNCTION get_childid
2448
2449
2450
2451	SUBROUTINE get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, sy_coord, sy_coord_b, &
2452	ny_coord, ny_coord_b, uz_coord, uz_coord_b )
2453
2454	IMPLICIT NONE
2455
2456	INTEGER,INTENT(IN) :: m !<
2457
2458	REAL(wp),INTENT(OUT) :: lx_coord, lx_coord_b !<
2459	REAL(wp),INTENT(OUT) :: rx_coord, rx_coord_b !<
2460	REAL(wp),INTENT(OUT) :: sy_coord, sy_coord_b !<
2461	REAL(wp),INTENT(OUT) :: ny_coord, ny_coord_b !<
2462	REAL(wp),INTENT(OUT) :: uz_coord, uz_coord_b !<
2463
2464
2465	lx_coord = childgrid(m)%lx_coord
2466	rx_coord = childgrid(m)%rx_coord
2467	sy_coord = childgrid(m)%sy_coord
2468	ny_coord = childgrid(m)%ny_coord
2469	uz_coord = childgrid(m)%uz_coord
2470
2471	lx_coord_b = childgrid(m)%lx_coord_b
2472	rx_coord_b = childgrid(m)%rx_coord_b
2473	sy_coord_b = childgrid(m)%sy_coord_b
2474	ny_coord_b = childgrid(m)%ny_coord_b
2475	uz_coord_b = childgrid(m)%uz_coord_b
2476
2477	END SUBROUTINE get_child_edges
2478
2479
2480
2481	SUBROUTINE get_child_gridspacing( m, dx, dy,dz )
2482
2483	IMPLICIT NONE
2484
2485	INTEGER, INTENT(IN) :: m !<
2486
2487	REAL(wp), INTENT(OUT) :: dx,dy !<
2488
2489	REAL(wp), INTENT(OUT), OPTIONAL :: dz !<
2490
2491
2492	dx = childgrid(m)%dx
2493	dy = childgrid(m)%dy
2494	IF ( PRESENT( dz ) ) THEN
2495	dz = childgrid(m)%dz
2496	ENDIF
2497
2498	END SUBROUTINE get_child_gridspacing
2499
2500
2501
2502	SUBROUTINE pmci_create_childs_parent_grid_arrays( name, is, ie, js, je, nzc, n )
2503
2504	IMPLICIT NONE
2505
2506	INTEGER(iwp), INTENT(IN) :: ie !< RENAME ie, is, je, js?
2507	INTEGER(iwp), INTENT(IN) :: is !<
2508	INTEGER(iwp), INTENT(IN) :: je !<
2509	INTEGER(iwp), INTENT(IN) :: js !<
2510	INTEGER(iwp), INTENT(IN) :: nzc !< nzc is pg%nz, but note that pg%nz is not the original nz of parent,
2511	!< but the highest parent-grid level needed for nesting.
2512	INTEGER(iwp), INTENT(IN), OPTIONAL :: n !< number of chemical species / salsa variables
2513
2514	CHARACTER(LEN=*), INTENT(IN) :: name !<
2515	!
2516	!-- Local variables:
2517	INTEGER(iwp) :: ierr !<
2518
2519	INTEGER(idp), POINTER,DIMENSION(:,:) :: i_2d !<
2520
2521	REAL(wp), POINTER,DIMENSION(:,:) :: p_2d !<
2522	REAL(wp), POINTER,DIMENSION(:,:,:) :: p_3d !<
2523
2524	NULLIFY( p_3d )
2525	NULLIFY( p_2d )
2526	NULLIFY( i_2d )
2527	!
2528	!-- List of array names, which can be coupled
2529	IF ( TRIM( name ) == "u" ) THEN
2530	IF ( .NOT. ALLOCATED( uc ) ) ALLOCATE( uc(0:nzc+1,js:je,is:ie) )
2531	p_3d => uc
2532	ELSEIF ( TRIM( name ) == "v" ) THEN
2533	IF ( .NOT. ALLOCATED( vc ) ) ALLOCATE( vc(0:nzc+1,js:je,is:ie) )
2534	p_3d => vc
2535	ELSEIF ( TRIM( name ) == "w" ) THEN
2536	IF ( .NOT. ALLOCATED( wc ) ) ALLOCATE( wc(0:nzc+1,js:je,is:ie) )
2537	p_3d => wc
2538	ELSEIF ( TRIM( name ) == "e" ) THEN
2539	IF ( .NOT. ALLOCATED( ec ) ) ALLOCATE( ec(0:nzc+1,js:je,is:ie) )
2540	p_3d => ec
2541	ELSEIF ( TRIM( name ) == "diss" ) THEN
2542	IF ( .NOT. ALLOCATED( dissc ) ) ALLOCATE( dissc(0:nzc+1,js:je,is:ie) )
2543	p_3d => dissc
2544	ELSEIF ( TRIM( name ) == "pt") THEN
2545	IF ( .NOT. ALLOCATED( ptc ) ) ALLOCATE( ptc(0:nzc+1,js:je,is:ie) )
2546	p_3d => ptc
2547	ELSEIF ( TRIM( name ) == "q") THEN
2548	IF ( .NOT. ALLOCATED( q_c ) ) ALLOCATE( q_c(0:nzc+1,js:je,is:ie) )
2549	p_3d => q_c
2550	ELSEIF ( TRIM( name ) == "qc") THEN
2551	IF ( .NOT. ALLOCATED( qcc ) ) ALLOCATE( qcc(0:nzc+1,js:je,is:ie) )
2552	p_3d => qcc
2553	ELSEIF ( TRIM( name ) == "qr") THEN
2554	IF ( .NOT. ALLOCATED( qrc ) ) ALLOCATE( qrc(0:nzc+1,js:je,is:ie) )
2555	p_3d => qrc
2556	ELSEIF ( TRIM( name ) == "nr") THEN
2557	IF ( .NOT. ALLOCATED( nrc ) ) ALLOCATE( nrc(0:nzc+1,js:je,is:ie) )
2558	p_3d => nrc
2559	ELSEIF ( TRIM( name ) == "nc") THEN
2560	IF ( .NOT. ALLOCATED( ncc ) ) ALLOCATE( ncc(0:nzc+1,js:je,is:ie) )
2561	p_3d => ncc
2562	ELSEIF ( TRIM( name ) == "s") THEN
2563	IF ( .NOT. ALLOCATED( sc ) ) ALLOCATE( sc(0:nzc+1,js:je,is:ie) )
2564	p_3d => sc
2565	ELSEIF ( TRIM( name ) == "nr_part") THEN
2566	IF ( .NOT. ALLOCATED( nr_partc ) ) ALLOCATE( nr_partc(js:je,is:ie) )
2567	i_2d => nr_partc
2568	ELSEIF ( TRIM( name ) == "part_adr") THEN
2569	IF ( .NOT. ALLOCATED( part_adrc ) ) ALLOCATE( part_adrc(js:je,is:ie) )
2570	i_2d => part_adrc
2571	ELSEIF ( TRIM( name(1:5) ) == "chem_" ) THEN
2572	IF ( .NOT. ALLOCATED( chem_spec_c ) ) ALLOCATE( chem_spec_c(0:nzc+1,js:je,is:ie,1:nspec) )
2573	p_3d => chem_spec_c(:,:,:,n)
2574	ELSEIF ( TRIM( name(1:3) ) == "an_" ) THEN
2575	IF ( .NOT. ALLOCATED( aerosol_number_c ) ) &
2576	ALLOCATE( aerosol_number_c(0:nzc+1,js:je,is:ie,1:nbins_aerosol) )
2577	p_3d => aerosol_number_c(:,:,:,n)
2578	ELSEIF ( TRIM( name(1:3) ) == "am_" ) THEN
2579	IF ( .NOT. ALLOCATED( aerosol_mass_c ) ) &
2580	ALLOCATE( aerosol_mass_c(0:nzc+1,js:je,is:ie,1:(nbins_aerosol*ncomponents_mass) ) )
2581	p_3d => aerosol_mass_c(:,:,:,n)
2582	ELSEIF ( TRIM( name(1:3) ) == "sg_" .AND. .NOT. salsa_gases_from_chem ) &
2583	THEN
2584	IF ( .NOT. ALLOCATED( salsa_gas_c ) ) &
2585	ALLOCATE( salsa_gas_c(0:nzc+1,js:je,is:ie,1:ngases_salsa) )
2586	p_3d => salsa_gas_c(:,:,:,n)
2587	!ELSEIF (trim(name) == "z0") then
2588	!IF (.not.allocated(z0c)) allocate(z0c(js:je, is:ie))
2589	!p_2d => z0c
2590	ENDIF
2591
2592	IF ( ASSOCIATED( p_3d ) ) THEN
2593	CALL pmc_c_set_dataarray( p_3d )
2594	ELSEIF ( ASSOCIATED( p_2d ) ) THEN
2595	CALL pmc_c_set_dataarray( p_2d )
2596	ELSEIF ( ASSOCIATED( i_2d ) ) THEN
2597	CALL pmc_c_set_dataarray( i_2d )
2598	ELSE
2599	!
2600	!-- Give only one message for the first child domain.
2601	IF ( myid == 0 .AND. cpl_id == 2 ) THEN
2602	message_string = 'pointer for array "' // TRIM( name ) // &
2603	'" can''t be associated'
2604	CALL message( 'pmci_create_childs_parent_grid_arrays', 'PA0170', 3, 2, 0, 6, 0 )
2605	ELSE
2606	!
2607	!-- Prevent others from continuing in case the abort is to come.
2608	CALL MPI_BARRIER( comm2d, ierr )
2609	ENDIF
2610
2611	ENDIF
2612
2613	END SUBROUTINE pmci_create_childs_parent_grid_arrays
2614
2615
2616	!
2617	! E N D O F S E T U P R O U T I N E S
2618	!
2619	SUBROUTINE pmci_parent_initialize
2620
2621	!
2622	!-- Send data for the children in order to let them create initial
2623	!-- conditions by interpolating the parent-domain fields.
2624	#if defined( __parallel )
2625	IMPLICIT NONE
2626
2627	INTEGER(iwp) :: child_id !<
2628	INTEGER(iwp) :: m !<
2629	REAL(wp) :: waittime !<
2630
2631
2632	DO m = 1, SIZE( pmc_parent_for_child ) - 1
2633	child_id = pmc_parent_for_child(m)
2634	CALL pmc_s_fillbuffer( child_id, waittime=waittime )
2635	ENDDO
2636
2637	#endif
2638	END SUBROUTINE pmci_parent_initialize
2639
2640
2641
2642	SUBROUTINE pmci_child_initialize
2643
2644	!
2645	!-- Create initial conditions for the current child domain by interpolating
2646	!-- the parent-domain fields.
2647	#if defined( __parallel )
2648	IMPLICIT NONE
2649
2650	INTEGER(iwp) :: ic !< Child-grid index in x-direction
2651	INTEGER(iwp) :: jc !< Child-grid index in y-direction
2652	INTEGER(iwp) :: kc !< Child-grid index in z-direction
2653	INTEGER(iwp) :: lb !< Running index for aerosol size bins
2654	INTEGER(iwp) :: lc !< Running index for aerosol mass bins
2655	INTEGER(iwp) :: lg !< Running index for salsa gases
2656	INTEGER(iwp) :: n !< Running index for chemical species
2657	REAL(wp) :: waittime !< Waiting time
2658
2659	!
2660	!-- Root model is never anyone's child
2661	IF ( cpl_id > 1 ) THEN
2662	!
2663	!-- Get data from the parent
2664	CALL pmc_c_getbuffer( waittime = waittime )
2665	!
2666	!-- The interpolation.
2667	CALL pmci_interp_1sto_all ( u, uc, kcto, iflu, ifuu, jflo, jfuo, kflo, kfuo, 'u' )
2668	CALL pmci_interp_1sto_all ( v, vc, kcto, iflo, ifuo, jflv, jfuv, kflo, kfuo, 'v' )
2669	CALL pmci_interp_1sto_all ( w, wc, kctw, iflo, ifuo, jflo, jfuo, kflw, kfuw, 'w' )
2670
2671	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
2672	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
2673	.NOT. constant_diffusion ) ) THEN
2674	CALL pmci_interp_1sto_all ( e, ec, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 'e' )
2675	ENDIF
2676
2677	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
2678	CALL pmci_interp_1sto_all ( diss, dissc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2679	ENDIF
2680
2681	IF ( .NOT. neutral ) THEN
2682	CALL pmci_interp_1sto_all ( pt, ptc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2683	ENDIF
2684
2685	IF ( humidity ) THEN
2686
2687	CALL pmci_interp_1sto_all ( q, q_c, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2688
2689	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
2690	CALL pmci_interp_1sto_all ( qc, qcc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2691	CALL pmci_interp_1sto_all ( nc, ncc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2692	ENDIF
2693
2694	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
2695	CALL pmci_interp_1sto_all ( qr, qrc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2696	CALL pmci_interp_1sto_all ( nr, nrc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2697	ENDIF
2698
2699	ENDIF
2700
2701	IF ( passive_scalar ) THEN
2702	CALL pmci_interp_1sto_all ( s, sc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2703	ENDIF
2704
2705	IF ( air_chemistry .AND. nest_chemistry ) THEN
2706	DO n = 1, nspec
2707	CALL pmci_interp_1sto_all ( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
2708	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2709	ENDDO
2710	ENDIF
2711
2712	IF ( salsa .AND. nest_salsa ) THEN
2713	DO lb = 1, nbins_aerosol
2714	CALL pmci_interp_1sto_all ( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
2715	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2716	ENDDO
2717	DO lc = 1, nbins_aerosol * ncomponents_mass
2718	CALL pmci_interp_1sto_all ( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
2719	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2720	ENDDO
2721	IF ( .NOT. salsa_gases_from_chem ) THEN
2722	DO lg = 1, ngases_salsa
2723	CALL pmci_interp_1sto_all ( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
2724	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' )
2725	ENDDO
2726	ENDIF
2727	ENDIF
2728
2729	IF ( topography /= 'flat' ) THEN
2730	!
2731	!-- Inside buildings set velocities back to zero.
2732	DO ic = nxlg, nxrg
2733	DO jc = nysg, nyng
2734	DO kc = nzb, nzt
2735	u(kc,jc,ic) = MERGE( u(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 1 ) )
2736	v(kc,jc,ic) = MERGE( v(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 2 ) )
2737	w(kc,jc,ic) = MERGE( w(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 3 ) )
2738	u_p(kc,jc,ic) = MERGE( u_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 1 ) )
2739	v_p(kc,jc,ic) = MERGE( v_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 2 ) )
2740	w_p(kc,jc,ic) = MERGE( w_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 3 ) )
2741	ENDDO
2742	ENDDO
2743	ENDDO
2744	ENDIF
2745	ENDIF
2746
2747
2748	CONTAINS
2749
2750
2751	SUBROUTINE pmci_interp_1sto_all( child_array, parent_array, kct, ifl, ifu, jfl, jfu, kfl, kfu, &
2752	var )
2753	!
2754	!-- Interpolation of the internal values for the child-domain initialization
2755	IMPLICIT NONE
2756
2757	INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node
2758
2759	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain
2760	!< parent cell - x direction
2761	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain
2762	!< parent cell - x direction
2763	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain
2764	!< parent cell - y direction
2765	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain
2766	!< parent cell - y direction
2767	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain
2768	!< parent cell - z direction
2769	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain
2770	!< parent cell - z direction
2771	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array
2772	REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array
2773
2774	CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's'
2775	!
2776	!-- Local variables:
2777	INTEGER(iwp) :: ic !< Running child-grid index in the x-direction
2778	INTEGER(iwp) :: icfirst !< Leftmost child-grid index initialized by the main loops (usually icfirst == icl_init)
2779	INTEGER(iwp) :: iclast !< Rightmost child-grid index initialized by the main loops (usually iclast == icr_init)
2780	INTEGER(iwp) :: icl_init !< Left child-grid index bound for initialization in the x-direction
2781	INTEGER(iwp) :: icr_init !< Right child-grid index bound for initialization in the x-direction
2782	INTEGER(iwp) :: jc !< Running child-grid index in the y-direction
2783	INTEGER(iwp) :: jcfirst !< Southmost child-grid index initialized by the main loops (usually jcfirst == jcs_init)
2784	INTEGER(iwp) :: jclast !< Northmost child-grid index initialized by the main loops (usually jclast == jcn_init)
2785	INTEGER(iwp) :: jcs_init !< South child-grid index bound for initialization in the y-direction
2786	INTEGER(iwp) :: jcn_init !< North child-grid index bound for initialization in the y-direction
2787	INTEGER(iwp) :: kc !< Running child-grid index in the z-direction
2788	INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction
2789	INTEGER(iwp) :: ipl_init !< Left parent-grid index bound for initialization in the x-direction
2790	INTEGER(iwp) :: ipr_init !< Right parent-grid index bound for initialization in the x-direction
2791	INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction
2792	INTEGER(iwp) :: jps_init !< South parent-grid index bound for initialization in the y-direction
2793	INTEGER(iwp) :: jpn_init !< North parent-grid index bound for initialization in the y-direction
2794	INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction
2795
2796
2797	ipl_init = ipl
2798	ipr_init = ipr
2799	jps_init = jps
2800	jpn_init = jpn
2801	icl_init = nxl
2802	icr_init = nxr
2803	jcs_init = nys
2804	jcn_init = nyn
2805
2806	IF ( nesting_mode /= 'vertical' ) THEN
2807	IF ( bc_dirichlet_l ) THEN
2808	ipl_init = ipl + 1
2809	icl_init = nxl - 1
2810	!
2811	!-- For u, nxl is a ghost node, but not for the other variables
2812	IF ( var == 'u' ) THEN
2813	ipl_init = ipl + 2
2814	icl_init = nxl
2815	ENDIF
2816	ENDIF
2817	IF ( bc_dirichlet_s ) THEN
2818	jps_init = jps + 1
2819	jcs_init = nys - 1
2820	!
2821	!-- For v, nys is a ghost node, but not for the other variables
2822	IF ( var == 'v' ) THEN
2823	jps_init = jps + 2
2824	jcs_init = nys
2825	ENDIF
2826	ENDIF
2827	IF ( bc_dirichlet_r ) THEN
2828	ipr_init = ipr - 1
2829	icr_init = nxr + 1
2830	ENDIF
2831	IF ( bc_dirichlet_n ) THEN
2832	jpn_init = jpn - 1
2833	jcn_init = nyn + 1
2834	ENDIF
2835	ENDIF
2836
2837	child_array(:,:,:) = 0.0_wp
2838
2839	IF ( var == 'u' ) THEN
2840
2841	icfirst = ifl(ipl_init)
2842	iclast = ifl(ipr_init+1) - 1
2843	jcfirst = jfl(jps_init)
2844	jclast = jfu(jpn_init)
2845	DO ip = ipl_init, ipr_init
2846	DO jp = jps_init, jpn_init
2847	DO kp = 0, kct + 1
2848
2849	DO ic = ifl(ip), ifl(ip+1)-1
2850	DO jc = jfl(jp), jfu(jp)
2851	DO kc = kfl(kp), MIN( kfu(kp), nzt+1 )
2852	child_array(kc,jc,ic) = parent_array(kp,jp,ip)
2853	ENDDO
2854	ENDDO
2855	ENDDO
2856
2857	ENDDO
2858	ENDDO
2859	ENDDO
2860
2861	ELSE IF ( var == 'v' ) THEN
2862
2863	icfirst = ifl(ipl_init)
2864	iclast = ifu(ipr_init)
2865	jcfirst = jfl(jps_init)
2866	jclast = jfl(jpn_init+1) - 1
2867	DO ip = ipl_init, ipr_init
2868	DO jp = jps_init, jpn_init
2869	DO kp = 0, kct + 1
2870
2871	DO ic = ifl(ip), ifu(ip)
2872	DO jc = jfl(jp), jfl(jp+1)-1
2873	DO kc = kfl(kp), MIN( kfu(kp), nzt+1 )
2874	child_array(kc,jc,ic) = parent_array(kp,jp,ip)
2875	ENDDO
2876	ENDDO
2877	ENDDO
2878
2879	ENDDO
2880	ENDDO
2881	ENDDO
2882
2883	ELSE IF ( var == 'w' ) THEN
2884
2885	icfirst = ifl(ipl_init)
2886	iclast = ifu(ipr_init)
2887	jcfirst = jfl(jps_init)
2888	jclast = jfu(jpn_init)
2889	DO ip = ipl_init, ipr_init
2890	DO jp = jps_init, jpn_init
2891	DO kp = 1, kct + 1
2892
2893	DO ic = ifl(ip), ifu(ip)
2894	DO jc = jfl(jp), jfu(jp)
2895	!
2896	!-- Because the kp-loop for w starts from kp=1 instead of 0
2897	child_array(nzb,jc,ic) = 0.0_wp
2898	DO kc = kfu(kp-1)+1, kfu(kp)
2899	child_array(kc,jc,ic) = parent_array(kp,jp,ip)
2900	ENDDO
2901	ENDDO
2902	ENDDO
2903
2904	ENDDO
2905	ENDDO
2906	ENDDO
2907
2908	ELSE ! scalars
2909
2910	icfirst = ifl(ipl_init)
2911	iclast = ifu(ipr_init)
2912	jcfirst = jfl(jps_init)
2913	jclast = jfu(jpn_init)
2914	DO ip = ipl_init, ipr_init
2915	DO jp = jps_init, jpn_init
2916	DO kp = 0, kct + 1
2917
2918	DO ic = ifl(ip), ifu(ip)
2919	DO jc = jfl(jp), jfu(jp)
2920	DO kc = kfl(kp), MIN( kfu(kp), nzt+1 )
2921	child_array(kc,jc,ic) = parent_array(kp,jp,ip)
2922	ENDDO
2923	ENDDO
2924	ENDDO
2925
2926	ENDDO
2927	ENDDO
2928	ENDDO
2929
2930	ENDIF ! var
2931	!
2932	!-- If the number of grid points in child subdomain in x- or y-direction
2933	!-- (nxr - nxl + 1 and/or nyn - nys + 1) is not integer divisible by the grid spacing
2934	!-- ratio in its direction (igsr and/or jgsr), the above loops will return with
2935	!-- unfilled gaps in the initial fields. These gaps, if present, are filled here.
2936	IF ( icfirst > icl_init ) THEN
2937	DO ic = icl_init, icfirst - 1
2938	child_array(:,:,ic) = child_array(:,:,icfirst)
2939	ENDDO
2940	ENDIF
2941	IF ( iclast < icr_init ) THEN
2942	DO ic = iclast + 1, icr_init
2943	child_array(:,:,ic) = child_array(:,:,iclast)
2944	ENDDO
2945	ENDIF
2946	IF ( jcfirst > jcs_init ) THEN
2947	DO jc = jcs_init, jcfirst - 1
2948	child_array(:,jc,:) = child_array(:,jcfirst,:)
2949	ENDDO
2950	ENDIF
2951	IF ( jclast < jcn_init ) THEN
2952	DO jc = jclast + 1, jcn_init
2953	child_array(:,jc,:) = child_array(:,jclast,:)
2954	ENDDO
2955	ENDIF
2956
2957	END SUBROUTINE pmci_interp_1sto_all
2958
2959	#endif
2960	END SUBROUTINE pmci_child_initialize
2961
2962
2963
2964	SUBROUTINE pmci_check_setting_mismatches
2965	!
2966	!-- Check for mismatches between settings of root and child variables
2967	!-- (e.g., all children have to follow the end_time settings of the root model).
2968	!-- The root model overwrites variables in the other models, so these variables
2969	!-- only need to be set once in file PARIN.
2970
2971	#if defined( __parallel )
2972
2973	USE control_parameters, &
2974	ONLY: dt_restart, end_time, message_string, restart_time, time_restart
2975
2976	IMPLICIT NONE
2977
2978	INTEGER :: ierr !< MPI error code
2979
2980	REAL(wp) :: dt_restart_root !<
2981	REAL(wp) :: end_time_root !<
2982	REAL(wp) :: restart_time_root !<
2983	REAL(wp) :: time_restart_root !<
2984
2985	!
2986	!-- Check the time to be simulated.
2987	!-- Here, and in the following, the root process communicates the respective
2988	!-- variable to all others, and its value will then be compared with the local
2989	!-- values.
2990	IF ( pmc_is_rootmodel() ) end_time_root = end_time
2991	CALL MPI_BCAST( end_time_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
2992
2993	IF ( .NOT. pmc_is_rootmodel() ) THEN
2994	IF ( end_time /= end_time_root ) THEN
2995	WRITE( message_string, * ) 'mismatch between root model and ', &
2996	'child settings:& end_time(root) = ', end_time_root, &
2997	'& end_time(child) = ', end_time, '& child value is set', &
2998	' to root value'
2999	CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
3000	0 )
3001	end_time = end_time_root
3002	ENDIF
3003	ENDIF
3004	!
3005	!-- Same for restart time
3006	IF ( pmc_is_rootmodel() ) restart_time_root = restart_time
3007	CALL MPI_BCAST( restart_time_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
3008
3009	IF ( .NOT. pmc_is_rootmodel() ) THEN
3010	IF ( restart_time /= restart_time_root ) THEN
3011	WRITE( message_string, * ) 'mismatch between root model and ', &
3012	'child settings: & restart_time(root) = ', restart_time_root, &
3013	'& restart_time(child) = ', restart_time, '& child ', &
3014	'value is set to root value'
3015	CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
3016	0 )
3017	restart_time = restart_time_root
3018	ENDIF
3019	ENDIF
3020	!
3021	!-- Same for dt_restart
3022	IF ( pmc_is_rootmodel() ) dt_restart_root = dt_restart
3023	CALL MPI_BCAST( dt_restart_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
3024
3025	IF ( .NOT. pmc_is_rootmodel() ) THEN
3026	IF ( dt_restart /= dt_restart_root ) THEN
3027	WRITE( message_string, * ) 'mismatch between root model and ', &
3028	'child settings: & dt_restart(root) = ', dt_restart_root, &
3029	'& dt_restart(child) = ', dt_restart, '& child ', &
3030	'value is set to root value'
3031	CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
3032	0 )
3033	dt_restart = dt_restart_root
3034	ENDIF
3035	ENDIF
3036	!
3037	!-- Same for time_restart
3038	IF ( pmc_is_rootmodel() ) time_restart_root = time_restart
3039	CALL MPI_BCAST( time_restart_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
3040
3041	IF ( .NOT. pmc_is_rootmodel() ) THEN
3042	IF ( time_restart /= time_restart_root ) THEN
3043	WRITE( message_string, * ) 'mismatch between root model and ', &
3044	'child settings: & time_restart(root) = ', time_restart_root, &
3045	'& time_restart(child) = ', time_restart, '& child ', &
3046	'value is set to root value'
3047	CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
3048	0 )
3049	time_restart = time_restart_root
3050	ENDIF
3051	ENDIF
3052
3053	#endif
3054
3055	END SUBROUTINE pmci_check_setting_mismatches
3056
3057
3058
3059	SUBROUTINE pmci_synchronize
3060
3061	#if defined( __parallel )
3062	!
3063	!-- Unify the time steps for each model and synchronize using
3064	!-- MPI_ALLREDUCE with the MPI_MIN operator over all processes using
3065	!-- the global communicator MPI_COMM_WORLD.
3066
3067	IMPLICIT NONE
3068
3069	INTEGER(iwp) :: ierr !< MPI error code
3070	REAL(wp) :: dtl !< Local time step of the current process
3071	REAL(wp) :: dtg !< Global time step defined as the global minimum of dtl of all processes
3072
3073	IF ( debug_output ) THEN
3074	WRITE( debug_string, * ) 'pmci_synchronize'
3075	CALL debug_message( debug_string, 'start' )
3076	ENDIF
3077
3078	dtl = dt_3d
3079	CALL MPI_ALLREDUCE( dtl, dtg, 1, MPI_REAL, MPI_MIN, MPI_COMM_WORLD, ierr )
3080	dt_3d = dtg
3081
3082	IF ( debug_output ) THEN
3083	WRITE( debug_string, * ) 'pmci_synchronize'
3084	CALL debug_message( debug_string, 'end' )
3085	ENDIF
3086
3087	#endif
3088	END SUBROUTINE pmci_synchronize
3089
3090
3091
3092	SUBROUTINE pmci_set_swaplevel( swaplevel )
3093
3094	!
3095	!-- After each Runge-Kutta sub-timestep, alternately set buffer one or buffer
3096	!-- two active
3097
3098	IMPLICIT NONE
3099
3100	INTEGER(iwp), INTENT(IN) :: swaplevel !< swaplevel (1 or 2) of PALM's timestep
3101
3102	INTEGER(iwp) :: child_id !< Child id of the child number m
3103	INTEGER(iwp) :: m !< Loop index over all children of the current parent
3104
3105	#if defined( __parallel )
3106	DO m = 1, SIZE( pmc_parent_for_child ) - 1
3107	child_id = pmc_parent_for_child(m)
3108	CALL pmc_s_set_active_data_array( child_id, swaplevel )
3109	ENDDO
3110	#endif
3111	END SUBROUTINE pmci_set_swaplevel
3112
3113
3114
3115	SUBROUTINE pmci_datatrans( local_nesting_mode )
3116	!
3117	!-- This subroutine controls the nesting according to the nestpar
3118	!-- parameter nesting_mode (two-way (default) or one-way) and the
3119	!-- order of anterpolations according to the nestpar parameter
3120	!-- nesting_datatransfer_mode (cascade, overlap or mixed (default)).
3121	!-- Although nesting_mode is a variable of this model, pass it as
3122	!-- an argument to allow for example to force one-way initialization
3123	!-- phase.
3124	!-- Note that interpolation ( parent_to_child ) must always be carried
3125	!-- out before anterpolation ( child_to_parent ).
3126
3127	IMPLICIT NONE
3128
3129	CHARACTER(LEN=*), INTENT(IN) :: local_nesting_mode !< Nesting mode: 'one-way', 'two-way' or 'vertical'
3130
3131	!
3132	!-- Debug location message
3133	IF ( debug_output ) THEN
3134	WRITE( debug_string, * ) 'pmci_datatrans'
3135	CALL debug_message( debug_string, 'start' )
3136	ENDIF
3137
3138	IF ( TRIM( local_nesting_mode ) == 'one-way' ) THEN
3139
3140	CALL pmci_child_datatrans( parent_to_child )
3141	CALL pmci_parent_datatrans( parent_to_child )
3142
3143	ELSE
3144
3145	IF ( nesting_datatransfer_mode == 'cascade' ) THEN
3146
3147	CALL pmci_child_datatrans( parent_to_child )
3148	CALL pmci_parent_datatrans( parent_to_child )
3149
3150	CALL pmci_parent_datatrans( child_to_parent )
3151	CALL pmci_child_datatrans( child_to_parent )
3152
3153	ELSEIF ( nesting_datatransfer_mode == 'overlap') THEN
3154
3155	CALL pmci_parent_datatrans( parent_to_child )
3156	CALL pmci_child_datatrans( parent_to_child )
3157
3158	CALL pmci_child_datatrans( child_to_parent )
3159	CALL pmci_parent_datatrans( child_to_parent )
3160
3161	ELSEIF ( TRIM( nesting_datatransfer_mode ) == 'mixed' ) THEN
3162
3163	CALL pmci_parent_datatrans( parent_to_child )
3164	CALL pmci_child_datatrans( parent_to_child )
3165
3166	CALL pmci_parent_datatrans( child_to_parent )
3167	CALL pmci_child_datatrans( child_to_parent )
3168
3169	ENDIF
3170
3171	ENDIF
3172	!
3173	!-- Debug location message
3174	IF ( debug_output ) THEN
3175	WRITE( debug_string, * ) 'pmci_datatrans'
3176	CALL debug_message( debug_string, 'end' )
3177	ENDIF
3178
3179	END SUBROUTINE pmci_datatrans
3180
3181
3182
3183	SUBROUTINE pmci_parent_datatrans( direction )
3184
3185	IMPLICIT NONE
3186
3187	INTEGER(iwp), INTENT(IN) :: direction !< Direction of the data transfer: 'parent_to_child' or 'child_to_parent'
3188
3189	#if defined( __parallel )
3190	INTEGER(iwp) :: child_id !< Child id of the child number m
3191	INTEGER(iwp) :: i !< Parent-grid index in x-direction
3192	INTEGER(iwp) :: j !< Parent-grid index in y-direction
3193	INTEGER(iwp) :: k !< Parent-grid index in z-direction
3194	INTEGER(iwp) :: m !< Loop index over all children of the current parent
3195
3196
3197	DO m = 1, SIZE( pmc_parent_for_child ) - 1
3198	child_id = pmc_parent_for_child(m)
3199	IF ( direction == parent_to_child ) THEN
3200	CALL cpu_log( log_point_s(71), 'pmc parent send', 'start' )
3201	CALL pmc_s_fillbuffer( child_id )
3202	CALL cpu_log( log_point_s(71), 'pmc parent send', 'stop' )
3203	ELSE
3204	!
3205	!-- Communication from child to parent
3206	CALL cpu_log( log_point_s(72), 'pmc parent recv', 'start' )
3207	child_id = pmc_parent_for_child(m)
3208	CALL pmc_s_getdata_from_buffer( child_id )
3209	CALL cpu_log( log_point_s(72), 'pmc parent recv', 'stop' )
3210	!
3211	!-- The anterpolated data is now available in u etc
3212	IF ( topography /= 'flat' ) THEN
3213	!
3214	!-- Inside buildings/topography reset velocities back to zero.
3215	!-- Scalars (pt, q, s, km, kh, p, sa, ...) are ignored at
3216	!-- present, maybe revise later.
3217	!-- Resetting of e is removed as unnecessary since e is not
3218	!-- anterpolated, and as incorrect since it overran the default
3219	!-- Neumann condition (bc_e_b).
3220	DO i = nxlg, nxrg
3221	DO j = nysg, nyng
3222	DO k = nzb, nzt+1
3223	u(k,j,i) = MERGE( u(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 1 ) )
3224	v(k,j,i) = MERGE( v(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 2 ) )
3225	w(k,j,i) = MERGE( w(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 3 ) )
3226	!
3227	!-- TO_DO: zero setting of temperature within topography creates
3228	!-- wrong results
3229	! pt(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp
3230	! IF ( humidity .OR. passive_scalar ) THEN
3231	! q(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp
3232	! ENDIF
3233	ENDDO
3234	ENDDO
3235	ENDDO
3236	ENDIF
3237	ENDIF
3238	ENDDO ! m
3239
3240	#endif
3241	END SUBROUTINE pmci_parent_datatrans
3242
3243
3244
3245	SUBROUTINE pmci_child_datatrans( direction )
3246
3247	IMPLICIT NONE
3248
3249	INTEGER(iwp), INTENT(IN) :: direction !< Transfer direction: parent_to_child or child_to_parent
3250
3251	#if defined( __parallel )
3252
3253	REAL(wp), DIMENSION(1) :: dtl !< Time step size
3254
3255
3256	dtl = dt_3d
3257	IF ( cpl_id > 1 ) THEN
3258
3259	IF ( direction == parent_to_child ) THEN
3260
3261	CALL cpu_log( log_point_s(73), 'pmc child recv', 'start' )
3262	CALL pmc_c_getbuffer( )
3263	CALL cpu_log( log_point_s(73), 'pmc child recv', 'stop' )
3264
3265	CALL cpu_log( log_point_s(75), 'pmc interpolation', 'start' )
3266	CALL pmci_interpolation
3267	CALL cpu_log( log_point_s(75), 'pmc interpolation', 'stop' )
3268
3269	ELSE
3270	!
3271	!-- direction == child_to_parent
3272	CALL cpu_log( log_point_s(76), 'pmc anterpolation', 'start' )
3273	CALL pmci_anterpolation
3274	CALL cpu_log( log_point_s(76), 'pmc anterpolation', 'stop' )
3275
3276	CALL cpu_log( log_point_s(74), 'pmc child send', 'start' )
3277	CALL pmc_c_putbuffer( )
3278	CALL cpu_log( log_point_s(74), 'pmc child send', 'stop' )
3279
3280	ENDIF
3281	ENDIF
3282
3283	CONTAINS
3284
3285
3286	SUBROUTINE pmci_interpolation
3287
3288	!
3289	!-- A wrapper routine for all interpolation actions
3290
3291	IMPLICIT NONE
3292
3293	INTEGER(iwp) :: ibgp !< Index running over the nbgp boundary ghost points in i-direction
3294	INTEGER(iwp) :: jbgp !< Index running over the nbgp boundary ghost points in j-direction
3295	INTEGER(iwp) :: lb !< Running index for aerosol size bins
3296	INTEGER(iwp) :: lc !< Running index for aerosol mass bins
3297	INTEGER(iwp) :: lg !< Running index for salsa gases
3298	INTEGER(iwp) :: n !< Running index for number of chemical species
3299
3300	!
3301	!-- In case of vertical nesting no interpolation is needed for the
3302	!-- horizontal boundaries
3303	IF ( nesting_mode /= 'vertical' ) THEN
3304	!
3305	!-- Left border pe:
3306	IF ( bc_dirichlet_l ) THEN
3307
3308	CALL pmci_interp_1sto_lr( u, uc, kcto, jflo, jfuo, kflo, kfuo, 'l', 'u' )
3309	CALL pmci_interp_1sto_lr( v, vc, kcto, jflv, jfuv, kflo, kfuo, 'l', 'v' )
3310	CALL pmci_interp_1sto_lr( w, wc, kctw, jflo, jfuo, kflw, kfuw, 'l', 'w' )
3311
3312	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
3313	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
3314	.NOT. constant_diffusion ) ) THEN
3315	! CALL pmci_interp_1sto_lr( e, ec, kcto, jflo, jfuo, kflo, kfuo, 'l', 'e' )
3316	!
3317	!-- Interpolation of e is replaced by the Neumann condition.
3318	DO ibgp = -nbgp, -1
3319	e(nzb:nzt,nys:nyn,ibgp) = e(nzb:nzt,nys:nyn,0)
3320	ENDDO
3321
3322	ENDIF
3323
3324	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
3325	CALL pmci_interp_1sto_lr( diss, dissc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3326	ENDIF
3327
3328	IF ( .NOT. neutral ) THEN
3329	CALL pmci_interp_1sto_lr( pt, ptc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3330	ENDIF
3331
3332	IF ( humidity ) THEN
3333
3334	CALL pmci_interp_1sto_lr( q, q_c, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3335
3336	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3337	CALL pmci_interp_1sto_lr( qc, qcc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3338	CALL pmci_interp_1sto_lr( nc, ncc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3339	ENDIF
3340
3341	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3342	CALL pmci_interp_1sto_lr( qr, qrc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3343	CALL pmci_interp_1sto_lr( nr, nrc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3344	ENDIF
3345
3346	ENDIF
3347
3348	IF ( passive_scalar ) THEN
3349	CALL pmci_interp_1sto_lr( s, sc, kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3350	ENDIF
3351
3352	IF ( air_chemistry .AND. nest_chemistry ) THEN
3353	DO n = 1, nspec
3354	CALL pmci_interp_1sto_lr( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3355	kcto, jflo, jfuo, kflo, kfuo, 'l', 's' )
3356	ENDDO
3357	ENDIF
3358
3359	IF ( salsa .AND. nest_salsa ) THEN
3360	DO lb = 1, nbins_aerosol
3361	CALL pmci_interp_1sto_lr( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3362	kcto, jflo, jfuo, kflo, kfuo, 'l', 's')
3363	ENDDO
3364	DO lc = 1, nbins_aerosol * ncomponents_mass
3365	CALL pmci_interp_1sto_lr( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3366	kcto, jflo, jfuo, kflo, kfuo, 'l', 's')
3367	ENDDO
3368	IF ( .NOT. salsa_gases_from_chem ) THEN
3369	DO lg = 1, ngases_salsa
3370	CALL pmci_interp_1sto_lr( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3371	kcto, jflo, jfuo, kflo, kfuo, 'l', 's')
3372	ENDDO
3373	ENDIF
3374	ENDIF
3375
3376	ENDIF
3377	!
3378	!-- Right border pe
3379	IF ( bc_dirichlet_r ) THEN
3380
3381	CALL pmci_interp_1sto_lr( u, uc, kcto, jflo, jfuo, kflo, kfuo, 'r', 'u' )
3382	CALL pmci_interp_1sto_lr( v, vc, kcto, jflv, jfuv, kflo, kfuo, 'r', 'v' )
3383	CALL pmci_interp_1sto_lr( w, wc, kctw, jflo, jfuo, kflw, kfuw, 'r', 'w' )
3384
3385	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
3386	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
3387	.NOT. constant_diffusion ) ) THEN
3388	! CALL pmci_interp_1sto_lr( e, ec, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 'r', 'e' )
3389	!
3390	!-- Interpolation of e is replaced by the Neumann condition.
3391	DO ibgp = nx+1, nx+nbgp
3392	e(nzb:nzt,nys:nyn,ibgp) = e(nzb:nzt,nys:nyn,nx)
3393	ENDDO
3394	ENDIF
3395
3396	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
3397	CALL pmci_interp_1sto_lr( diss, dissc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3398	ENDIF
3399
3400	IF ( .NOT. neutral ) THEN
3401	CALL pmci_interp_1sto_lr( pt, ptc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3402	ENDIF
3403
3404	IF ( humidity ) THEN
3405	CALL pmci_interp_1sto_lr( q, q_c, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3406
3407	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3408	CALL pmci_interp_1sto_lr( qc, qcc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3409	CALL pmci_interp_1sto_lr( nc, ncc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3410	ENDIF
3411
3412	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3413	CALL pmci_interp_1sto_lr( qr, qrc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3414	CALL pmci_interp_1sto_lr( nr, nrc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3415	ENDIF
3416
3417	ENDIF
3418
3419	IF ( passive_scalar ) THEN
3420	CALL pmci_interp_1sto_lr( s, sc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3421	ENDIF
3422
3423	IF ( air_chemistry .AND. nest_chemistry ) THEN
3424	DO n = 1, nspec
3425	CALL pmci_interp_1sto_lr( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3426	kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3427	ENDDO
3428	ENDIF
3429
3430	IF ( salsa .AND. nest_salsa ) THEN
3431	DO lb = 1, nbins_aerosol
3432	CALL pmci_interp_1sto_lr( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3433	kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3434	ENDDO
3435	DO lc = 1, nbins_aerosol * ncomponents_mass
3436	CALL pmci_interp_1sto_lr( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3437	kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3438	ENDDO
3439	IF ( .NOT. salsa_gases_from_chem ) THEN
3440	DO lg = 1, ngases_salsa
3441	CALL pmci_interp_1sto_lr( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3442	kcto, jflo, jfuo, kflo, kfuo, 'r', 's' )
3443	ENDDO
3444	ENDIF
3445	ENDIF
3446
3447	ENDIF
3448	!
3449	!-- South border pe
3450	IF ( bc_dirichlet_s ) THEN
3451
3452	CALL pmci_interp_1sto_sn( v, vc, kcto, iflo, ifuo, kflo, kfuo, 's', 'v' )
3453	CALL pmci_interp_1sto_sn( w, wc, kctw, iflo, ifuo, kflw, kfuw, 's', 'w' )
3454	CALL pmci_interp_1sto_sn( u, uc, kcto, iflu, ifuu, kflo, kfuo, 's', 'u' )
3455
3456	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
3457	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
3458	.NOT. constant_diffusion ) ) THEN
3459	! CALL pmci_interp_1sto_sn( e, ec, kcto, iflo, ifuo, kflo, kfuo, 's', 'e' )
3460	!
3461	!-- Interpolation of e is replaced by the Neumann condition.
3462	DO jbgp = -nbgp, -1
3463	e(nzb:nzt,jbgp,nxl:nxr) = e(nzb:nzt,0,nxl:nxr)
3464	ENDDO
3465	ENDIF
3466
3467	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
3468	CALL pmci_interp_1sto_sn( diss, dissc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3469	ENDIF
3470
3471	IF ( .NOT. neutral ) THEN
3472	CALL pmci_interp_1sto_sn( pt, ptc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3473	ENDIF
3474
3475	IF ( humidity ) THEN
3476	CALL pmci_interp_1sto_sn( q, q_c, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3477
3478	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3479	CALL pmci_interp_1sto_sn( qc, qcc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3480	CALL pmci_interp_1sto_sn( nc, ncc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3481	ENDIF
3482
3483	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3484	CALL pmci_interp_1sto_sn( qr, qrc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3485	CALL pmci_interp_1sto_sn( nr, nrc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3486	ENDIF
3487
3488	ENDIF
3489
3490	IF ( passive_scalar ) THEN
3491	CALL pmci_interp_1sto_sn( s, sc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3492	ENDIF
3493
3494	IF ( air_chemistry .AND. nest_chemistry ) THEN
3495	DO n = 1, nspec
3496	CALL pmci_interp_1sto_sn( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3497	kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3498	ENDDO
3499	ENDIF
3500
3501	IF ( salsa .AND. nest_salsa ) THEN
3502	DO lb = 1, nbins_aerosol
3503	CALL pmci_interp_1sto_sn( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3504	kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3505	ENDDO
3506	DO lc = 1, nbins_aerosol * ncomponents_mass
3507	CALL pmci_interp_1sto_sn( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3508	kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3509	ENDDO
3510	IF ( .NOT. salsa_gases_from_chem ) THEN
3511	DO lg = 1, ngases_salsa
3512	CALL pmci_interp_1sto_sn( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3513	kcto, iflo, ifuo, kflo, kfuo, 's', 's' )
3514	ENDDO
3515	ENDIF
3516	ENDIF
3517
3518	ENDIF
3519	!
3520	!-- North border pe
3521	IF ( bc_dirichlet_n ) THEN
3522
3523	CALL pmci_interp_1sto_sn( v, vc, kcto, iflo, ifuo, kflo, kfuo, 'n', 'v' )
3524	CALL pmci_interp_1sto_sn( w, wc, kctw, iflo, ifuo, kflw, kfuw, 'n', 'w' )
3525	CALL pmci_interp_1sto_sn( u, uc, kcto, iflu, ifuu, kflo, kfuo, 'n', 'u' )
3526
3527	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
3528	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
3529	.NOT. constant_diffusion ) ) THEN
3530	! CALL pmci_interp_1sto_sn( e, ec, kcto, iflo, ifuo, kflo, kfuo, 'n', 'e' )
3531	!
3532	!-- Interpolation of e is replaced by the Neumann condition.
3533	DO jbgp = ny+1, ny+nbgp
3534	e(nzb:nzt,jbgp,nxl:nxr) = e(nzb:nzt,ny,nxl:nxr)
3535	ENDDO
3536	ENDIF
3537
3538	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
3539	CALL pmci_interp_1sto_sn( diss, dissc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3540	ENDIF
3541
3542	IF ( .NOT. neutral ) THEN
3543	CALL pmci_interp_1sto_sn( pt, ptc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3544	ENDIF
3545
3546	IF ( humidity ) THEN
3547	CALL pmci_interp_1sto_sn( q, q_c, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3548
3549	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3550	CALL pmci_interp_1sto_sn( qc, qcc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3551	CALL pmci_interp_1sto_sn( nc, ncc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3552	ENDIF
3553
3554	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3555	CALL pmci_interp_1sto_sn( qr, qrc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3556	CALL pmci_interp_1sto_sn( nr, nrc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3557	ENDIF
3558
3559	ENDIF
3560
3561	IF ( passive_scalar ) THEN
3562	CALL pmci_interp_1sto_sn( s, sc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3563	ENDIF
3564
3565	IF ( air_chemistry .AND. nest_chemistry ) THEN
3566	DO n = 1, nspec
3567	CALL pmci_interp_1sto_sn( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3568	kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3569	ENDDO
3570	ENDIF
3571
3572	IF ( salsa .AND. nest_salsa ) THEN
3573	DO lb = 1, nbins_aerosol
3574	CALL pmci_interp_1sto_sn( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3575	kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3576	ENDDO
3577	DO lc = 1, nbins_aerosol * ncomponents_mass
3578	CALL pmci_interp_1sto_sn( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3579	kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3580	ENDDO
3581	IF ( .NOT. salsa_gases_from_chem ) THEN
3582	DO lg = 1, ngases_salsa
3583	CALL pmci_interp_1sto_sn( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3584	kcto, iflo, ifuo, kflo, kfuo, 'n', 's' )
3585	ENDDO
3586	ENDIF
3587	ENDIF
3588
3589	ENDIF
3590
3591	ENDIF ! IF ( nesting_mode /= 'vertical' )
3592	!
3593	!-- All PEs are top-border PEs
3594	CALL pmci_interp_1sto_t( w, wc, kctw, iflo, ifuo, jflo, jfuo, 'w' )
3595	CALL pmci_interp_1sto_t( u, uc, kcto, iflu, ifuu, jflo, jfuo, 'u' )
3596	CALL pmci_interp_1sto_t( v, vc, kcto, iflo, ifuo, jflv, jfuv, 'v' )
3597
3598	IF ( ( rans_mode_parent .AND. rans_mode ) .OR. &
3599	( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. &
3600	.NOT. constant_diffusion ) ) THEN
3601	! CALL pmci_interp_1sto_t( e, ec, kcto, iflo, ifuo, jflo, jfuo, 'e' )
3602	!
3603	!-- Interpolation of e is replaced by the Neumann condition.
3604	e(nzt+1,nys:nyn,nxl:nxr) = e(nzt,nys:nyn,nxl:nxr)
3605	ENDIF
3606
3607	IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN
3608	CALL pmci_interp_1sto_t( diss, dissc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3609	ENDIF
3610
3611	IF ( .NOT. neutral ) THEN
3612	CALL pmci_interp_1sto_t( pt, ptc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3613	ENDIF
3614
3615	IF ( humidity ) THEN
3616	CALL pmci_interp_1sto_t( q, q_c, kcto, iflo, ifuo, jflo, jfuo, 's' )
3617	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3618	CALL pmci_interp_1sto_t( qc, qcc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3619	CALL pmci_interp_1sto_t( nc, ncc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3620	ENDIF
3621	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3622	CALL pmci_interp_1sto_t( qr, qrc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3623	CALL pmci_interp_1sto_t( nr, nrc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3624	ENDIF
3625	ENDIF
3626
3627	IF ( passive_scalar ) THEN
3628	CALL pmci_interp_1sto_t( s, sc, kcto, iflo, ifuo, jflo, jfuo, 's' )
3629	ENDIF
3630
3631	IF ( air_chemistry .AND. nest_chemistry ) THEN
3632	DO n = 1, nspec
3633	CALL pmci_interp_1sto_t( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3634	kcto, iflo, ifuo, jflo, jfuo, 's' )
3635	ENDDO
3636	ENDIF
3637
3638	IF ( salsa .AND. nest_salsa ) THEN
3639	DO lb = 1, nbins_aerosol
3640	CALL pmci_interp_1sto_t( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3641	kcto, iflo, ifuo, jflo, jfuo, 's' )
3642	ENDDO
3643	DO lc = 1, nbins_aerosol * ncomponents_mass
3644	CALL pmci_interp_1sto_t( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3645	kcto, iflo, ifuo, jflo, jfuo, 's' )
3646	ENDDO
3647	IF ( .NOT. salsa_gases_from_chem ) THEN
3648	DO lg = 1, ngases_salsa
3649	CALL pmci_interp_1sto_t( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3650	kcto, iflo, ifuo, jflo, jfuo, 's' )
3651	ENDDO
3652	ENDIF
3653	ENDIF
3654
3655	END SUBROUTINE pmci_interpolation
3656
3657
3658
3659	SUBROUTINE pmci_anterpolation
3660
3661	!
3662	!-- A wrapper routine for all anterpolation actions.
3663	!-- Note that TKE is not anterpolated.
3664	IMPLICIT NONE
3665	INTEGER(iwp) :: lb !< Running index for aerosol size bins
3666	INTEGER(iwp) :: lc !< Running index for aerosol mass bins
3667	INTEGER(iwp) :: lg !< Running index for salsa gases
3668	INTEGER(iwp) :: n !< Running index for number of chemical species
3669
3670
3671	CALL pmci_anterp_tophat( u, uc, kcto, iflu, ifuu, jflo, jfuo, kflo, kfuo, ijkfc_u, 'u' )
3672	CALL pmci_anterp_tophat( v, vc, kcto, iflo, ifuo, jflv, jfuv, kflo, kfuo, ijkfc_v, 'v' )
3673	CALL pmci_anterp_tophat( w, wc, kctw, iflo, ifuo, jflo, jfuo, kflw, kfuw, ijkfc_w, 'w' )
3674	!
3675	!-- Anterpolation of TKE and dissipation rate if parent and child are in
3676	!-- RANS mode.
3677	IF ( rans_mode_parent .AND. rans_mode ) THEN
3678	CALL pmci_anterp_tophat( e, ec, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 'e' )
3679	!
3680	!-- Anterpolation of dissipation rate only if TKE-e closure is applied.
3681	IF ( rans_tke_e ) THEN
3682	CALL pmci_anterp_tophat( diss, dissc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, &
3683	ijkfc_s, 'diss' )
3684	ENDIF
3685
3686	ENDIF
3687
3688	IF ( .NOT. neutral ) THEN
3689	CALL pmci_anterp_tophat( pt, ptc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 'pt' )
3690	ENDIF
3691
3692	IF ( humidity ) THEN
3693
3694	CALL pmci_anterp_tophat( q, q_c, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 'q' )
3695
3696	IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN
3697
3698	CALL pmci_anterp_tophat( qc, qcc, kcto, iflo, ifuo, jflo, jfuo, &
3699	kflo, kfuo, ijkfc_s, 'qc' )
3700
3701	CALL pmci_anterp_tophat( nc, ncc, kcto, iflo, ifuo, jflo, jfuo, &
3702	kflo, kfuo, ijkfc_s, 'nc' )
3703
3704	ENDIF
3705
3706	IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN
3707
3708	CALL pmci_anterp_tophat( qr, qrc, kcto, iflo, ifuo, jflo, jfuo, &
3709	kflo, kfuo, ijkfc_s, 'qr' )
3710
3711	CALL pmci_anterp_tophat( nr, nrc, kcto, iflo, ifuo, jflo, jfuo, &
3712	kflo, kfuo, ijkfc_s, 'nr' )
3713
3714	ENDIF
3715
3716	ENDIF
3717
3718	IF ( passive_scalar ) THEN
3719	CALL pmci_anterp_tophat( s, sc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' )
3720	ENDIF
3721
3722	IF ( air_chemistry .AND. nest_chemistry ) THEN
3723	DO n = 1, nspec
3724	CALL pmci_anterp_tophat( chem_species(n)%conc, chem_spec_c(:,:,:,n), &
3725	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' )
3726	ENDDO
3727	ENDIF
3728
3729	IF ( salsa .AND. nest_salsa ) THEN
3730	DO lb = 1, nbins_aerosol
3731	CALL pmci_anterp_tophat( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), &
3732	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' )
3733	ENDDO
3734	DO lc = 1, nbins_aerosol * ncomponents_mass
3735	CALL pmci_anterp_tophat( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), &
3736	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' )
3737	ENDDO
3738	IF ( .NOT. salsa_gases_from_chem ) THEN
3739	DO lg = 1, ngases_salsa
3740	CALL pmci_anterp_tophat( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), &
3741	kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' )
3742	ENDDO
3743	ENDIF
3744	ENDIF
3745
3746	END SUBROUTINE pmci_anterpolation
3747
3748
3749
3750	SUBROUTINE pmci_interp_1sto_lr( child_array, parent_array, kct, jfl, jfu, kfl, kfu, edge, var )
3751	!
3752	!-- Interpolation of ghost-node values used as the child-domain boundary
3753	!-- conditions. This subroutine handles the left and right boundaries.
3754	IMPLICIT NONE
3755
3756	INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node
3757
3758	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain
3759	!< parent cell - y direction
3760	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain
3761	!< parent cell - y direction
3762	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain
3763	!< parent cell - z direction
3764	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain
3765	!< parent cell - z direction
3766
3767	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array
3768	REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array
3769
3770	CHARACTER(LEN=1), INTENT(IN) :: edge !< Edge symbol: 'l' or 'r'
3771	CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's'
3772	!
3773	!-- Local variables:
3774	INTEGER(iwp) :: icb !< Fixed child-grid index in the x-direction pointing to the node just behind the
3775	!< boundary-value node
3776	INTEGER(iwp) :: icbc !< Fixed child-grid index in the x-direction pointing to the boundary-value nodes
3777	INTEGER(iwp) :: icbgp !< Index running over the redundant boundary ghost points in the x-direction
3778	INTEGER(iwp) :: ierr !< MPI error code
3779	INTEGER(iwp) :: ipbeg !< Parent-grid index in the x-direction pointing to the starting point of workarr_lr
3780	!< in the parent-grid array
3781	INTEGER(iwp) :: ipw !< Reduced parent-grid index in the x-direction for workarr_lr pointing to
3782	!< the boundary ghost node
3783	INTEGER(iwp) :: ipwp !< Reduced parent-grid index in the x-direction for workarr_lr pointing to
3784	!< the first prognostic node
3785	INTEGER(iwp) :: jc !< Running child-grid index in the y-direction
3786	INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction
3787	INTEGER(iwp) :: kc !< Running child-grid index in the z-direction
3788	INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction
3789
3790	REAL(wp) :: cb !< Interpolation coefficient for the boundary ghost node
3791	REAL(wp) :: cp !< Interpolation coefficient for the first prognostic node
3792	REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in x direction from the parent-grid data
3793	REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in x direction from the parent-grid data
3794
3795	!
3796	!-- Check which edge is to be handled
3797	IF ( edge == 'l' ) THEN
3798	!
3799	!-- For u, nxl is a ghost node, but not for the other variables
3800	IF ( var == 'u' ) THEN
3801	icbc = nxl
3802	icb = icbc - 1
3803	ipw = 2
3804	ipwp = ipw ! This is redundant when var == 'u'
3805	ipbeg = ipl
3806	ELSE
3807	icbc = nxl - 1
3808	icb = icbc - 1
3809	ipw = 1
3810	ipwp = 2
3811	ipbeg = ipl
3812	ENDIF
3813	ELSEIF ( edge == 'r' ) THEN
3814	IF ( var == 'u' ) THEN
3815	icbc = nxr + 1
3816	icb = icbc + 1
3817	ipw = 1
3818	ipwp = ipw ! This is redundant when var == 'u'
3819	ipbeg = ipr - 2
3820	ELSE
3821	icbc = nxr + 1
3822	icb = icbc + 1
3823	ipw = 1
3824	ipwp = 0
3825	ipbeg = ipr - 2
3826	ENDIF
3827	ENDIF
3828	!
3829	!-- Interpolation coefficients
3830	IF ( interpolation_scheme_lrsn == 1 ) THEN
3831	cb = 1.0_wp ! 1st-order upwind
3832	ELSE IF ( interpolation_scheme_lrsn == 2 ) THEN
3833	cb = 0.5_wp ! 2nd-order central
3834	ELSE
3835	cb = 0.5_wp ! 2nd-order central (default)
3836	ENDIF
3837	cp = 1.0_wp - cb
3838	!
3839	!-- Substitute the necessary parent-grid data to the work array workarr_lr.
3840	workarr_lr = 0.0_wp
3841	IF ( pdims(2) > 1 ) THEN
3842	#if defined( __parallel )
3843	IF ( bc_dirichlet_s ) THEN
3844	workarr_lr(0:pg%nz+1,jpsw:jpnw-1,0:2) = parent_array(0:pg%nz+1,jpsw:jpnw-1,ipbeg:ipbeg+2)
3845	ELSE IF ( bc_dirichlet_n ) THEN
3846	workarr_lr(0:pg%nz+1,jpsw+1:jpnw,0:2) = parent_array(0:pg%nz+1,jpsw+1:jpnw,ipbeg:ipbeg+2)
3847	ELSE
3848	workarr_lr(0:pg%nz+1,jpsw+1:jpnw-1,0:2) &
3849	= parent_array(0:pg%nz+1,jpsw+1:jpnw-1,ipbeg:ipbeg+2)
3850	ENDIF
3851	!
3852	!-- South-north exchange if more than one PE subdomain in the y-direction.
3853	!-- Note that in case of 3-D nesting the south (psouth == MPI_PROC_NULL)
3854	!-- and north (pnorth == MPI_PROC_NULL) boundaries are not exchanged
3855	!-- because the nest domain is not cyclic.
3856	!-- From south to north
3857	CALL MPI_SENDRECV( workarr_lr(0,jpsw+1,0), 1, workarr_lr_exchange_type, psouth, 0, &
3858	workarr_lr(0,jpnw,0), 1, workarr_lr_exchange_type, pnorth, 0, comm2d, &
3859	status, ierr )
3860	!
3861	!-- From north to south
3862	CALL MPI_SENDRECV( workarr_lr(0,jpnw-1,0), 1, workarr_lr_exchange_type, pnorth, 1, &
3863	workarr_lr(0,jpsw,0), 1, workarr_lr_exchange_type, psouth, 1, comm2d, &
3864	status, ierr )
3865	#endif
3866	ELSE
3867	workarr_lr(0:pg%nz+1,jpsw:jpnw,0:2) = parent_array(0:pg%nz+1,jpsw:jpnw,ipbeg:ipbeg+2)
3868	ENDIF
3869
3870	IF ( var == 'u' ) THEN
3871
3872	DO jp = jpsw, jpnw
3873	DO kp = 0, kct
3874
3875	DO jc = jfl(jp), jfu(jp)
3876	DO kc = kfl(kp), kfu(kp)
3877	child_array(kc,jc,icbc) = workarr_lr(kp,jp,ipw)
3878	ENDDO
3879	ENDDO
3880
3881	ENDDO
3882	ENDDO
3883
3884	ELSE IF ( var == 'v' ) THEN
3885
3886	DO jp = jpsw, jpnw-1
3887	DO kp = 0, kct
3888	!
3889	!-- First interpolate to the flux point
3890	c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp)
3891	c_interp_2 = cb * workarr_lr(kp,jp+1,ipw) + cp * workarr_lr(kp,jp+1,ipwp)
3892	!
3893	!-- Use averages of the neighbouring matching grid-line values
3894	DO jc = jfl(jp), jfl(jp+1)
3895	child_array(kfl(kp):kfu(kp),jc,icbc) = 0.5_wp * ( c_interp_1 + c_interp_2 )
3896	ENDDO
3897	!
3898	!-- Then set the values along the matching grid-lines
3899	IF ( MOD( jfl(jp), jgsr ) == 0 ) THEN
3900	child_array(kfl(kp):kfu(kp),jfl(jp),icbc) = c_interp_1
3901	ENDIF
3902	ENDDO
3903	ENDDO
3904	!
3905	!-- Finally, set the values along the last matching grid-line
3906	IF ( MOD( jfl(jpnw), jgsr ) == 0 ) THEN
3907	DO kp = 0, kct
3908	c_interp_1 = cb * workarr_lr(kp,jpnw,ipw) + cp * workarr_lr(kp,jpnw,ipwp)
3909	child_array(kfl(kp):kfu(kp),jfl(jpnw),icbc) = c_interp_1
3910	ENDDO
3911	ENDIF
3912	!
3913	!-- A gap may still remain in some cases if the subdomain size is not
3914	!-- divisible by the grid-spacing ratio. In such a case, fill the
3915	!-- gap. Note however, this operation may produce some additional
3916	!-- momentum conservation error.
3917	IF ( jfl(jpnw) < nyn ) THEN
3918	DO kp = 0, kct
3919	DO jc = jfl(jpnw) + 1, nyn
3920	child_array(kfl(kp):kfu(kp),jc,icbc) = child_array(kfl(kp):kfu(kp),jfl(jpnw),icbc)
3921	ENDDO
3922	ENDDO
3923	ENDIF
3924
3925	ELSE IF ( var == 'w' ) THEN
3926
3927	DO jp = jpsw, jpnw
3928	DO kp = 0, kct + 1 ! It is important to go up to kct+1
3929	!
3930	!-- Interpolate to the flux point
3931	c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp)
3932	!
3933	!-- First substitute only the matching-node values
3934	child_array(kfu(kp),jfl(jp):jfu(jp),icbc) = c_interp_1
3935
3936	ENDDO
3937	ENDDO
3938
3939	DO jp = jpsw, jpnw
3940	DO kp = 1, kct + 1 ! It is important to go up to kct+1
3941	!
3942	!-- Then fill up the nodes in between with the averages
3943	DO kc = kfu(kp-1) + 1, kfu(kp) - 1
3944	child_array(kc,jfl(jp):jfu(jp),icbc) = &
3945	0.5_wp * ( child_array(kfu(kp-1),jfl(jp):jfu(jp),icbc) &
3946	+ child_array(kfu(kp),jfl(jp):jfu(jp),icbc) )
3947	ENDDO
3948
3949	ENDDO
3950	ENDDO
3951
3952	ELSE ! any scalar
3953
3954	DO jp = jpsw, jpnw
3955	DO kp = 0, kct
3956	!
3957	!-- Interpolate to the flux point
3958	c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp)
3959	DO jc = jfl(jp), jfu(jp)
3960	DO kc = kfl(kp), kfu(kp)
3961	child_array(kc,jc,icbc) = c_interp_1
3962	ENDDO
3963	ENDDO
3964
3965	ENDDO
3966	ENDDO
3967
3968	ENDIF ! var
3969	!
3970	!-- Fill up also the redundant 2nd and 3rd ghost-node layers
3971	IF ( edge == 'l' ) THEN
3972	DO icbgp = -nbgp, icb
3973	child_array(0:nzt+1,nysg:nyng,icbgp) = child_array(0:nzt+1,nysg:nyng,icbc)
3974	ENDDO
3975	ELSEIF ( edge == 'r' ) THEN
3976	DO icbgp = icb, nx+nbgp
3977	child_array(0:nzt+1,nysg:nyng,icbgp) = child_array(0:nzt+1,nysg:nyng,icbc)
3978	ENDDO
3979	ENDIF
3980
3981	END SUBROUTINE pmci_interp_1sto_lr
3982
3983
3984
3985	SUBROUTINE pmci_interp_1sto_sn( child_array, parent_array, kct, ifl, ifu, kfl, kfu, edge, var )
3986	!
3987	!-- Interpolation of ghost-node values used as the child-domain boundary
3988	!-- conditions. This subroutine handles the south and north boundaries.
3989	IMPLICIT NONE
3990
3991	INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node
3992
3993	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain
3994	!< parent cell - x direction
3995	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain
3996	!< parent cell - x direction
3997	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain
3998	!< parent cell - z direction
3999	INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain
4000	!< parent cell - z direction
4001
4002	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array
4003	REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array
4004
4005	CHARACTER(LEN=1), INTENT(IN) :: edge !< Edge symbol: 's' or 'n'
4006	CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's'
4007	!
4008	!-- Local variables:
4009	INTEGER(iwp) :: ic !< Running child-grid index in the x-direction
4010	INTEGER(iwp) :: ierr !< MPI error code
4011	INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction
4012	INTEGER(iwp) :: jcb !< Fixed child-grid index in the y-direction pointing to the node just behind the
4013	!< boundary-value node
4014	INTEGER(iwp) :: jcbc !< Fixed child-grid index in the y-direction pointing to the boundary-value nodes
4015	INTEGER(iwp) :: jcbgp !< Index running over the redundant boundary ghost points in y-direction
4016	INTEGER(iwp) :: kc !< Running child-grid index in the z-direction
4017	INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction
4018	INTEGER(iwp) :: jpbeg !< Parent-grid index in the y-direction pointing to the starting point of workarr_sn
4019	!< in the parent-grid array
4020	INTEGER(iwp) :: jpw !< Reduced parent-grid index in the y-direction for workarr_sn pointing to
4021	!< the boundary ghost node
4022	INTEGER(iwp) :: jpwp !< Reduced parent-grid index in the y-direction for workarr_sn pointing to
4023	!< the first prognostic node
4024	REAL(wp) :: cb !< Interpolation coefficient for the boundary ghost node
4025	REAL(wp) :: cp !< Interpolation coefficient for the first prognostic node
4026	REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in x direction from the parent-grid data
4027	REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in x direction from the parent-grid data
4028
4029	!
4030	!-- Check which edge is to be handled: south or north
4031	IF ( edge == 's' ) THEN
4032	!
4033	!-- For v, nys is a ghost node, but not for the other variables
4034	IF ( var == 'v' ) THEN
4035	jcbc = nys
4036	jcb = jcbc - 1
4037	jpw = 2
4038	jpwp = 2 ! This is redundant when var == 'v'
4039	jpbeg = jps
4040	ELSE
4041	jcbc = nys - 1
4042	jcb = jcbc - 1
4043	jpw = 1
4044	jpwp = 2
4045	jpbeg = jps
4046	ENDIF
4047	ELSEIF ( edge == 'n' ) THEN
4048	IF ( var == 'v' ) THEN
4049	jcbc = nyn + 1
4050	jcb = jcbc + 1
4051	jpw = 1
4052	jpwp = 0 ! This is redundant when var == 'v'
4053	jpbeg = jpn - 2
4054	ELSE
4055	jcbc = nyn + 1
4056	jcb = jcbc + 1
4057	jpw = 1
4058	jpwp = 0
4059	jpbeg = jpn - 2
4060	ENDIF
4061	ENDIF
4062	!
4063	!-- Interpolation coefficients
4064	IF ( interpolation_scheme_lrsn == 1 ) THEN
4065	cb = 1.0_wp ! 1st-order upwind
4066	ELSE IF ( interpolation_scheme_lrsn == 2 ) THEN
4067	cb = 0.5_wp ! 2nd-order central
4068	ELSE
4069	cb = 0.5_wp ! 2nd-order central (default)
4070	ENDIF
4071	cp = 1.0_wp - cb
4072	!
4073	!-- Substitute the necessary parent-grid data to the work array workarr_sn.
4074	workarr_sn = 0.0_wp
4075	IF ( pdims(1) > 1 ) THEN
4076	#if defined( __parallel )
4077	IF ( bc_dirichlet_l ) THEN
4078	workarr_sn(0:pg%nz+1,0:2,iplw:iprw-1) = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw:iprw-1)
4079	ELSE IF ( bc_dirichlet_r ) THEN
4080	workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw) = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw)
4081	ELSE
4082	workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw-1) &
4083	= parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw-1)
4084	ENDIF
4085	!
4086	!-- Left-right exchange if more than one PE subdomain in the x-direction.
4087	!-- Note that in case of 3-D nesting the left (pleft == MPI_PROC_NULL) and
4088	!-- right (pright == MPI_PROC_NULL) boundaries are not exchanged because
4089	!-- the nest domain is not cyclic.
4090	!-- From left to right
4091	CALL MPI_SENDRECV( workarr_sn(0,0,iplw+1), 1, workarr_sn_exchange_type, pleft, 0, &
4092	workarr_sn(0,0,iprw), 1, workarr_sn_exchange_type, pright, 0, comm2d, &
4093	status, ierr )
4094	!
4095	!-- From right to left
4096	CALL MPI_SENDRECV( workarr_sn(0,0,iprw-1), 1, workarr_sn_exchange_type, pright, 1, &
4097	workarr_sn(0,0,iplw), 1, workarr_sn_exchange_type, pleft, 1, comm2d, &
4098	status, ierr )
4099	#endif
4100	ELSE
4101	workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw-1) &
4102	= parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw-1)
4103	ENDIF
4104
4105	IF ( var == 'v' ) THEN
4106
4107	DO ip = iplw, iprw
4108	DO kp = 0, kct
4109
4110	DO ic = ifl(ip), ifu(ip)
4111	DO kc = kfl(kp), kfu(kp)
4112	child_array(kc,jcbc,ic) = workarr_sn(kp,jpw,ip)
4113	ENDDO
4114	ENDDO
4115
4116	ENDDO
4117	ENDDO
4118
4119	ELSE IF ( var == 'u' ) THEN
4120
4121	DO ip = iplw, iprw - 1
4122	DO kp = 0, kct
4123	!
4124	!-- First interpolate to the flux point
4125	c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip)
4126	c_interp_2 = cb * workarr_sn(kp,jpw,ip+1) + cp * workarr_sn(kp,jpwp,ip+1)
4127	!
4128	!-- Use averages of the neighbouring matching grid-line values
4129	DO ic = ifl(ip), ifl(ip+1)
4130	child_array(kfl(kp):kfu(kp),jcbc,ic) = 0.5_wp * ( c_interp_1 + c_interp_2 )
4131	ENDDO
4132	!
4133	!-- Then set the values along the matching grid-lines
4134	IF ( MOD( ifl(ip), igsr ) == 0 ) THEN
4135	child_array(kfl(kp):kfu(kp),jcbc,ifl(ip)) = c_interp_1
4136	ENDIF
4137
4138	ENDDO
4139	ENDDO
4140	!
4141	!-- Finally, set the values along the last matching grid-line
4142	IF ( MOD( ifl(iprw), igsr ) == 0 ) THEN
4143	DO kp = 0, kct
4144	c_interp_1 = cb * workarr_sn(kp,jpw,iprw) + cp * workarr_sn(kp,jpwp,iprw)
4145	child_array(kfl(kp):kfu(kp),jcbc,ifl(iprw)) = c_interp_1
4146	ENDDO
4147	ENDIF
4148	!
4149	!-- A gap may still remain in some cases if the subdomain size is not
4150	!-- divisible by the grid-spacing ratio. In such a case, fill the
4151	!-- gap. Note however, this operation may produce some additional
4152	!-- momentum conservation error.
4153	IF ( ifl(iprw) < nxr ) THEN
4154	DO kp = 0, kct
4155	DO ic = ifl(iprw) + 1, nxr
4156	child_array(kfl(kp):kfu(kp),jcbc,ic) = child_array(kfl(kp):kfu(kp),jcbc,ifl(iprw))
4157	ENDDO
4158	ENDDO
4159	ENDIF
4160
4161	ELSE IF ( var == 'w' ) THEN
4162
4163	DO ip = iplw, iprw
4164	DO kp = 0, kct + 1 ! It is important to go up to kct+1
4165	!
4166	!-- Interpolate to the flux point
4167	c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip)
4168	!
4169	!-- First substitute only the matching-node values
4170	child_array(kfu(kp),jcbc,ifl(ip):ifu(ip)) = c_interp_1
4171
4172	ENDDO
4173	ENDDO
4174
4175	DO ip = iplw, iprw
4176	DO kp = 1, kct + 1 ! It is important to go up to kct + 1
4177	!
4178	!-- Then fill up the nodes in between with the averages
4179	DO kc = kfu(kp-1) + 1, kfu(kp) - 1
4180	child_array(kc,jcbc,ifl(ip):ifu(ip)) = &
4181	0.5_wp * ( child_array(kfu(kp-1),jcbc,ifl(ip):ifu(ip)) &
4182	+ child_array(kfu(kp),jcbc,ifl(ip):ifu(ip)) )
4183	ENDDO
4184
4185	ENDDO
4186	ENDDO
4187
4188	ELSE ! Any scalar
4189
4190	DO ip = iplw, iprw
4191	DO kp = 0, kct
4192	!
4193	!-- Interpolate to the flux point
4194	c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip)
4195	DO ic = ifl(ip), ifu(ip)
4196	DO kc = kfl(kp), kfu(kp)
4197	child_array(kc,jcbc,ic) = c_interp_1
4198	ENDDO
4199	ENDDO
4200
4201	ENDDO
4202	ENDDO
4203
4204	ENDIF ! var
4205	!
4206	!-- Fill up also the redundant 2nd and 3rd ghost-node layers
4207	IF ( edge == 's' ) THEN
4208	DO jcbgp = -nbgp, jcb
4209	child_array(0:nzt+1,jcbgp,nxlg:nxrg) = child_array(0:nzt+1,jcbc,nxlg:nxrg)
4210	ENDDO
4211	ELSEIF ( edge == 'n' ) THEN
4212	DO jcbgp = jcb, ny+nbgp
4213	child_array(0:nzt+1,jcbgp,nxlg:nxrg) = child_array(0:nzt+1,jcbc,nxlg:nxrg)
4214	ENDDO
4215	ENDIF
4216
4217	END SUBROUTINE pmci_interp_1sto_sn
4218
4219
4220
4221	SUBROUTINE pmci_interp_1sto_t( child_array, parent_array, kct, ifl, ifu, jfl, jfu, var )
4222	!
4223	!-- Interpolation of ghost-node values used as the child-domain boundary
4224	!-- conditions. This subroutine handles the top boundary.
4225	IMPLICIT NONE
4226
4227	INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node
4228
4229	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain
4230	!< parent cell - x direction
4231	INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain
4232	!< parent cell - x direction
4233	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain
4234	!< parent cell - y direction
4235	INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain
4236	!< parent cell - y direction
4237
4238	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array
4239	REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array
4240
4241	CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's'
4242	!
4243	!-- Local variables:
4244	INTEGER(iwp) :: ic !< Running child-grid index in the x-direction
4245	INTEGER(iwp) :: ierr !< MPI error code
4246	INTEGER(iwp) :: iplc !< Lower parent-grid index limit in the x-direction for copying parent-grid
4247	!< array data to workarr_t
4248	INTEGER(iwp) :: iprc !< Upper parent-grid index limit in the x-direction for copying parent-grid
4249	!< array data to workarr_t
4250	INTEGER(iwp) :: jc !< Running child-grid index in the y-direction
4251	INTEGER(iwp) :: jpsc !< Lower parent-grid index limit in the y-direction for copying parent-grid
4252	!< array data to workarr_t
4253	INTEGER(iwp) :: jpnc !< Upper parent-grid-index limit in the y-direction for copying parent-grid
4254	!< array data to workarr_t
4255	INTEGER(iwp) :: kc !< Vertical child-grid index fixed to the boundary-value level
4256	INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction
4257	INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction
4258	INTEGER(iwp) :: kpw !< Reduced parent-grid index in the z-direction for workarr_t pointing to
4259	!< the boundary ghost node
4260	REAL(wp) :: c31 !< Interpolation coefficient for the 3rd-order WS scheme
4261	REAL(wp) :: c32 !< Interpolation coefficient for the 3rd-order WS scheme
4262	REAL(wp) :: c33 !< Interpolation coefficient for the 3rd-order WS scheme
4263	REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in z direction from the parent-grid data
4264	REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in z direction from the parent-grid data
4265
4266
4267	IF ( var == 'w' ) THEN
4268	kc = nzt
4269	ELSE
4270	kc = nzt + 1
4271	ENDIF
4272	kpw = 1
4273	!
4274	!-- Interpolation coefficients
4275	IF ( interpolation_scheme_t == 1 ) THEN
4276	c31 = 0.0_wp ! 1st-order upwind
4277	c32 = 1.0_wp
4278	c33 = 0.0_wp
4279	ELSE IF ( interpolation_scheme_t == 2 ) THEN
4280	c31 = 0.5_wp ! 2nd-order central
4281	c32 = 0.5_wp
4282	c33 = 0.0_wp
4283	ELSE
4284	c31 = 2.0_wp / 6.0_wp ! 3rd-order WS upwind biased (default)
4285	c32 = 5.0_wp / 6.0_wp
4286	c33 = -1.0_wp / 6.0_wp
4287	ENDIF
4288	!
4289	!-- Substitute the necessary parent-grid data to the work array.
4290	!-- Note that the dimension of workarr_t is (0:2,jpsw:jpnw,iplw:iprw),
4291	!-- And the jc?w and ic?w-index bounds depend on the location of the PE-
4292	!-- subdomain relative to the side boundaries.
4293	iplc = iplw + 1
4294	iprc = iprw - 1
4295	jpsc = jpsw + 1
4296	jpnc = jpnw - 1
4297	IF ( bc_dirichlet_l ) THEN
4298	iplc = iplw
4299	ENDIF
4300	IF ( bc_dirichlet_r ) THEN
4301	iprc = iprw
4302	ENDIF
4303	IF ( bc_dirichlet_s ) THEN
4304	jpsc = jpsw
4305	ENDIF
4306	IF ( bc_dirichlet_n ) THEN
4307	jpnc = jpnw
4308	ENDIF
4309	workarr_t = 0.0_wp
4310	workarr_t(0:2,jpsc:jpnc,iplc:iprc) = parent_array(kct:kct+2,jpsc:jpnc,iplc:iprc)
4311	!
4312	!-- Left-right exchange if more than one PE subdomain in the x-direction.
4313	!-- Note that in case of 3-D nesting the lef