Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |