Home

Context Navigation

source: palm/trunk/SOURCE/pmc_particle_interface.f90 @ 2834

Last change on this file since 2834 was 2809, checked in by schwenkel, 7 years ago
Bugfix for gfortran: Replace the function C_SIZEOF with STORAGE_SIZE
Property svn:keywords set to `Id`
File size: 42.8 KB

Line
1	MODULE pmc_particle_interface
2
3	!------------------------------------------------------------------------------!
4	! This file is part of PALM.
5	!
6	! PALM is free software: you can redistribute it and/or modify it under the
7	! terms of the GNU General Public License as published by the Free Software
8	! Foundation, either version 3 of the License, or (at your option) any later
9	! version.
10	!
11	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
12	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
13	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
14	!
15	! You should have received a copy of the GNU General Public License along with
16	! PALM. If not, see <http://www.gnu.org/licenses/>.
17	!
18	! Copyright 1997-2018 Leibniz Universitaet Hannover
19	!------------------------------------------------------------------------------!
20	!
21	! Current revisions:
22	! ------------------
23	!
24	!
25	! Former revisions:
26	! -----------------!
27	! $Id: pmc_particle_interface.f90 2809 2018-02-15 09:55:58Z raasch $
28	! Bugfix for gfortran: Replace the function C_SIZEOF with STORAGE_SIZE.
29	!
30	!
31	! Initial Version (by Klaus Ketelsen)
32	!
33	!
34	! Description:
35	! ------------
36	! Introduce particle transfer in nested models.
37	! Interface to palm lpm model to handel particle transfer between parent and
38	! child model.
39	!------------------------------------------------------------------------------!
40
41	USE, INTRINSIC :: ISO_C_BINDING
42	USE MPI
43
44	USE kinds
45
46	USE pegrid, &
47	ONLY: myidx,myidy
48
49	USE indices, &
50	ONLY: nx, ny, nxl, nxr, nys, nyn, nzb, nzt, nxlg, nxrg, nysg, nyng, nbgp
51
52	USE grid_variables, &
53	ONLY: dx, dy
54
55	USE arrays_3d, &
56	ONLY: zw
57
58	USE control_parameters, &
59	ONLY: message_string
60
61	USE particle_attributes, &
62	ONLY: prt_count, particles, grid_particles, &
63	particle_type, number_of_particles, zero_particle, &
64	ibc_par_t, ibc_par_lr, ibc_par_ns, alloc_factor
65
66	USE lpm_pack_and_sort_mod
67
68	USE pmc_general, &
69	ONLY: pedef
70
71	USE pmc_parent, &
72	ONLY: children, pmc_s_fillbuffer, pmc_s_getdata_from_buffer, &
73	pmc_s_get_child_npes
74
75	USE pmc_child, &
76	ONLY: me, pmc_c_getbuffer, pmc_c_putbuffer
77
78	USE pmc_interface, &
79	ONLY: cpl_id, get_number_of_childs, nr_part, part_adr, nested_run, &
80	get_childid, get_child_edges, nr_partc, part_adrc, &
81	coarse_bound, coord_x, coord_y, cg, get_child_gridspacing, &
82	lower_left_coord_x, lower_left_coord_y
83
84	USE pmc_handle_communicator, &
85	ONLY: pmc_parent_for_child
86
87	USE pmc_mpi_wrapper, &
88	ONLY: pmc_send_to_parent, pmc_recv_from_child
89
90	USE lpm_exchange_horiz_mod, &
91	ONLY: realloc_particles_array
92
93	IMPLICIT NONE
94	PRIVATE
95	SAVE
96
97	TYPE coarse_particle_def
98	INTEGER(iwp) :: nr_particle
99
100	TYPE(particle_type),ALLOCATABLE,DIMENSION(:) :: parent_particles
101	END TYPE coarse_particle_def
102
103	INTEGER(iwp),PARAMETER :: Min_particles_per_column = 100 !<
104	INTEGER(iwp),PARAMETER :: max_nr_particle_in_rma_win = 100000 !<
105
106	INTEGER(iwp) :: nr_fine_in_coarse !< Number of fine grid cells in coarse grid (one direction)
107	INTEGER(iwp) :: particle_win_child !<
108
109	INTEGER(iwp),ALLOCATABLE,DIMENSION(:) :: particle_win_parent !<
110
111	TYPE(C_PTR), ALLOCATABLE,DIMENSION(:) :: buf_ptr !<
112
113	TYPE(particle_type), DIMENSION(:),POINTER :: particle_in_win !<
114
115	TYPE(coarse_particle_def),ALLOCATABLE,DIMENSION(:,:) :: coarse_particles !<
116
117
118
119	INTERFACE pmcp_g_init
120	MODULE PROCEDURE pmcp_g_init
121	END INTERFACE pmcp_g_init
122
123	INTERFACE pmcp_g_alloc_win
124	MODULE PROCEDURE pmcp_g_alloc_win
125	END INTERFACE pmcp_g_alloc_win
126
127	INTERFACE pmcp_c_get_particle_from_parent
128	MODULE PROCEDURE pmcp_c_get_particle_from_parent
129	END INTERFACE pmcp_c_get_particle_from_parent
130
131	INTERFACE pmcp_c_send_particle_to_parent
132	MODULE PROCEDURE pmcp_c_send_particle_to_parent
133	END INTERFACE pmcp_c_send_particle_to_parent
134
135	INTERFACE pmcp_p_fill_particle_win
136	MODULE PROCEDURE pmcp_p_fill_particle_win
137	END INTERFACE pmcp_p_fill_particle_win
138
139	INTERFACE pmcp_p_empty_particle_win
140	MODULE PROCEDURE pmcp_p_empty_particle_win
141	END INTERFACE pmcp_p_empty_particle_win
142
143	INTERFACE pmcp_g_print_number_of_particles
144	MODULE PROCEDURE pmcp_g_print_number_of_particles
145	END INTERFACE pmcp_g_print_number_of_particles
146
147	INTERFACE pmcp_p_delete_particles_in_fine_grid_area
148	MODULE PROCEDURE pmcp_p_delete_particles_in_fine_grid_area
149	END INTERFACE pmcp_p_delete_particles_in_fine_grid_area
150
151	PUBLIC pmcp_g_init, pmcp_g_alloc_win, pmcp_c_get_particle_from_parent
152	PUBLIC pmcp_c_send_particle_to_parent, pmcp_p_fill_particle_win, pmcp_g_print_number_of_particles
153	PUBLIC pmcp_p_empty_particle_win, pmcp_p_delete_particles_in_fine_grid_area
154
155	CONTAINS
156
157	!------------------------------------------------------------------------------!
158	! Description:
159	! ------------
160	!> general routine:
161	!> Initializing actions of the particle interface
162	!> check particle boundary conditions for the child models
163	!------------------------------------------------------------------------------!
164	SUBROUTINE pmcp_g_init
165
166	IMPLICIT NONE
167
168	INTEGER(iwp) :: nr_childs !< Number of child models of the current model
169
170	nr_childs = get_number_of_childs()
171	!
172	!-- Check if the current model has child models
173	IF ( nr_childs > 0 ) THEN
174	ALLOCATE( nr_part(nysg:nyng, nxlg:nxrg) )
175	ALLOCATE( part_adr(nysg:nyng, nxlg:nxrg) )
176	nr_part = 0
177	part_adr = 0
178	ENDIF
179
180	!
181	!-- Set the boundary conditions to nested for all non root (i.e child) models
182	IF ( cpl_id > 1 ) THEN
183
184	IF ( ibc_par_t /= 3 ) THEN
185	ibc_par_t = 3
186	message_string = 'In Child model: bc_par_t is automatically set to nested '
187	CALL message( 'pmcp_g_init ', 'PA0477', 0, 1, 0, 6, 0 ) ! PA number has to be adjusted
188	ENDIF
189
190	IF ( ibc_par_lr /= 3 ) THEN
191	ibc_par_lr = 3
192	message_string = 'In Child model: bc_par_lr is automatically set to nested '
193	CALL message( 'pmcp_g_init ', 'PA0478', 0, 1, 0, 6, 0 ) ! PA number has to be adjusted
194	ENDIF
195
196	IF ( ibc_par_ns /= 3 ) THEN
197	ibc_par_ns = 3
198	message_string = 'In Child model: bc_par_ns is automatically set to nested '
199	CALL message( 'pmcp_g_init ', 'PA0479', 0, 1, 0, 6, 0 ) ! PA number has to be adjusted
200	ENDIF
201
202	ENDIF
203
204	END SUBROUTINE pmcp_g_init
205	!------------------------------------------------------------------------------!
206	! Description:
207	! ------------
208	!> general routine:
209	!> allocate the MPI windows
210	!------------------------------------------------------------------------------!
211	SUBROUTINE pmcp_g_alloc_win
212
213	IMPLICIT NONE
214
215	INTEGER(iwp) :: m !< loop index
216	INTEGER(iwp) :: ierr !< error code
217	INTEGER(iwp) :: icl !< left boundary in coarse(parent) index space
218	INTEGER(iwp) :: icr !< right boundary in coarse(parent) index space
219	INTEGER(iwp) :: jcs !< south boundary in coarse(parent) index space
220	INTEGER(iwp) :: jcn !< north boundary in coarse(parent) index space
221	INTEGER(iwp) :: child_id !< Id of a child model
222	INTEGER(iwp) :: nr_childs !< Number of child models of the current model
223
224	INTEGER(KIND=MPI_ADDRESS_KIND) :: parsize !<
225	INTEGER(KIND=MPI_ADDRESS_KIND) :: winsize !<
226	TYPE(C_PTR), SAVE :: ptr !<
227
228	TYPE(particle_type),DIMENSION(:),POINTER :: win_buffer !<
229
230	INTEGER(iwp),DIMENSION(1) :: buf_shape !<
231
232	!
233	!-- If the model has a parent model prepare the structures for transfer
234	IF ( cpl_id > 1 ) THEN
235
236	parsize = STORAGE_SIZE(zero_particle)/8
237
238	CALL MPI_ALLOC_MEM( parsize , MPI_INFO_NULL, ptr, ierr )
239	buf_shape(1) = 1
240	CALL C_F_POINTER( ptr, win_buffer, buf_shape )
241	CALL MPI_WIN_CREATE( win_buffer, parsize, parsize, &
242	MPI_INFO_NULL, me%intra_comm, particle_win_child, &
243	ierr )
244
245	!
246	!-- Child domain boundaries in the parent index space
247	icl = coarse_bound(1)
248	icr = coarse_bound(2)
249	jcs = coarse_bound(3)
250	jcn = coarse_bound(4)
251
252	! ALLOCATE( coarse_particles(jcs:jcn, icl:icr) )
253
254	coarse_particles(:,:)%nr_particle = 0
255	ENDIF
256
257	!
258	!-- If the model has child models prepare the structures for transfer
259	nr_childs = get_number_of_childs()
260	IF ( nr_childs > 0 ) THEN
261	ALLOCATE( particle_win_parent(nr_childs) )
262	ALLOCATE( buf_ptr(nr_childs) )
263	DO m = 1, nr_childs
264	child_id = get_childid(m)
265	parsize = STORAGE_SIZE(zero_particle)/8
266
267	winsize = max_nr_particle_in_rma_win * parsize
268	CALL MPI_ALLOC_MEM( winsize , MPI_INFO_NULL, buf_ptr(m), ierr )
269	buf_shape(1) = max_nr_particle_in_rma_win
270	CALL C_F_POINTER( buf_ptr(m), win_buffer, buf_shape )
271	CALL MPI_WIN_CREATE( win_buffer, winsize, parsize, &
272	MPI_INFO_NULL, children(child_id)%intra_comm, &
273	particle_win_parent(m), ierr )
274	ENDDO
275	ENDIF
276
277	END SUBROUTINE pmcp_g_alloc_win
278	!------------------------------------------------------------------------------!
279	! Description:
280	! ------------
281	!> child routine:
282	!> Read/get particles out of the parent MPI window
283	!------------------------------------------------------------------------------!
284	SUBROUTINE pmcp_c_get_particle_from_parent
285
286	IMPLICIT NONE
287
288	INTEGER(iwp) :: i !< x grid index
289	INTEGER(iwp) :: icl !< left boundary in coarse(parent) index space
290	INTEGER(iwp) :: ierr !< error code
291	INTEGER(iwp) :: ij !< combined xy index for the buffer array
292	INTEGER(iwp) :: ip !< loop index (child PEs)
293	INTEGER(iwp) :: j !< y grid index
294	INTEGER(iwp) :: jcs !< south boundary in coarse(parent) index space
295	INTEGER(iwp) :: nr !< number of particles to receive from a parent box
296
297	INTEGER(KIND=MPI_ADDRESS_KIND) :: parsize !<
298	INTEGER(KIND=MPI_ADDRESS_KIND) :: target_disp !<
299
300	TYPE(pedef), POINTER :: ape !< TO_DO Klaus: give a description and better name of the variable
301
302	IF ( cpl_id > 1 ) THEN
303
304	CALL pmc_c_getbuffer( particle_transfer = .TRUE. ) !Get number of particle/column and offset in RMA window xx
305
306	!
307	!-- Wait for buffer to fill.
308	!
309	!-- The parent side (in pmc_s_fillbuffer) is filling the buffer in the MPI RMA window
310	!-- When the filling is complete, a MPI_BARRIER is called.
311	!-- The child is not allowd to access the parent-buffer before it is completely filled
312	!-- Synchronization is done implicitely in pmc_c_getbuffer and pmc_s_fillbuffer on the parent side
313
314	icl = coarse_bound(1)
315	jcs = coarse_bound(3)
316
317	DO ip = 1, me%inter_npes
318
319	ape => me%pes(ip)
320
321	DO ij = 1, ape%nrele
322	j = ape%locind(ij)%j + jcs - 1
323	i = ape%locind(ij)%i + icl - 1
324	nr = nr_partc(j,i)
325	IF ( nr > 0 ) THEN
326
327	CALL check_and_alloc_coarse_particle (i, j, nr)
328	parsize = STORAGE_SIZE(zero_particle)/8
329	target_disp = part_adrc(j,i) - 1
330	CALL MPI_WIN_LOCK( MPI_LOCK_SHARED , ip - 1, 0, particle_win_child, ierr )
331	CALL MPI_GET( coarse_particles(j,i)%parent_particles, &
332	nr * parsize, MPI_BYTE, ip - 1, &
333	target_disp, nr * parsize, &
334	MPI_BYTE, particle_win_child, ierr )
335	CALL MPI_WIN_UNLOCK( ip - 1, particle_win_child, ierr )
336	ENDIF
337	coarse_particles(j,i)%nr_particle = nr
338	ENDDO
339	ENDDO
340
341	CALL c_copy_particle_to_child_grid
342	ENDIF
343
344	END SUBROUTINE pmcp_c_get_particle_from_parent
345	!------------------------------------------------------------------------------!
346	! Description:
347	! ------------
348	!> child routine:
349	!> Write/put particles into the parent MPI window
350	!------------------------------------------------------------------------------!
351	SUBROUTINE pmcp_c_send_particle_to_parent
352
353	IMPLICIT NONE
354
355	INTEGER(iwp) :: disp_offset !<
356	INTEGER(iwp) :: i !< x loop index
357	INTEGER(iwp) :: icl !< left boundary in coarse(parent) index space
358	INTEGER(iwp) :: icr !< right boundary in coarse(parent) index space
359	INTEGER(iwp) :: ierr !< error code
360	INTEGER(iwp) :: ij !< combined xy index for the buffer array
361	INTEGER(iwp) :: ip !< loop index (child PEs)
362	INTEGER(iwp) :: j !< y loop index
363	INTEGER(iwp) :: jcs !< south boundary in coarse(parent) index space
364	INTEGER(iwp) :: jcn !< north boundary in coarse(parent) index space
365	INTEGER(iwp) :: max_nr_particle_per_pe !< maximum number of particles per PE (depending on grid apect ratio)
366	INTEGER(iwp) :: n !< shorter variable name for nr_fine_in_coarse
367	INTEGER(iwp) :: nr !< shorter variabel name for nr_partc
368	INTEGER(iwp) :: pe_offset !< offset index of the current PE
369
370	INTEGER(KIND=MPI_ADDRESS_KIND) :: parsize !<
371	INTEGER(KIND=MPI_ADDRESS_KIND) :: target_disp !<
372
373	REAL(wp) :: eps=0.00001 !< used in calculations to avoid rounding errors
374	REAL(wp) :: xx !< number of fine grid cells inside a coarse grid cell in x-direction
375	REAL(wp) :: yy !< number of fine grid cells inside a coarse grid cell in y-direction
376
377	TYPE(pedef), POINTER :: ape !< TO_DO Klaus: give a description and better name of the variable
378
379	! TYPE(particle_type) :: dummy_part !< dummy particle (needed for size calculations)
380
381
382	IF ( cpl_id > 1 ) THEN
383	CALL c_copy_particle_to_coarse_grid
384
385	!
386	!-- Child domain boundaries in the parent index space
387
388	icl = coarse_bound(1)
389	icr = coarse_bound(2)
390	jcs = coarse_bound(3)
391	jcn = coarse_bound(4)
392
393	nr_partc = 0
394
395	DO i = icl, icr
396	DO j = jcs, jcn
397	nr_partc(j,i) = coarse_particles(j,i)%nr_particle
398	ENDDO
399	ENDDO
400	part_adrc = 0
401
402	!
403	!-- compute number of fine grid cells in coarse grid (one direction)
404	xx = (cg%dx + eps) / dx ! +eps to avoid rounding error
405	yy = (cg%dy + eps) / dy
406	nr_fine_in_coarse = MAX( INT(xx), INT(yy) )
407
408	IF ( MOD( coord_x(0), cg%dx ) /= 0.0 .OR. MOD( coord_y(0), cg%dy ) /= 0.0 ) THEN
409	nr_fine_in_coarse = nr_fine_in_coarse + 1
410	ENDIF
411
412	!
413	!-- Assign a number to my child PE to select different areas in the RMA window on server side
414	!-- With this number a square of child PEs is defined which share the same coarse grid cells
415
416	n = nr_fine_in_coarse ! local variable n to make folloing statements shorter
417	pe_offset = MOD( myidx, n ) * n + MOD( myidy, n)
418	max_nr_particle_per_pe = max_nr_particle_in_rma_win / ( n * n )
419	disp_offset = pe_offset * max_nr_particle_per_pe
420	parsize = STORAGE_SIZE(zero_particle)/8
421	DO ip = 1, me%inter_npes
422
423	ape => me%pes(ip)
424
425	target_disp = disp_offset
426	DO ij = 1, ape%nrele
427	j = ape%locind(ij)%j + jcs - 1
428	i = ape%locind(ij)%i + icl - 1
429	nr = nr_partc(j,i)
430	IF( nr > 0 ) THEN
431	IF ( target_disp + nr - disp_offset >= max_nr_particle_per_pe ) THEN
432	WRITE(9,*) 'RMA window too small on child ', &
433	target_disp + nr - disp_offset, &
434	max_nr_particle_per_pe, max_nr_particle_in_rma_win
435	message_string = 'RMA window too small on child'
436	CALL message( 'pmci_create_child_arrays', 'PA0480', 3, 2, 0, 6, 0 ) ! PA number has to be adjusted
437	ENDIF
438	CALL MPI_WIN_LOCK( MPI_LOCK_SHARED , ip - 1, 0, particle_win_child, ierr )
439	CALL MPI_PUT( coarse_particles(j,i)%parent_particles, &
440	nr * parsize, MPI_BYTE, ip - 1, &
441	target_disp, nr * parsize , &
442	MPI_BYTE, particle_win_child, ierr )
443	CALL MPI_WIN_UNLOCK( ip - 1, particle_win_child, ierr )
444	part_adrc(j,i) = target_disp + 1
445	target_disp = target_disp + nr
446	ENDIF
447	ENDDO
448	ENDDO
449
450	CALL pmc_c_putbuffer ( particle_transfer = .TRUE. ) !Send new number of particle/column and offset to parent
451
452	ENDIF
453
454	END SUBROUTINE pmcp_c_send_particle_to_parent
455	!------------------------------------------------------------------------------!
456	! Description:
457	! ------------
458	!> parent routine:
459	!> write particles into the MPI window
460	!------------------------------------------------------------------------------!
461	SUBROUTINE pmcp_p_fill_particle_win
462
463	IMPLICIT NONE
464
465	LOGICAL :: active_particle !< Particles located in the fine/child grid area are marked as active (to be transferred)
466	LOGICAL,SAVE :: lfirst = .TRUE. !<
467
468	INTEGER(iwp) :: child_id !< id of the child model
469	INTEGER(iwp) :: i !< x grid box index
470	INTEGER(iwp) :: ierr !< error code
471	INTEGER(iwp) :: ij !< combined xy index for the buffer array
472	INTEGER(iwp) :: ip !< loop index (child PEs)
473	INTEGER(iwp) :: j !< y grid box index
474	INTEGER(iwp) :: k !< z grid box index
475	INTEGER(iwp) :: m !< loop index (number of childs)
476	INTEGER(iwp) :: n !< loop index (number of particles)
477	INTEGER(iwp) :: nr_part_col !< Number of particles to transfer per column
478	INTEGER(iwp) :: number_of_particles !<
479	INTEGER(iwp) :: pindex !<
480	INTEGER(iwp) :: tot_particle_count !< Total number of particles per child
481
482	REAL(wp) :: dx_child !< child grid spacing
483	REAL(wp) :: dy_child !< child grid spacing
484	REAL(wp) :: dz_child !< child grid spacing
485	REAL(wp) :: ny_coord !< north coordinate of child grid
486	REAL(wp) :: ny_coord_b !< north coordinate of child grid boundary
487	REAL(wp) :: lx_coord !< left coordinate of child grid
488	REAL(wp) :: lx_coord_b !< left coordinate of child grid boundary
489	REAL(wp) :: rx_coord !< right coordinate of child grid
490	REAL(wp) :: rx_coord_b !< right coordinate of child grid boundary
491	REAL(wp) :: sy_coord !< south coordinate of child grid
492	REAL(wp) :: sy_coord_b !< south coordinate of child grid boundary
493	REAL(wp) :: uz_coord !< upper coordinate of child grid
494	REAL(wp) :: uz_coord_b !< upper coordinate of child grid boundary
495	REAL(wp) :: x !< particle position
496	REAL(wp) :: xo !< origin of particle
497	REAL(wp) :: y !< particle position
498	REAL(wp) :: yo !< origin of particle
499	REAL(wp) :: z !< particle position
500
501	INTEGER(iwp),DIMENSION(1) :: buf_shape !<
502
503	TYPE(pedef), POINTER :: ape !< TO_DO Klaus: give a description and better name of the variable
504
505	DO m = 1, get_number_of_childs()
506
507	child_id = get_childid(m)
508
509	CALL get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, &
510	sy_coord, sy_coord_b, ny_coord, ny_coord_b, &
511	uz_coord, uz_coord_b)
512
513	CALL get_child_gridspacing( m, dx_child, dy_child, dz_child )
514
515	IF ( lfirst ) THEN
516	WRITE(9,'(a,5f10.2)') 'edges ',lx_coord,rx_coord,sy_coord,ny_coord,uz_coord
517	WRITE(9,'(a,5f10.2)') 'edges boundary ',lx_coord_b,rx_coord_b,sy_coord_b,ny_coord_b,uz_coord_b
518	WRITE(9,'(a,5f10.2)') 'child spacing ',dx_child, dy_child, dz_child,lower_left_coord_x,lower_left_coord_y
519	ENDIF
520	!
521	!-- reset values for every child
522	tot_particle_count = 0
523	nr_part = 0
524	part_adr = 0
525	pindex = 1
526
527	buf_shape(1) = max_nr_particle_in_rma_win
528	CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
529
530	DO ip = 1, children(child_id)%inter_npes
531
532	ape => children(child_id)%pes(ip)
533
534	nr_part_col = 0
535
536	DO ij = 1, ape%nrele
537
538	!
539	!-- Inside the PMC adressing of 3d arrays starts with 1
540	i = ape%locind(ij)%i + nxl - nbgp - 1
541	j = ape%locind(ij)%j + nys - nbgp - 1
542	nr_part_col = 0 ! Number of particles to transfer per column
543	part_adr(j,i) = pindex
544
545	DO k = nzb + 1, nzt
546	number_of_particles = prt_count(k,j,i)
547
548	IF ( number_of_particles <= 0 ) CYCLE
549
550	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
551
552	! Select particles within boundary area
553
554	DO n = 1, number_of_particles
555	x = particles(n)%x
556	y = particles(n)%y
557	z = particles(n)%z
558	!
559	!-- check if the particle is located in the fine grid area
560	active_particle = ((x > lx_coord .AND. x < rx_coord) .AND. &
561	(y > sy_coord .AND. y < ny_coord) .AND. &
562	(z > 0.000000001 .AND. z < uz_coord))
563	IF ( active_particle .AND. particles(n)%particle_mask ) THEN
564
565	particle_in_win(pindex) = particles(n)
566	!
567	!-- Change particle positions and origin relative to global origin
568	particle_in_win(pindex)%x = particle_in_win(pindex)%x + lower_left_coord_x
569	particle_in_win(pindex)%y = particle_in_win(pindex)%y + lower_left_coord_y
570	particle_in_win(pindex)%origin_x = particle_in_win(pindex)%origin_x + lower_left_coord_x
571	particle_in_win(pindex)%origin_y = particle_in_win(pindex)%origin_y + lower_left_coord_y
572
573	tot_particle_count = tot_particle_count + 1
574	nr_part_col = nr_part_col + 1
575	pindex = pindex + 1
576	IF ( pindex > max_nr_particle_in_rma_win ) THEN
577	WRITE(9,*) 'RMA window too small on parent ',pindex, max_nr_particle_in_rma_win
578	message_string = 'RMA window too small on parent'
579	CALL message( 'pmci_create_child_arrays', 'PA0481', 3, 2, 0, 6, 0 ) ! PA number has to be adjusted
580	ENDIF
581	END IF
582	ENDDO
583	ENDDO
584	nr_part(j,i) = nr_part_col
585	ENDDO
586	ENDDO
587
588	CALL pmc_s_fillbuffer( child_id, particle_transfer = .TRUE. )
589	ENDDO
590
591	lfirst = .FALSE.
592
593	END SUBROUTINE pmcp_p_fill_particle_win
594
595	!------------------------------------------------------------------------------!
596	! Description:
597	! ------------
598	!> parent routine:
599	!> delete particles from the MPI window
600	!------------------------------------------------------------------------------!
601	SUBROUTINE pmcp_p_empty_particle_win
602	IMPLICIT NONE
603
604	INTEGER(iwp) :: child_id !< model id of the child
605	INTEGER(iwp) :: ip !< loop index (child PEs)
606	INTEGER(iwp) :: m !< loop index (number of childs)
607
608	INTEGER(iwp),DIMENSION(1) :: buf_shape !<
609
610	DO m = 1, get_number_of_childs()
611
612	child_id = get_childid(m)
613
614	buf_shape(1) = max_nr_particle_in_rma_win
615	CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
616
617	!
618	!-- In some cells of the coarse grid, there are contributions from more than one child process
619	!-- Therfore p_copy_particle_to_org_grid is done for one child process per call
620
621	DO ip = 1, pmc_s_get_child_npes( child_id )
622	nr_part = 0
623	part_adr = 0
624
625	CALL pmc_s_getdata_from_buffer( child_id, particle_transfer = .TRUE., child_process_nr = ip )
626
627	CALL p_copy_particle_to_org_grid( m, child_id )
628	ENDDO
629
630	ENDDO
631
632	END SUBROUTINE pmcp_p_empty_particle_win
633
634	!------------------------------------------------------------------------------!
635	! Description:
636	! ------------
637	!> parent routine:
638	!> After the transfer mark all parent particles that are still inside on of the
639	!> child areas for deletion.
640	!------------------------------------------------------------------------------!
641	SUBROUTINE pmcp_p_delete_particles_in_fine_grid_area
642
643	IMPLICIT NONE
644
645	LOGICAL :: to_delete !< particles outside of model domain are marked as to_delete
646
647	INTEGER(iwp) :: i !< loop index (x grid)
648	INTEGER(iwp) :: j !< loop index (y grid)
649	INTEGER(iwp) :: k !< loop index (z grid)
650	INTEGER(iwp) :: m !< loop index (number of particles)
651	INTEGER(iwp) :: n !< loop index (number of childs)
652
653	REAL(wp) :: dx_child !< child grid spacing
654	REAL(wp) :: dy_child !< child grid spacing
655	REAL(wp) :: dz_child !< child grid spacing
656	REAL(wp) :: ny_coord !< north coordinate of child grid
657	REAL(wp) :: ny_coord_b !< north coordinate of child grid boundary
658	REAL(wp) :: lx_coord !< left coordinate of child grid
659	REAL(wp) :: lx_coord_b !< left coordinate of child grid boundary
660	REAL(wp) :: rx_coord !< right coordinate of child grid
661	REAL(wp) :: rx_coord_b !< right coordinate of child grid boundary
662	REAL(wp) :: sy_coord !< south coordinate of child grid
663	REAL(wp) :: sy_coord_b !< south coordinate of child grid boundary
664	REAL(wp) :: uz_coord !< upper coordinate of child grid
665	REAL(wp) :: uz_coord_b !< upper coordinate of child grid boundary
666	REAL(wp) :: x !< particle position
667	REAL(wp) :: y !< particle position
668	REAL(wp) :: z !< particle position
669
670	DO m = 1, get_number_of_childs()
671	CALL get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, &
672	sy_coord, sy_coord_b, ny_coord, ny_coord_b, &
673	uz_coord, uz_coord_b )
674
675
676	CALL get_child_gridspacing( m, dx_child, dy_child, dz_child )
677
678	DO i = nxl, nxr
679	DO j = nys, nyn
680	DO k = nzb, nzt
681	number_of_particles = prt_count(k,j,i)
682
683	IF ( number_of_particles == 0 ) CYCLE
684
685	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
686
687	DO n = 1, number_of_particles
688	x = particles(n)%x
689	y = particles(n)%y
690	z = particles(n)%z
691
692	to_delete = ((x > lx_coord .AND. x < rx_coord) .AND. &
693	(y > sy_coord .AND. y < ny_coord) .AND. &
694	(z > 0.000000001 .AND. z < uz_coord))
695
696	IF ( to_delete ) THEN
697	particles(n)%particle_mask = .FALSE.
698	ENDIF
699	ENDDO
700	ENDDO
701	ENDDO
702	ENDDO
703
704	ENDDO
705
706	END SUBROUTINE pmcp_p_delete_particles_in_fine_grid_area
707	!------------------------------------------------------------------------------!
708	! Description:
709	! ------------
710	!> general routine:
711	!> print the total number of of the current model and its child models into a file
712	!------------------------------------------------------------------------------!
713	SUBROUTINE pmcp_g_print_number_of_particles (my_time,local_nr_particles)
714
715	USE pegrid, &
716	ONLY: myid
717
718	IMPLICIT NONE
719
720	INTEGER(iwp),INTENT(IN) :: local_nr_particles !<
721
722	REAL(wp),INTENT(IN) :: my_time !<
723
724	LOGICAL, SAVE :: is_file_open=.FALSE. !<
725
726	INTEGER(iwp) :: child_id !<
727	INTEGER(iwp) :: child_nr_particles !< total number of particles in all child models
728	INTEGER(iwp) :: ierr !< error code
729	INTEGER(iwp) :: m !< loop index (number of childs)
730
731	CHARACTER(LEN=16) :: fname !< filename
732
733	INTEGER(iwp),DIMENSION(2) :: ivalr !< integer value to be received
734	INTEGER(iwp),DIMENSION(2) :: ivals !< integer value to be send
735
736
737	child_nr_particles = 0
738	IF ( myid == 0 ) THEN
739	IF ( cpl_id > 1 ) THEN
740	ivals(1) = local_nr_particles
741	CALL pmc_send_to_parent( ivals, 1, 0, 400, ierr )
742	ENDIF
743	DO m = 1, SIZE( pmc_parent_for_child ) - 1
744
745	child_id = pmc_parent_for_child(m)
746	CALL pmc_recv_from_child( child_id, ivalr, 1, 0, 400, ierr )
747	child_nr_particles = child_nr_particles + ivalr(1)
748	ENDDO
749
750	IF ( SIZE( pmc_parent_for_child ) > 1 ) THEN
751	IF ( .NOT. is_file_open ) THEN !kk muss noch auf file_open umgestellt werden
752	WRITE(fname,'(a,i2.2)') 'nr_particles_',cpl_id
753	OPEN (333,file = fname)
754	is_file_open = .true.
755	ENDIF
756	WRITE(333,'(f12.2,3i10)') my_time,local_nr_particles + child_nr_particles,local_nr_particles,child_nr_particles
757	ENDIF
758	ENDIF
759
760	END SUBROUTINE pmcp_g_print_number_of_particles
761
762	!------------------------------------------------------------------------------!
763	!------------------------------------------------------------------------------!
764	! Private subroutines
765	!------------------------------------------------------------------------------!
766	!------------------------------------------------------------------------------!
767
768	!------------------------------------------------------------------------------!
769	! Description:
770	! ------------
771	!> child routine
772	!> update the size of the local buffer (coarse_particles)
773	!------------------------------------------------------------------------------!
774	SUBROUTINE check_and_alloc_coarse_particle (ic,jc,nr,with_copy)
775
776	IMPLICIT NONE
777
778	INTEGER(iwp),INTENT(IN) :: ic !< coarse x grid index
779	INTEGER(iwp),INTENT(IN) :: jc !< coarse y grid index
780	INTEGER(iwp),INTENT(IN) :: nr !< number of particles to be transferred/stored in local buffer
781
782	LOGICAL,INTENT(IN),OPTIONAL :: with_copy !< copy particles in buffer? or reallocate empty buffer
783
784	LOGICAL :: with_copy_lo !< local variable of with copy
785
786	INTEGER(iwp) :: new_size !< new size of the local buffer
787	INTEGER(iwp) :: old_size !< old size of the local buffer
788
789	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: tmp_particles_d !<
790
791	with_copy_lo = .FALSE.
792	IF ( PRESENT( with_copy ) ) with_copy_lo = with_copy
793
794	IF ( .NOT. ALLOCATED( coarse_particles(jc,ic)%parent_particles ) ) THEN
795	new_size = MAX( nr, Min_particles_per_column )
796	ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
797	ELSEIF ( nr > SIZE( coarse_particles(jc,ic)%parent_particles ) ) THEN
798
799	old_size = SIZE( coarse_particles(jc,ic)%parent_particles )
800	new_size = old_size * ( 1.0_wp + alloc_factor / 100.0_wp )
801	new_size = MAX( nr, new_size, old_size + Min_particles_per_column )
802
803	!
804	!-- Copy existing particles to new particle buffer
805	IF ( with_copy_lo ) THEN
806	ALLOCATE( tmp_particles_d(old_size) )
807	tmp_particles_d(1:old_size) = coarse_particles(jc,ic)%parent_particles
808
809	DEALLOCATE( coarse_particles(jc,ic)%parent_particles )
810	ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
811
812	coarse_particles(jc,ic)%parent_particles(1:old_size) = tmp_particles_d(1:old_size)
813	coarse_particles(jc,ic)%parent_particles(old_size+1:new_size) = zero_particle
814
815	DEALLOCATE( tmp_particles_d )
816	!
817	!-- allocate or reallocate an empty buffer
818	ELSE
819	DEALLOCATE( coarse_particles(jc,ic)%parent_particles )
820	ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
821	ENDIF
822	ENDIF
823
824	END SUBROUTINE check_and_alloc_coarse_particle
825
826	!------------------------------------------------------------------------------!
827	! Description:
828	! ------------
829	!> child routine:
830	!> copy/sort particles out of the local buffer into the respective grid boxes
831	!------------------------------------------------------------------------------!
832	SUBROUTINE c_copy_particle_to_child_grid
833
834	IMPLICIT NONE
835
836	INTEGER(iwp) :: ic !< coarse x grid index
837	INTEGER(iwp) :: icl !< left boundary in coarse(parent) index space
838	INTEGER(iwp) :: icr !< right boundary in coarse(parent) index space
839	INTEGER(iwp) :: ip !< x grid index
840	INTEGER(iwp) :: jc !< coarse y grid index
841	INTEGER(iwp) :: jcn !< north boundary in coarse(parent) index space
842	INTEGER(iwp) :: jcs !< south boundary in coarse(parent) index space
843	INTEGER(iwp) :: jp !< y grid index
844	INTEGER(iwp) :: kp !< z grid index
845	INTEGER(iwp) :: n !< loop index (number of particles)
846	INTEGER(iwp) :: nr !< short variable for name number or particles
847
848	REAL(wp) :: xc !< child x coordinate
849	REAL(wp) :: xoc !< child x origin
850	REAL(wp) :: yc !< child y coordinate
851	REAL(wp) :: yoc !< child y origin
852	REAL(wp) :: zc !< child z coordinate
853
854	!
855	!-- Child domain boundaries in the parent index space
856	icl = coarse_bound(1)
857	icr = coarse_bound(2)
858	jcs = coarse_bound(3)
859	jcn = coarse_bound(4)
860
861	DO ic = icl, icr
862	DO jc = jcs, jcn
863	nr = coarse_particles(jc,ic)%nr_particle
864
865	IF ( nr > 0 ) THEN
866	DO n = 1, nr
867	xc = coarse_particles(jc,ic)%parent_particles(n)%x-lower_left_coord_x
868	yc = coarse_particles(jc,ic)%parent_particles(n)%y-lower_left_coord_y
869	zc = coarse_particles(jc,ic)%parent_particles(n)%z
870	xoc = coarse_particles(jc,ic)%parent_particles(n)%origin_x-lower_left_coord_x
871	yoc = coarse_particles(jc,ic)%parent_particles(n)%origin_y-lower_left_coord_y
872	ip = xc / dx
873	jp = yc / dy
874	kp = nzt
875	DO WHILE ( zw(kp-1) > zc .AND. kp > nzb + 1 ) ! kk search loop has to be optimzed !!!
876	kp = kp - 1
877	ENDDO
878
879	IF ( ip >= nxl .AND. ip <= nxr .AND. jp >= nys .AND. jp <= nyn ) THEN
880	prt_count(kp,jp,ip) = prt_count(kp,jp,ip) + 1
881	IF ( prt_count(kp,jp,ip) > SIZE( grid_particles(kp,jp,ip)%particles ) ) THEN
882	CALL realloc_particles_array( ip, jp, kp )
883	ENDIF
884	coarse_particles(jc,ic)%parent_particles(n)%x = xc ! Adjust coordinates to child grid
885	coarse_particles(jc,ic)%parent_particles(n)%y = yc
886	coarse_particles(jc,ic)%parent_particles(n)%origin_x = xoc ! Adjust origins to child grid
887	coarse_particles(jc,ic)%parent_particles(n)%origin_y = yoc
888	grid_particles(kp,jp,ip)%particles(prt_count(kp,jp,ip)) = coarse_particles(jc,ic)%parent_particles(n)
889	ENDIF
890	ENDDO
891	ENDIF
892	ENDDO
893	ENDDO
894
895	END SUBROUTINE c_copy_particle_to_child_grid
896	!------------------------------------------------------------------------------!
897	! Description:
898	! ------------
899	!> child routine:
900	!> copy particles which left the model area into the local buffer
901	!------------------------------------------------------------------------------!
902	SUBROUTINE c_copy_particle_to_coarse_grid
903
904	IMPLICIT NONE
905
906	LOGICAL :: boundary_particle !<
907
908	INTEGER(iwp) :: i !< loop index (x grid)
909	INTEGER(iwp) :: ic !< loop index (coarse x grid)
910	INTEGER(iwp) :: icl !< left boundary in coarse(parent) index space
911	INTEGER(iwp) :: icr !< left boundary in coarse(parent) index space
912	INTEGER(iwp) :: ierr !< error code
913	INTEGER(iwp) :: j !< loop index (y grid)
914	INTEGER(iwp) :: jc !< loop index (coarse y grid)
915	INTEGER(iwp) :: jcs !< south boundary in coarse(parent) index space
916	INTEGER(iwp) :: jcn !< north boundary in coarse(parent) index space
917	INTEGER(iwp) :: k !< loop index (z grid)
918	INTEGER(iwp) :: n !< loop index (number of particles)
919
920	REAL(iwp) :: x !< x coordinate
921	REAL(iwp) :: xo !< x origin
922	REAL(iwp) :: x_left !< absolute left boundary
923	REAL(iwp) :: x_right !< absolute right boundary
924	REAL(iwp) :: y !< left boundary in coarse(parent) index space
925	REAL(iwp) :: yo !< y origin
926	REAL(iwp) :: y_north !< absolute north boundary
927	REAL(iwp) :: y_south !< absoulte south boundary
928	REAL(iwp) :: z !< z coordinate
929
930	!
931	!-- Child domain boundaries in the parent index space
932
933	icl = coarse_bound(1)
934	icr = coarse_bound(2)
935	jcs = coarse_bound(3)
936	jcn = coarse_bound(4)
937
938	!
939	!-- absolute coordinates
940	x_left = coord_x(0)
941	x_right = coord_x(nx) + dx
942	y_south = coord_y(0)
943	y_north = coord_y(ny) + dy
944
945	! Clear Particle Buffer
946
947	DO ic = icl, icr
948	DO jc = jcs, jcn
949	coarse_particles(jc,ic)%nr_particle = 0
950	ENDDO
951	ENDDO
952
953	!
954	!-- Determine particles which leave the inner area in east or west dirextion
955	!-- Compute only first (nxl) and last (nxr) loop iterration.
956	DO i = nxl, nxr
957	DO j = nys, nyn
958	DO k = nzb + 1, nzt
959	number_of_particles = prt_count(k,j,i)
960	IF ( number_of_particles <= 0 ) CYCLE
961	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
962	DO n = 1, number_of_particles
963	IF ( particles(n)%particle_mask ) THEN
964	x = particles(n)%x+ lower_left_coord_x
965	y = particles(n)%y+ lower_left_coord_y
966	xo = particles(n)%origin_x + lower_left_coord_x
967	yo = particles(n)%origin_y + lower_left_coord_y
968	z = particles(n)%z
969
970	boundary_particle = .FALSE.
971	boundary_particle = boundary_particle .OR. (x < x_left .OR. x > x_right)
972	boundary_particle = boundary_particle .OR. (y < y_south .OR. y > y_north)
973	boundary_particle = boundary_particle .OR. (z > zw(nzt))
974
975	IF ( boundary_particle ) THEN
976	ic = x / cg%dx !TODO anpassen auf Mehrfachnesting
977	jc = y / cg%dy
978
979	IF ( ic >= icl .AND. ic <= icr .AND. jc >= jcs .AND. jc <= jcn ) THEN
980	coarse_particles(jc,ic)%nr_particle = coarse_particles(jc,ic)%nr_particle + 1
981	CALL check_and_alloc_coarse_particle( ic, jc, coarse_particles(jc,ic)%nr_particle, with_copy=.TRUE. )
982
983	coarse_particles(jc,ic)%parent_particles(coarse_particles(jc,ic)%nr_particle) = particles(n)
984	coarse_particles(jc,ic)%parent_particles(coarse_particles(jc,ic)%nr_particle)%x = x !adapt to global coordinates
985	coarse_particles(jc,ic)%parent_particles(coarse_particles(jc,ic)%nr_particle)%y = y
986	coarse_particles(jc,ic)%parent_particles(coarse_particles(jc,ic)%nr_particle)%origin_x = xo
987	coarse_particles(jc,ic)%parent_particles(coarse_particles(jc,ic)%nr_particle)%origin_y = yo
988	!-- Mark particle as deleted after copying it to the transfer buffer
989	grid_particles(k,j,i)%particles(n)%particle_mask = .FALSE.
990	ELSE
991	WRITE(9,'(a,10i6)') 'This should not happen ',i,j,k,ic,jc,icl,icr,jcs,jcn
992	CALL MPI_Abort( MPI_COMM_WORLD, 9999, ierr )
993	ENDIF
994	ENDIF
995	ENDIF
996	ENDDO
997	ENDDO
998	ENDDO
999	ENDDO
1000
1001	!
1002	!- Pack particles (eliminate those marked for deletion),
1003	!- determine new number of particles
1004	CALL lpm_sort_in_subboxes
1005
1006	END SUBROUTINE c_copy_particle_to_coarse_grid
1007	!------------------------------------------------------------------------------!
1008	! Description:
1009	! ------------
1010	!> parent routine:
1011	!> copy/sort particles from the MPI window into the respective grid boxes
1012	!------------------------------------------------------------------------------!
1013	SUBROUTINE p_copy_particle_to_org_grid( m, child_id )
1014
1015	IMPLICIT NONE
1016
1017	INTEGER(iwp),INTENT(IN) :: child_id !<
1018	INTEGER(iwp),INTENT(IN) :: m !<
1019
1020	INTEGER(iwp) :: i !< loop index (x grid)
1021	INTEGER(iwp) :: j !< loop index (y grid)
1022	INTEGER(iwp) :: k !< loop index (z grid)
1023	INTEGER(iwp) :: n !< loop index (nr part)
1024	INTEGER(iwp) :: nr !< short variable name for nr_part
1025	INTEGER(iwp) :: pindex !< short variable name part_adr
1026
1027	REAL(iwp) :: x !< x coordinate
1028	REAL(iwp) :: xo !< x origin
1029	REAL(iwp) :: y !< y coordinate
1030	REAL(iwp) :: yo !< y origin
1031	REAL(iwp) :: z !< z coordinate
1032
1033	INTEGER(iwp),DIMENSION(1) :: buf_shape !<
1034
1035	buf_shape(1) = max_nr_particle_in_rma_win
1036	CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
1037
1038	DO i = nxl, nxr
1039	DO j = nys, nyn
1040	nr = nr_part(j,i)
1041	IF ( nr > 0 ) THEN
1042	pindex = part_adr(j,i)
1043	DO n = 1, nr
1044	x = particle_in_win(pindex)%x-lower_left_coord_x
1045	y = particle_in_win(pindex)%y-lower_left_coord_y
1046	z = particle_in_win(pindex)%z
1047	xo = particle_in_win(pindex)%origin_x-lower_left_coord_x
1048	yo = particle_in_win(pindex)%origin_y-lower_left_coord_y
1049	k = nzt + 1
1050	DO WHILE ( zw(k-1) > z .AND. k > nzb + 1 ) ! kk search loop has to be optimzed !!!
1051	k = k - 1
1052	END DO
1053
1054	prt_count(k,j,i) = prt_count(k,j,i) + 1
1055	IF ( prt_count(k,j,i) > SIZE( grid_particles(k,j,i)%particles ) ) THEN
1056	CALL realloc_particles_array( i, j, k )
1057	ENDIF
1058	grid_particles(k,j,i)%particles(prt_count(k,j,i)) = particle_in_win(pindex)
1059
1060	!
1061	!-- Update particle positions and origins relative to parent domain
1062	grid_particles(k,j,i)%particles(prt_count(k,j,i))%x = x
1063	grid_particles(k,j,i)%particles(prt_count(k,j,i))%y = y
1064	grid_particles(k,j,i)%particles(prt_count(k,j,i))%origin_x = xo
1065	grid_particles(k,j,i)%particles(prt_count(k,j,i))%origin_y = yo
1066	pindex = pindex + 1
1067	ENDDO
1068	ENDIF
1069	ENDDO
1070	ENDDO
1071
1072	END SUBROUTINE p_copy_particle_to_org_grid
1073
1074
1075	END MODULE pmc_particle_interface

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |