Home

Context Navigation

source: palm/trunk/SOURCE/pmc_interface_mod.f90 @ 3909

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |