Home

Context Navigation

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

Last change on this file since 3943 was 3943, checked in by maronga, 5 years ago
bugfixes in urban surface model; output of greenz roof transpiration added/corrected; minor formatting improvements
Property svn:keywords set to `Id`
File size: 204.1 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |