Home

Context Navigation

source: palm/trunk/SOURCE/lpm_exchange_horiz.f90 @ 1989

Last change on this file since 1989 was 1937, checked in by suehring, 8 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 46.0 KB

Rev	Line
[1873]	1	!> @file lpm_exchange_horiz.f90
[1036]	2	!--------------------------------------------------------------------------------!
	3	! This file is part of PALM.
	4	!
	5	! PALM is free software: you can redistribute it and/or modify it under the terms
	6	! of the GNU General Public License as published by the Free Software Foundation,
	7	! either version 3 of the License, or (at your option) any later version.
	8	!
	9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
	10	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
	11	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
	12	!
	13	! You should have received a copy of the GNU General Public License along with
	14	! PALM. If not, see <http://www.gnu.org/licenses/>.
	15	!
[1818]	16	! Copyright 1997-2016 Leibniz Universitaet Hannover
[1036]	17	!--------------------------------------------------------------------------------!
	18	!
[849]	19	! Current revisions:
	20	! ------------------
[1930]	21	!
[1937]	22	!
[1930]	23	! Former revisions:
	24	! -----------------
	25	! $Id: lpm_exchange_horiz.f90 1937 2016-06-13 13:40:04Z gronemeier $
	26	!
[1937]	27	! 1936 2016-06-13 13:37:44Z suehring
	28	! Deallocation of unused memory
	29	!
[1930]	30	! 1929 2016-06-09 16:25:25Z suehring
[1929]	31	! Bugfixes:
	32	! - reallocation of new particles
	33	! ( did not work for small number of min_nr_particle )
	34	! - dynamical reallocation of north-south exchange arrays ( particles got lost )
	35	! - north-south exchange ( nr_move_north, nr_move_south were overwritten by zero )
	36	! - horizontal particle boundary conditions in serial mode
	37	!
	38	! Remove unused variables
	39	! Descriptions in variable declaration blocks added
[1321]	40	!
[1874]	41	! 1873 2016-04-18 14:50:06Z maronga
	42	! Module renamed (removed _mod)
	43	!
	44	!
[1851]	45	! 1850 2016-04-08 13:29:27Z maronga
	46	! Module renamed
	47	!
	48	!
[1823]	49	! 1822 2016-04-07 07:49:42Z hoffmann
	50	! Tails removed. Unused variables removed.
	51	!
[1784]	52	! 1783 2016-03-06 18:36:17Z raasch
	53	! new netcdf-module included
	54	!
[1692]	55	! 1691 2015-10-26 16:17:44Z maronga
	56	! Formatting corrections.
	57	!
[1686]	58	! 1685 2015-10-08 07:32:13Z raasch
	59	! bugfix concerning vertical index offset in case of ocean
	60	!
[1683]	61	! 1682 2015-10-07 23:56:08Z knoop
	62	! Code annotations made doxygen readable
	63	!
[1360]	64	! 1359 2014-04-11 17:15:14Z hoffmann
	65	! New particle structure integrated.
	66	! Kind definition added to all floating point numbers.
	67	!
[1329]	68	! 1327 2014-03-21 11:00:16Z raasch
	69	! -netcdf output queries
	70	!
[1321]	71	! 1320 2014-03-20 08:40:49Z raasch
[1320]	72	! ONLY-attribute added to USE-statements,
	73	! kind-parameters added to all INTEGER and REAL declaration statements,
	74	! kinds are defined in new module kinds,
	75	! comment fields (!:) to be used for variable explanations added to
	76	! all variable declaration statements
[849]	77	!
[1319]	78	! 1318 2014-03-17 13:35:16Z raasch
	79	! module interfaces removed
	80	!
[1037]	81	! 1036 2012-10-22 13:43:42Z raasch
	82	! code put under GPL (PALM 3.9)
	83	!
[852]	84	! 851 2012-03-15 14:32:58Z raasch
	85	! Bugfix: resetting of particle_mask and tail mask moved from end of this
	86	! routine to lpm
	87	!
[850]	88	! 849 2012-03-15 10:35:09Z raasch
	89	! initial revision (former part of advec_particles)
[849]	90	!
[850]	91	!
[849]	92	! Description:
	93	! ------------
[1822]	94	! Exchange of particles between the subdomains.
[849]	95	!------------------------------------------------------------------------------!
[1682]	96	MODULE lpm_exchange_horiz_mod
	97
[849]	98
[1320]	99	USE control_parameters, &
[1783]	100	ONLY: dz, message_string, simulated_time
[1320]	101
	102	USE cpulog, &
	103	ONLY: cpu_log, log_point_s
	104
	105	USE grid_variables, &
	106	ONLY: ddx, ddy, dx, dy
	107
	108	USE indices, &
[1359]	109	ONLY: nx, nxl, nxr, ny, nyn, nys, nzb, nzt
[1320]	110
	111	USE kinds
	112
[1359]	113	USE lpm_pack_arrays_mod, &
	114	ONLY: lpm_pack_arrays
	115
[1783]	116	USE netcdf_interface, &
	117	ONLY: netcdf_data_format
	118
[1320]	119	USE particle_attributes, &
[1822]	120	ONLY: alloc_factor, deleted_particles, grid_particles, &
	121	ibc_par_lr, ibc_par_ns, min_nr_particle, &
	122	mpi_particle_type, number_of_particles, &
	123	offset_ocean_nzt, offset_ocean_nzt_m1, particles, &
	124	particle_type, prt_count, trlp_count_sum, &
[1359]	125	trlp_count_recv_sum, trnp_count_sum, trnp_count_recv_sum, &
	126	trrp_count_sum, trrp_count_recv_sum, trsp_count_sum, &
[1822]	127	trsp_count_recv_sum, zero_particle
[1320]	128
[849]	129	USE pegrid
	130
	131	IMPLICIT NONE
	132
[1682]	133	INTEGER(iwp), PARAMETER :: NR_2_direction_move = 10000 !<
	134	INTEGER(iwp) :: nr_move_north !<
	135	INTEGER(iwp) :: nr_move_south !<
[1359]	136
[1929]	137	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: move_also_north
	138	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: move_also_south
[1359]	139
	140	SAVE
	141
	142	PRIVATE
[1936]	143	PUBLIC lpm_exchange_horiz, lpm_move_particle, realloc_particles_array, &
	144	dealloc_particles_array
[1359]	145
	146	INTERFACE lpm_exchange_horiz
	147	MODULE PROCEDURE lpm_exchange_horiz
	148	END INTERFACE lpm_exchange_horiz
	149
	150	INTERFACE lpm_move_particle
	151	MODULE PROCEDURE lpm_move_particle
	152	END INTERFACE lpm_move_particle
	153
	154	INTERFACE realloc_particles_array
	155	MODULE PROCEDURE realloc_particles_array
	156	END INTERFACE realloc_particles_array
	157
[1936]	158	INTERFACE dealloc_particles_array
	159	MODULE PROCEDURE dealloc_particles_array
	160	END INTERFACE dealloc_particles_array
[1359]	161	CONTAINS
	162
[1682]	163	!------------------------------------------------------------------------------!
	164	! Description:
	165	! ------------
	166	!> Exchange between subdomains.
	167	!> As soon as one particle has moved beyond the boundary of the domain, it
	168	!> is included in the relevant transfer arrays and marked for subsequent
	169	!> deletion on this PE.
	170	!> First sweep for crossings in x direction. Find out first the number of
	171	!> particles to be transferred and allocate temporary arrays needed to store
	172	!> them.
	173	!> For a one-dimensional decomposition along y, no transfer is necessary,
	174	!> because the particle remains on the PE, but the particle coordinate has to
	175	!> be adjusted.
	176	!------------------------------------------------------------------------------!
[1359]	177	SUBROUTINE lpm_exchange_horiz
	178
	179	IMPLICIT NONE
	180
[1929]	181	INTEGER(iwp) :: i !< grid index (x) of particle positition
	182	INTEGER(iwp) :: ip !< index variable along x
	183	INTEGER(iwp) :: j !< grid index (y) of particle positition
	184	INTEGER(iwp) :: jp !< index variable along y
	185	INTEGER(iwp) :: kp !< index variable along z
	186	INTEGER(iwp) :: n !< particle index variable
	187	INTEGER(iwp) :: trlp_count !< number of particles send to left PE
	188	INTEGER(iwp) :: trlp_count_recv !< number of particles receive from right PE
	189	INTEGER(iwp) :: trnp_count !< number of particles send to north PE
	190	INTEGER(iwp) :: trnp_count_recv !< number of particles receive from south PE
	191	INTEGER(iwp) :: trrp_count !< number of particles send to right PE
	192	INTEGER(iwp) :: trrp_count_recv !< number of particles receive from left PE
	193	INTEGER(iwp) :: trsp_count !< number of particles send to south PE
	194	INTEGER(iwp) :: trsp_count_recv !< number of particles receive from north PE
[849]	195
[1929]	196	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: rvlp !< particles received from right PE
	197	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: rvnp !< particles received from south PE
	198	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: rvrp !< particles received from left PE
	199	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: rvsp !< particles received from north PE
	200	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: trlp !< particles send to left PE
	201	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: trnp !< particles send to north PE
	202	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: trrp !< particles send to right PE
	203	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: trsp !< particles send to south PE
[849]	204
[1359]	205	CALL cpu_log( log_point_s(23), 'lpm_exchange_horiz', 'start' )
	206
[849]	207	#if defined( __parallel )
	208
[1822]	209	!
	210	!-- Exchange between subdomains.
	211	!-- As soon as one particle has moved beyond the boundary of the domain, it
	212	!-- is included in the relevant transfer arrays and marked for subsequent
	213	!-- deletion on this PE.
	214	!-- First sweep for crossings in x direction. Find out first the number of
	215	!-- particles to be transferred and allocate temporary arrays needed to store
	216	!-- them.
	217	!-- For a one-dimensional decomposition along y, no transfer is necessary,
	218	!-- because the particle remains on the PE, but the particle coordinate has to
	219	!-- be adjusted.
[849]	220	trlp_count = 0
	221	trrp_count = 0
	222
	223	trlp_count_recv = 0
	224	trrp_count_recv = 0
	225
	226	IF ( pdims(1) /= 1 ) THEN
	227	!
[1359]	228	!-- First calculate the storage necessary for sending and receiving the data.
	229	!-- Compute only first (nxl) and last (nxr) loop iterration.
	230	DO ip = nxl, nxr, nxr - nxl
	231	DO jp = nys, nyn
	232	DO kp = nzb+1, nzt
[849]	233
[1359]	234	number_of_particles = prt_count(kp,jp,ip)
	235	IF ( number_of_particles <= 0 ) CYCLE
	236	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	237	DO n = 1, number_of_particles
	238	IF ( particles(n)%particle_mask ) THEN
	239	i = ( particles(n)%x + 0.5_wp * dx ) * ddx
[1929]	240	!
	241	!-- Above calculation does not work for indices less than zero
[1359]	242	IF ( particles(n)%x < -0.5_wp * dx ) i = -1
	243
	244	IF ( i < nxl ) THEN
	245	trlp_count = trlp_count + 1
	246	ELSEIF ( i > nxr ) THEN
	247	trrp_count = trrp_count + 1
	248	ENDIF
	249	ENDIF
	250	ENDDO
	251
	252	ENDDO
	253	ENDDO
[849]	254	ENDDO
	255
	256	IF ( trlp_count == 0 ) trlp_count = 1
	257	IF ( trrp_count == 0 ) trrp_count = 1
	258
	259	ALLOCATE( trlp(trlp_count), trrp(trrp_count) )
	260
[1359]	261	trlp = zero_particle
	262	trrp = zero_particle
[849]	263
	264	trlp_count = 0
	265	trrp_count = 0
	266
	267	ENDIF
[1359]	268	!
	269	!-- Compute only first (nxl) and last (nxr) loop iterration
[1929]	270	DO ip = nxl, nxr, nxr-nxl
[1359]	271	DO jp = nys, nyn
	272	DO kp = nzb+1, nzt
	273	number_of_particles = prt_count(kp,jp,ip)
	274	IF ( number_of_particles <= 0 ) CYCLE
	275	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	276	DO n = 1, number_of_particles
	277	!
	278	!-- Only those particles that have not been marked as 'deleted' may
	279	!-- be moved.
	280	IF ( particles(n)%particle_mask ) THEN
[849]	281
[1359]	282	i = ( particles(n)%x + 0.5_wp * dx ) * ddx
[1929]	283	!
	284	!-- Above calculation does not work for indices less than zero
[1359]	285	IF ( particles(n)%x < - 0.5_wp * dx ) i = -1
[849]	286
[1359]	287	IF ( i < nxl ) THEN
	288	IF ( i < 0 ) THEN
[1929]	289	!
	290	!-- Apply boundary condition along x
[1359]	291	IF ( ibc_par_lr == 0 ) THEN
[1929]	292	!
	293	!-- Cyclic condition
[1359]	294	IF ( pdims(1) == 1 ) THEN
	295	particles(n)%x = ( nx + 1 ) * dx + particles(n)%x
	296	particles(n)%origin_x = ( nx + 1 ) * dx + &
	297	particles(n)%origin_x
	298	ELSE
	299	trlp_count = trlp_count + 1
	300	trlp(trlp_count) = particles(n)
	301	trlp(trlp_count)%x = ( nx + 1 ) * dx + trlp(trlp_count)%x
	302	trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x + &
	303	( nx + 1 ) * dx
	304	particles(n)%particle_mask = .FALSE.
	305	deleted_particles = deleted_particles + 1
[849]	306
[1359]	307	IF ( trlp(trlp_count)%x >= (nx + 0.5_wp)* dx - 1.0E-12_wp ) THEN
	308	trlp(trlp_count)%x = trlp(trlp_count)%x - 1.0E-10_wp
	309	!++ why is 1 subtracted in next statement???
	310	trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x - 1
	311	ENDIF
[849]	312
[1359]	313	ENDIF
[849]	314
[1359]	315	ELSEIF ( ibc_par_lr == 1 ) THEN
[1929]	316	!
	317	!-- Particle absorption
[1359]	318	particles(n)%particle_mask = .FALSE.
	319	deleted_particles = deleted_particles + 1
[849]	320
[1359]	321	ELSEIF ( ibc_par_lr == 2 ) THEN
[1929]	322	!
	323	!-- Particle reflection
[1359]	324	particles(n)%x = -particles(n)%x
	325	particles(n)%speed_x = -particles(n)%speed_x
[849]	326
[1359]	327	ENDIF
	328	ELSE
[1929]	329	!
	330	!-- Store particle data in the transfer array, which will be
	331	!-- send to the neighbouring PE
[1359]	332	trlp_count = trlp_count + 1
	333	trlp(trlp_count) = particles(n)
	334	particles(n)%particle_mask = .FALSE.
	335	deleted_particles = deleted_particles + 1
[849]	336
[1359]	337	ENDIF
[849]	338
[1359]	339	ELSEIF ( i > nxr ) THEN
	340	IF ( i > nx ) THEN
[1929]	341	!
	342	!-- Apply boundary condition along x
[1359]	343	IF ( ibc_par_lr == 0 ) THEN
[1929]	344	!
	345	!-- Cyclic condition
[1359]	346	IF ( pdims(1) == 1 ) THEN
	347	particles(n)%x = particles(n)%x - ( nx + 1 ) * dx
	348	particles(n)%origin_x = particles(n)%origin_x - &
	349	( nx + 1 ) * dx
	350	ELSE
	351	trrp_count = trrp_count + 1
	352	trrp(trrp_count) = particles(n)
	353	trrp(trrp_count)%x = trrp(trrp_count)%x - ( nx + 1 ) * dx
	354	trrp(trrp_count)%origin_x = trrp(trrp_count)%origin_x - &
	355	( nx + 1 ) * dx
	356	particles(n)%particle_mask = .FALSE.
	357	deleted_particles = deleted_particles + 1
[849]	358
[1359]	359	ENDIF
[849]	360
[1359]	361	ELSEIF ( ibc_par_lr == 1 ) THEN
[1929]	362	!
	363	!-- Particle absorption
[1359]	364	particles(n)%particle_mask = .FALSE.
	365	deleted_particles = deleted_particles + 1
[849]	366
[1359]	367	ELSEIF ( ibc_par_lr == 2 ) THEN
[1929]	368	!
	369	!-- Particle reflection
[1359]	370	particles(n)%x = 2 * ( nx * dx ) - particles(n)%x
	371	particles(n)%speed_x = -particles(n)%speed_x
[849]	372
[1359]	373	ENDIF
	374	ELSE
[1929]	375	!
	376	!-- Store particle data in the transfer array, which will be send
	377	!-- to the neighbouring PE
[1359]	378	trrp_count = trrp_count + 1
	379	trrp(trrp_count) = particles(n)
	380	particles(n)%particle_mask = .FALSE.
	381	deleted_particles = deleted_particles + 1
[849]	382
[1359]	383	ENDIF
[849]	384
[1359]	385	ENDIF
	386	ENDIF
[1929]	387
[1359]	388	ENDDO
	389	ENDDO
	390	ENDDO
[849]	391	ENDDO
	392
	393	!
[1929]	394	!-- Allocate arrays required for north-south exchange, as these
	395	!-- are used directly after particles are exchange along x-direction.
	396	ALLOCATE( move_also_north(1:NR_2_direction_move) )
	397	ALLOCATE( move_also_south(1:NR_2_direction_move) )
	398
	399	nr_move_north = 0
	400	nr_move_south = 0
	401	!
[849]	402	!-- Send left boundary, receive right boundary (but first exchange how many
	403	!-- and check, if particle storage must be extended)
	404	IF ( pdims(1) /= 1 ) THEN
	405
	406	CALL MPI_SENDRECV( trlp_count, 1, MPI_INTEGER, pleft, 0, &
	407	trrp_count_recv, 1, MPI_INTEGER, pright, 0, &
	408	comm2d, status, ierr )
	409
[1359]	410	ALLOCATE(rvrp(MAX(1,trrp_count_recv)))
[849]	411
[1359]	412	CALL MPI_SENDRECV( trlp(1)%radius, max(1,trlp_count), mpi_particle_type,&
	413	pleft, 1, rvrp(1)%radius, &
	414	max(1,trrp_count_recv), mpi_particle_type, pright, 1,&
[849]	415	comm2d, status, ierr )
	416
[1359]	417	IF ( trrp_count_recv > 0 ) CALL Add_particles_to_gridcell(rvrp(1:trrp_count_recv))
	418
	419	DEALLOCATE(rvrp)
	420
[849]	421	!
	422	!-- Send right boundary, receive left boundary
	423	CALL MPI_SENDRECV( trrp_count, 1, MPI_INTEGER, pright, 0, &
	424	trlp_count_recv, 1, MPI_INTEGER, pleft, 0, &
	425	comm2d, status, ierr )
	426
[1359]	427	ALLOCATE(rvlp(MAX(1,trlp_count_recv)))
[849]	428
[1359]	429	CALL MPI_SENDRECV( trrp(1)%radius, max(1,trrp_count), mpi_particle_type,&
	430	pright, 1, rvlp(1)%radius, &
	431	max(1,trlp_count_recv), mpi_particle_type, pleft, 1, &
[849]	432	comm2d, status, ierr )
	433
[1359]	434	IF ( trlp_count_recv > 0 ) CALL Add_particles_to_gridcell(rvlp(1:trlp_count_recv))
	435
[1822]	436	DEALLOCATE( rvlp )
[1359]	437	DEALLOCATE( trlp, trrp )
[849]	438
	439	ENDIF
	440
	441	!
	442	!-- Check whether particles have crossed the boundaries in y direction. Note
	443	!-- that this case can also apply to particles that have just been received
	444	!-- from the adjacent right or left PE.
	445	!-- Find out first the number of particles to be transferred and allocate
	446	!-- temporary arrays needed to store them.
[1929]	447	!-- For a one-dimensional decomposition along y, no transfer is necessary,
[849]	448	!-- because the particle remains on the PE.
[1359]	449	trsp_count = nr_move_south
	450	trnp_count = nr_move_north
[849]	451
	452	trsp_count_recv = 0
	453	trnp_count_recv = 0
	454
	455	IF ( pdims(2) /= 1 ) THEN
	456	!
	457	!-- First calculate the storage necessary for sending and receiving the
	458	!-- data
[1359]	459	DO ip = nxl, nxr
[1929]	460	DO jp = nys, nyn, nyn-nys !compute only first (nys) and last (nyn) loop iterration
[1359]	461	DO kp = nzb+1, nzt
	462	number_of_particles = prt_count(kp,jp,ip)
	463	IF ( number_of_particles <= 0 ) CYCLE
	464	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	465	DO n = 1, number_of_particles
	466	IF ( particles(n)%particle_mask ) THEN
	467	j = ( particles(n)%y + 0.5_wp * dy ) * ddy
[849]	468	!
[1359]	469	!-- Above calculation does not work for indices less than zero
	470	IF ( particles(n)%y < -0.5_wp * dy ) j = -1
[849]	471
[1359]	472	IF ( j < nys ) THEN
	473	trsp_count = trsp_count + 1
	474	ELSEIF ( j > nyn ) THEN
	475	trnp_count = trnp_count + 1
	476	ENDIF
	477	ENDIF
	478	ENDDO
	479	ENDDO
	480	ENDDO
[849]	481	ENDDO
	482
	483	IF ( trsp_count == 0 ) trsp_count = 1
	484	IF ( trnp_count == 0 ) trnp_count = 1
	485
	486	ALLOCATE( trsp(trsp_count), trnp(trnp_count) )
	487
[1359]	488	trsp = zero_particle
	489	trnp = zero_particle
[849]	490
[1359]	491	trsp_count = nr_move_south
	492	trnp_count = nr_move_north
[1929]	493
[1359]	494	trsp(1:nr_move_south) = move_also_south(1:nr_move_south)
	495	trnp(1:nr_move_north) = move_also_north(1:nr_move_north)
	496
[849]	497	ENDIF
	498
[1359]	499	DO ip = nxl, nxr
	500	DO jp = nys, nyn, nyn-nys ! compute only first (nys) and last (nyn) loop iterration
	501	DO kp = nzb+1, nzt
	502	number_of_particles = prt_count(kp,jp,ip)
	503	IF ( number_of_particles <= 0 ) CYCLE
	504	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	505	DO n = 1, number_of_particles
[849]	506	!
[1359]	507	!-- Only those particles that have not been marked as 'deleted' may
	508	!-- be moved.
	509	IF ( particles(n)%particle_mask ) THEN
[1929]	510
[1359]	511	j = ( particles(n)%y + 0.5_wp * dy ) * ddy
[849]	512	!
[1359]	513	!-- Above calculation does not work for indices less than zero
	514	IF ( particles(n)%y < -0.5_wp * dy ) j = -1
[849]	515
[1359]	516	IF ( j < nys ) THEN
	517	IF ( j < 0 ) THEN
[849]	518	!
[1359]	519	!-- Apply boundary condition along y
	520	IF ( ibc_par_ns == 0 ) THEN
[849]	521	!
[1359]	522	!-- Cyclic condition
	523	IF ( pdims(2) == 1 ) THEN
	524	particles(n)%y = ( ny + 1 ) * dy + particles(n)%y
	525	particles(n)%origin_y = ( ny + 1 ) * dy + &
[1929]	526	particles(n)%origin_y
[1359]	527	ELSE
[1929]	528	trsp_count = trsp_count + 1
	529	trsp(trsp_count) = particles(n)
[1359]	530	trsp(trsp_count)%y = ( ny + 1 ) * dy + &
[1929]	531	trsp(trsp_count)%y
[1359]	532	trsp(trsp_count)%origin_y = trsp(trsp_count)%origin_y &
[1929]	533	+ ( ny + 1 ) * dy
[1359]	534	particles(n)%particle_mask = .FALSE.
	535	deleted_particles = deleted_particles + 1
[849]	536
[1359]	537	IF ( trsp(trsp_count)%y >= (ny+0.5_wp)* dy - 1.0E-12_wp ) THEN
	538	trsp(trsp_count)%y = trsp(trsp_count)%y - 1.0E-10_wp
	539	!++ why is 1 subtracted in next statement???
	540	trsp(trsp_count)%origin_y = &
[1929]	541	trsp(trsp_count)%origin_y - 1
[1359]	542	ENDIF
[849]	543
[1359]	544	ENDIF
[849]	545
[1359]	546	ELSEIF ( ibc_par_ns == 1 ) THEN
[849]	547	!
[1359]	548	!-- Particle absorption
	549	particles(n)%particle_mask = .FALSE.
[1929]	550	deleted_particles = deleted_particles + 1
[849]	551
[1359]	552	ELSEIF ( ibc_par_ns == 2 ) THEN
[849]	553	!
[1359]	554	!-- Particle reflection
	555	particles(n)%y = -particles(n)%y
	556	particles(n)%speed_y = -particles(n)%speed_y
[849]	557
[1359]	558	ENDIF
	559	ELSE
[849]	560	!
[1359]	561	!-- Store particle data in the transfer array, which will
	562	!-- be send to the neighbouring PE
	563	trsp_count = trsp_count + 1
	564	trsp(trsp_count) = particles(n)
	565	particles(n)%particle_mask = .FALSE.
	566	deleted_particles = deleted_particles + 1
[849]	567
[1359]	568	ENDIF
[849]	569
[1359]	570	ELSEIF ( j > nyn ) THEN
	571	IF ( j > ny ) THEN
[849]	572	!
[1929]	573	!-- Apply boundary condition along y
[1359]	574	IF ( ibc_par_ns == 0 ) THEN
[849]	575	!
[1359]	576	!-- Cyclic condition
	577	IF ( pdims(2) == 1 ) THEN
[1929]	578	particles(n)%y = particles(n)%y - ( ny + 1 ) * dy
[1359]	579	particles(n)%origin_y = &
[1929]	580	particles(n)%origin_y - ( ny + 1 ) * dy
[1359]	581	ELSE
[1929]	582	trnp_count = trnp_count + 1
	583	trnp(trnp_count) = particles(n)
[1359]	584	trnp(trnp_count)%y = &
[1929]	585	trnp(trnp_count)%y - ( ny + 1 ) * dy
[1359]	586	trnp(trnp_count)%origin_y = &
[1929]	587	trnp(trnp_count)%origin_y - ( ny + 1 ) * dy
[1359]	588	particles(n)%particle_mask = .FALSE.
[1929]	589	deleted_particles = deleted_particles + 1
[1359]	590	ENDIF
[849]	591
[1359]	592	ELSEIF ( ibc_par_ns == 1 ) THEN
[849]	593	!
[1359]	594	!-- Particle absorption
	595	particles(n)%particle_mask = .FALSE.
	596	deleted_particles = deleted_particles + 1
[849]	597
[1359]	598	ELSEIF ( ibc_par_ns == 2 ) THEN
[849]	599	!
[1359]	600	!-- Particle reflection
	601	particles(n)%y = 2 * ( ny * dy ) - particles(n)%y
	602	particles(n)%speed_y = -particles(n)%speed_y
[849]	603
[1359]	604	ENDIF
	605	ELSE
[849]	606	!
[1359]	607	!-- Store particle data in the transfer array, which will
	608	!-- be send to the neighbouring PE
	609	trnp_count = trnp_count + 1
	610	trnp(trnp_count) = particles(n)
	611	particles(n)%particle_mask = .FALSE.
	612	deleted_particles = deleted_particles + 1
[849]	613
[1359]	614	ENDIF
	615
	616	ENDIF
[849]	617	ENDIF
[1359]	618	ENDDO
	619	ENDDO
	620	ENDDO
[849]	621	ENDDO
	622
	623	!
	624	!-- Send front boundary, receive back boundary (but first exchange how many
	625	!-- and check, if particle storage must be extended)
	626	IF ( pdims(2) /= 1 ) THEN
	627
	628	CALL MPI_SENDRECV( trsp_count, 1, MPI_INTEGER, psouth, 0, &
	629	trnp_count_recv, 1, MPI_INTEGER, pnorth, 0, &
	630	comm2d, status, ierr )
	631
[1359]	632	ALLOCATE(rvnp(MAX(1,trnp_count_recv)))
[1929]	633
[1359]	634	CALL MPI_SENDRECV( trsp(1)%radius, trsp_count, mpi_particle_type, &
	635	psouth, 1, rvnp(1)%radius, &
[849]	636	trnp_count_recv, mpi_particle_type, pnorth, 1, &
	637	comm2d, status, ierr )
	638
[1359]	639	IF ( trnp_count_recv > 0 ) CALL Add_particles_to_gridcell(rvnp(1:trnp_count_recv))
	640
	641	DEALLOCATE(rvnp)
	642
[849]	643	!
	644	!-- Send back boundary, receive front boundary
	645	CALL MPI_SENDRECV( trnp_count, 1, MPI_INTEGER, pnorth, 0, &
	646	trsp_count_recv, 1, MPI_INTEGER, psouth, 0, &
	647	comm2d, status, ierr )
	648
[1359]	649	ALLOCATE(rvsp(MAX(1,trsp_count_recv)))
[849]	650
[1359]	651	CALL MPI_SENDRECV( trnp(1)%radius, trnp_count, mpi_particle_type, &
	652	pnorth, 1, rvsp(1)%radius, &
[849]	653	trsp_count_recv, mpi_particle_type, psouth, 1, &
	654	comm2d, status, ierr )
	655
[1359]	656	IF ( trsp_count_recv > 0 ) CALL Add_particles_to_gridcell(rvsp(1:trsp_count_recv))
	657
	658	DEALLOCATE(rvsp)
	659
[849]	660	number_of_particles = number_of_particles + trsp_count_recv
	661
	662	DEALLOCATE( trsp, trnp )
	663
	664	ENDIF
	665
[1929]	666	DEALLOCATE( move_also_north )
	667	DEALLOCATE( move_also_south )
	668
[849]	669	#else
	670
[1929]	671	DO ip = nxl, nxr, nxr-nxl
	672	DO jp = nys, nyn
	673	DO kp = nzb+1, nzt
	674	number_of_particles = prt_count(kp,jp,ip)
	675	IF ( number_of_particles <= 0 ) CYCLE
	676	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	677	DO n = 1, number_of_particles
[849]	678	!
[1929]	679	!-- Apply boundary conditions
[849]	680
[1929]	681	IF ( particles(n)%x < -0.5_wp * dx ) THEN
[849]	682
[1929]	683	IF ( ibc_par_lr == 0 ) THEN
[849]	684	!
[1929]	685	!-- Cyclic boundary. Relevant coordinate has to be changed.
	686	particles(n)%x = ( nx + 1 ) * dx + particles(n)%x
[1822]	687
[1929]	688	ELSEIF ( ibc_par_lr == 1 ) THEN
[849]	689	!
[1929]	690	!-- Particle absorption
	691	particles(n)%particle_mask = .FALSE.
	692	deleted_particles = deleted_particles + 1
[1822]	693
[1929]	694	ELSEIF ( ibc_par_lr == 2 ) THEN
[849]	695	!
[1929]	696	!-- Particle reflection
	697	particles(n)%x = -dx - particles(n)%x
	698	particles(n)%speed_x = -particles(n)%speed_x
	699	ENDIF
[849]	700
[1929]	701	ELSEIF ( particles(n)%x >= ( nx + 0.5_wp ) * dx ) THEN
[849]	702
[1929]	703	IF ( ibc_par_lr == 0 ) THEN
[849]	704	!
[1929]	705	!-- Cyclic boundary. Relevant coordinate has to be changed.
	706	particles(n)%x = particles(n)%x - ( nx + 1 ) * dx
[1822]	707
[1929]	708	ELSEIF ( ibc_par_lr == 1 ) THEN
[849]	709	!
[1929]	710	!-- Particle absorption
	711	particles(n)%particle_mask = .FALSE.
	712	deleted_particles = deleted_particles + 1
[1822]	713
[1929]	714	ELSEIF ( ibc_par_lr == 2 ) THEN
[849]	715	!
[1929]	716	!-- Particle reflection
	717	particles(n)%x = ( nx + 1 ) * dx - particles(n)%x
	718	particles(n)%speed_x = -particles(n)%speed_x
	719	ENDIF
[849]	720
[1929]	721	ENDIF
	722	ENDDO
	723	ENDDO
	724	ENDDO
	725	ENDDO
[849]	726
[1929]	727	DO ip = nxl, nxr
	728	DO jp = nys, nyn, nyn-nys
	729	DO kp = nzb+1, nzt
	730	number_of_particles = prt_count(kp,jp,ip)
	731	IF ( number_of_particles <= 0 ) CYCLE
	732	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	733	DO n = 1, number_of_particles
[849]	734
[1929]	735	IF ( particles(n)%y < -0.5_wp * dy ) THEN
	736
	737	IF ( ibc_par_ns == 0 ) THEN
[849]	738	!
[1929]	739	!-- Cyclic boundary. Relevant coordinate has to be changed.
	740	particles(n)%y = ( ny + 1 ) * dy + particles(n)%y
[1822]	741
[1929]	742	ELSEIF ( ibc_par_ns == 1 ) THEN
[849]	743	!
[1929]	744	!-- Particle absorption
	745	particles(n)%particle_mask = .FALSE.
	746	deleted_particles = deleted_particles + 1
[1822]	747
[1929]	748	ELSEIF ( ibc_par_ns == 2 ) THEN
[849]	749	!
[1929]	750	!-- Particle reflection
	751	particles(n)%y = -dy - particles(n)%y
	752	particles(n)%speed_y = -particles(n)%speed_y
	753	ENDIF
[849]	754
[1929]	755	ELSEIF ( particles(n)%y >= ( ny + 0.5_wp ) * dy ) THEN
[849]	756
[1929]	757	IF ( ibc_par_ns == 0 ) THEN
[849]	758	!
[1929]	759	!-- Cyclic boundary. Relevant coordinate has to be changed.
	760	particles(n)%y = particles(n)%y - ( ny + 1 ) * dy
[1822]	761
[1929]	762	ELSEIF ( ibc_par_ns == 1 ) THEN
[849]	763	!
[1929]	764	!-- Particle absorption
	765	particles(n)%particle_mask = .FALSE.
	766	deleted_particles = deleted_particles + 1
[1822]	767
[1929]	768	ELSEIF ( ibc_par_ns == 2 ) THEN
[849]	769	!
[1929]	770	!-- Particle reflection
	771	particles(n)%y = ( ny + 1 ) * dy - particles(n)%y
	772	particles(n)%speed_y = -particles(n)%speed_y
	773	ENDIF
[849]	774
[1929]	775	ENDIF
	776
	777	ENDDO
	778	ENDDO
	779	ENDDO
[849]	780	ENDDO
	781	#endif
	782
	783	!
	784	!-- Accumulate the number of particles transferred between the subdomains
	785	#if defined( __parallel )
	786	trlp_count_sum = trlp_count_sum + trlp_count
	787	trlp_count_recv_sum = trlp_count_recv_sum + trlp_count_recv
	788	trrp_count_sum = trrp_count_sum + trrp_count
	789	trrp_count_recv_sum = trrp_count_recv_sum + trrp_count_recv
	790	trsp_count_sum = trsp_count_sum + trsp_count
	791	trsp_count_recv_sum = trsp_count_recv_sum + trsp_count_recv
	792	trnp_count_sum = trnp_count_sum + trnp_count
	793	trnp_count_recv_sum = trnp_count_recv_sum + trnp_count_recv
	794	#endif
	795
[1359]	796	CALL cpu_log( log_point_s(23), 'lpm_exchange_horiz', 'stop' )
[849]	797
	798	END SUBROUTINE lpm_exchange_horiz
[1359]	799
[1682]	800	!------------------------------------------------------------------------------!
	801	! Description:
	802	! ------------
	803	!> If a particle moves from one processor to another, this subroutine moves
	804	!> the corresponding elements from the particle arrays of the old grid cells
	805	!> to the particle arrays of the new grid cells.
	806	!------------------------------------------------------------------------------!
[1359]	807	SUBROUTINE Add_particles_to_gridcell (particle_array)
	808
	809	IMPLICIT NONE
	810
[1929]	811	INTEGER(iwp) :: ip !< grid index (x) of particle
	812	INTEGER(iwp) :: jp !< grid index (x) of particle
	813	INTEGER(iwp) :: kp !< grid index (x) of particle
	814	INTEGER(iwp) :: n !< index variable of particle
	815	INTEGER(iwp) :: pindex !< dummy argument for new number of particles per grid box
[1359]	816
[1682]	817	LOGICAL :: pack_done !<
[1359]	818
[1929]	819	TYPE(particle_type), DIMENSION(:), INTENT(IN) :: particle_array !< new particles in a grid box
	820	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: temp_ns !< temporary particle array for reallocation
[1359]	821
	822	pack_done = .FALSE.
	823
[1929]	824	DO n = 1, SIZE(particle_array)
[1359]	825
[1929]	826	IF ( .NOT. particle_array(n)%particle_mask ) CYCLE
	827
[1359]	828	ip = ( particle_array(n)%x + 0.5_wp * dx ) * ddx
	829	jp = ( particle_array(n)%y + 0.5_wp * dy ) * ddy
[1929]	830	kp = particle_array(n)%z / dz + 1 + offset_ocean_nzt
[1359]	831
	832	IF ( ip >= nxl .AND. ip <= nxr .AND. jp >= nys .AND. jp <= nyn &
	833	.AND. kp >= nzb+1 .AND. kp <= nzt) THEN ! particle stays on processor
	834	number_of_particles = prt_count(kp,jp,ip)
	835	particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
	836
	837	pindex = prt_count(kp,jp,ip)+1
	838	IF( pindex > SIZE(grid_particles(kp,jp,ip)%particles) ) THEN
	839	IF ( pack_done ) THEN
	840	CALL realloc_particles_array (ip,jp,kp)
	841	ELSE
	842	CALL lpm_pack_arrays
	843	prt_count(kp,jp,ip) = number_of_particles
	844	pindex = prt_count(kp,jp,ip)+1
	845	IF ( pindex > SIZE(grid_particles(kp,jp,ip)%particles) ) THEN
	846	CALL realloc_particles_array (ip,jp,kp)
	847	ENDIF
	848	pack_done = .TRUE.
	849	ENDIF
	850	ENDIF
	851	grid_particles(kp,jp,ip)%particles(pindex) = particle_array(n)
	852	prt_count(kp,jp,ip) = pindex
	853	ELSE
[1929]	854	IF ( jp <= nys - 1 ) THEN
[1359]	855	nr_move_south = nr_move_south+1
[1929]	856	!
	857	!-- Before particle information is swapped to exchange-array, check
	858	!-- if enough memory is allocated. If required, reallocate exchange
	859	!-- array.
	860	IF ( nr_move_south > SIZE(move_also_south) ) THEN
	861	!
	862	!-- At first, allocate further temporary array to swap particle
	863	!-- information.
	864	ALLOCATE( temp_ns(SIZE(move_also_south)+NR_2_direction_move) )
	865	temp_ns(1:nr_move_south-1) = move_also_south(1:nr_move_south-1)
	866	DEALLOCATE( move_also_south )
	867	ALLOCATE( move_also_south(SIZE(temp_ns)) )
	868	move_also_south(1:nr_move_south-1) = temp_ns(1:nr_move_south-1)
	869	DEALLOCATE( temp_ns )
	870
	871	ENDIF
	872
[1359]	873	move_also_south(nr_move_south) = particle_array(n)
[1929]	874
[1359]	875	IF ( jp == -1 ) THEN
	876	move_also_south(nr_move_south)%y = &
	877	move_also_south(nr_move_south)%y + ( ny + 1 ) * dy
	878	move_also_south(nr_move_south)%origin_y = &
	879	move_also_south(nr_move_south)%origin_y + ( ny + 1 ) * dy
	880	ENDIF
[1929]	881	ELSEIF ( jp >= nyn+1 ) THEN
[1359]	882	nr_move_north = nr_move_north+1
[1929]	883	!
	884	!-- Before particle information is swapped to exchange-array, check
	885	!-- if enough memory is allocated. If required, reallocate exchange
	886	!-- array.
	887	IF ( nr_move_north > SIZE(move_also_north) ) THEN
	888	!
	889	!-- At first, allocate further temporary array to swap particle
	890	!-- information.
	891	ALLOCATE( temp_ns(SIZE(move_also_north)+NR_2_direction_move) )
	892	temp_ns(1:nr_move_north-1) = move_also_south(1:nr_move_north-1)
	893	DEALLOCATE( move_also_north )
	894	ALLOCATE( move_also_north(SIZE(temp_ns)) )
	895	move_also_north(1:nr_move_north-1) = temp_ns(1:nr_move_north-1)
	896	DEALLOCATE( temp_ns )
	897
	898	ENDIF
	899
[1359]	900	move_also_north(nr_move_north) = particle_array(n)
	901	IF ( jp == ny+1 ) THEN
	902	move_also_north(nr_move_north)%y = &
	903	move_also_north(nr_move_north)%y - ( ny + 1 ) * dy
	904	move_also_north(nr_move_north)%origin_y = &
	905	move_also_north(nr_move_north)%origin_y - ( ny + 1 ) * dy
	906	ENDIF
	907	ELSE
	908	WRITE(0,'(a,8i7)') 'particle out of range ',myid,ip,jp,kp,nxl,nxr,nys,nyn
	909	ENDIF
	910	ENDIF
	911	ENDDO
	912
	913	RETURN
	914
	915	END SUBROUTINE Add_particles_to_gridcell
	916
	917
	918
	919
[1682]	920	!------------------------------------------------------------------------------!
	921	! Description:
	922	! ------------
	923	!> If a particle moves from one grid cell to another (on the current
	924	!> processor!), this subroutine moves the corresponding element from the
	925	!> particle array of the old grid cell to the particle array of the new grid
	926	!> cell.
	927	!------------------------------------------------------------------------------!
[1359]	928	SUBROUTINE lpm_move_particle
	929
	930	IMPLICIT NONE
	931
[1929]	932	INTEGER(iwp) :: i !< grid index (x) of particle position
	933	INTEGER(iwp) :: ip !< index variable along x
	934	INTEGER(iwp) :: j !< grid index (y) of particle position
	935	INTEGER(iwp) :: jp !< index variable along y
	936	INTEGER(iwp) :: k !< grid index (z) of particle position
	937	INTEGER(iwp) :: kp !< index variable along z
	938	INTEGER(iwp) :: n !< index variable for particle array
	939	INTEGER(iwp) :: np_old_cell !< number of particles per grid box before moving
	940	INTEGER(iwp) :: n_start !< start index
	941	INTEGER(iwp) :: pindex !< dummy argument for number of new particle per grid box
[1359]	942
[1682]	943	LOGICAL :: pack_done !<
[1359]	944
[1929]	945	TYPE(particle_type), DIMENSION(:), POINTER :: particles_old_cell !< particles before moving
[1359]	946
	947	CALL cpu_log( log_point_s(41), 'lpm_move_particle', 'start' )
	948
	949	DO ip = nxl, nxr
	950	DO jp = nys, nyn
	951	DO kp = nzb+1, nzt
	952
	953	np_old_cell = prt_count(kp,jp,ip)
	954	IF ( np_old_cell <= 0 ) CYCLE
	955	particles_old_cell => grid_particles(kp,jp,ip)%particles(1:np_old_cell)
	956	n_start = -1
	957
	958	DO n = 1, np_old_cell
	959	i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx
	960	j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy
	961	k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt
	962	!
	963	!-- Check, if particle has moved to another grid cell.
	964	IF ( i /= ip .OR. j /= jp .OR. k /= kp ) THEN
	965	!
	966	!-- The particle has moved to another grid cell. Now check, if
	967	!-- particle stays on the same processor.
	968	IF ( i >= nxl .AND. i <= nxr .AND. j >= nys .AND. &
	969	j <= nyn .AND. k >= nzb+1 .AND. k <= nzt) THEN
	970	!
	971	!-- If the particle stays on the same processor, the particle
	972	!-- will be added to the particle array of the new processor.
	973	number_of_particles = prt_count(k,j,i)
	974	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
	975
	976	pindex = prt_count(k,j,i)+1
	977	IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) &
	978	THEN
	979	n_start = n
	980	EXIT
	981	ENDIF
	982
	983	grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n)
	984	prt_count(k,j,i) = pindex
	985
	986	particles_old_cell(n)%particle_mask = .FALSE.
	987	ENDIF
	988	ENDIF
	989	ENDDO
	990
[1691]	991	IF ( n_start >= 0 ) THEN
[1359]	992	pack_done = .FALSE.
	993	DO n = n_start, np_old_cell
	994	i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx
	995	j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy
	996	k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt
	997	IF ( i /= ip .OR. j /= jp .OR. k /= kp ) THEN
	998	!
	999	!-- Particle is in different box
	1000	IF ( i >= nxl .AND. i <= nxr .AND. j >= nys .AND. &
	1001	j <= nyn .AND. k >= nzb+1 .AND. k <= nzt) THEN
	1002	!
	1003	!-- Particle stays on processor
	1004	number_of_particles = prt_count(k,j,i)
	1005	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
	1006
	1007	pindex = prt_count(k,j,i)+1
	1008	IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) &
	1009	THEN
	1010	IF ( pack_done ) THEN
	1011	CALL realloc_particles_array(i,j,k)
	1012	ELSE
	1013	CALL lpm_pack_arrays
	1014	prt_count(k,j,i) = number_of_particles
	1015	!
	1016	!-- If number of particles in the new grid box
	1017	!-- exceeds its allocated memory, the particle array
	1018	!-- will be reallocated
	1019	IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) THEN
	1020	CALL realloc_particles_array(i,j,k)
	1021	ENDIF
	1022
	1023	pack_done = .TRUE.
	1024	ENDIF
	1025	ENDIF
	1026
	1027	grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n)
	1028	prt_count(k,j,i) = pindex
	1029
	1030	particles_old_cell(n)%particle_mask = .FALSE.
	1031	ENDIF
	1032	ENDIF
	1033	ENDDO
	1034	ENDIF
	1035	ENDDO
	1036	ENDDO
	1037	ENDDO
	1038
	1039	CALL cpu_log( log_point_s(41), 'lpm_move_particle', 'stop' )
	1040
	1041	RETURN
	1042
	1043	END SUBROUTINE lpm_move_particle
	1044
[1682]	1045	!------------------------------------------------------------------------------!
	1046	! Description:
	1047	! ------------
[1929]	1048	!> If the allocated memory for the particle array do not suffice to add arriving
	1049	!> particles from neighbour grid cells, this subrouting reallocates the
	1050	!> particle array to assure enough memory is available.
[1682]	1051	!------------------------------------------------------------------------------!
[1359]	1052	SUBROUTINE realloc_particles_array (i,j,k,size_in)
	1053
	1054	IMPLICIT NONE
	1055
[1682]	1056	INTEGER(iwp), INTENT(in) :: i !<
	1057	INTEGER(iwp), INTENT(in) :: j !<
	1058	INTEGER(iwp), INTENT(in) :: k !<
[1822]	1059	INTEGER(iwp), INTENT(in), OPTIONAL :: size_in !<
[1359]	1060
[1682]	1061	INTEGER(iwp) :: old_size !<
	1062	INTEGER(iwp) :: new_size !<
	1063	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: tmp_particles_d !<
	1064	TYPE(particle_type), DIMENSION(500) :: tmp_particles_s !<
[1359]	1065
	1066	old_size = SIZE(grid_particles(k,j,i)%particles)
	1067
	1068	IF ( PRESENT(size_in) ) THEN
	1069	new_size = size_in
	1070	ELSE
[1822]	1071	new_size = old_size * ( 1.0_wp + alloc_factor / 100.0_wp )
[1359]	1072	ENDIF
	1073
[1929]	1074	new_size = MAX( new_size, min_nr_particle, old_size + 1 )
[1359]	1075
	1076	IF ( old_size <= 500 ) THEN
	1077
	1078	tmp_particles_s(1:old_size) = grid_particles(k,j,i)%particles(1:old_size)
	1079
	1080	DEALLOCATE(grid_particles(k,j,i)%particles)
	1081	ALLOCATE(grid_particles(k,j,i)%particles(new_size))
	1082
	1083	grid_particles(k,j,i)%particles(1:old_size) = tmp_particles_s(1:old_size)
	1084	grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
	1085
	1086	ELSE
	1087
	1088	ALLOCATE(tmp_particles_d(new_size))
	1089	tmp_particles_d(1:old_size) = grid_particles(k,j,i)%particles
	1090
	1091	DEALLOCATE(grid_particles(k,j,i)%particles)
	1092	ALLOCATE(grid_particles(k,j,i)%particles(new_size))
	1093
	1094	grid_particles(k,j,i)%particles(1:old_size) = tmp_particles_d(1:old_size)
	1095	grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
	1096
	1097	DEALLOCATE(tmp_particles_d)
	1098
	1099	ENDIF
	1100	particles => grid_particles(k,j,i)%particles(1:number_of_particles)
	1101
	1102	RETURN
	1103	END SUBROUTINE realloc_particles_array
	1104
[1936]	1105
	1106
	1107
	1108
	1109	SUBROUTINE dealloc_particles_array
	1110
	1111	IMPLICIT NONE
	1112
	1113	INTEGER(iwp) :: i
	1114	INTEGER(iwp) :: j
	1115	INTEGER(iwp) :: k
	1116	INTEGER(iwp) :: old_size !<
	1117	INTEGER(iwp) :: new_size !<
	1118
	1119	LOGICAL :: dealloc
	1120
	1121	TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: tmp_particles_d !<
	1122	TYPE(particle_type), DIMENSION(500) :: tmp_particles_s !<
	1123
	1124	DO i = nxl, nxr
	1125	DO j = nys, nyn
	1126	DO k = nzb+1, nzt
	1127	!
	1128	!-- Determine number of active particles
	1129	number_of_particles = prt_count(k,j,i)
	1130	!
	1131	!-- Determine allocated memory size
	1132	old_size = SIZE( grid_particles(k,j,i)%particles )
	1133	!
	1134	!-- Check for large unused memory
	1135	dealloc = ( ( number_of_particles < min_nr_particle .AND. &
	1136	old_size > min_nr_particle ) .OR. &
	1137	( number_of_particles > min_nr_particle .AND. &
	1138	old_size - number_of_particles * &
	1139	( 1.0_wp + 0.01_wp * alloc_factor ) > 0.0_wp ) )
	1140
	1141
	1142	IF ( dealloc ) THEN
	1143	IF ( number_of_particles < min_nr_particle ) THEN
	1144	new_size = min_nr_particle
	1145	ELSE
	1146	new_size = INT( number_of_particles * ( 1.0_wp + 0.01_wp * alloc_factor ) )
	1147	ENDIF
	1148
	1149	IF ( number_of_particles <= 500 ) THEN
	1150
	1151	tmp_particles_s(1:number_of_particles) = grid_particles(k,j,i)%particles(1:number_of_particles)
	1152
	1153	DEALLOCATE(grid_particles(k,j,i)%particles)
	1154	ALLOCATE(grid_particles(k,j,i)%particles(new_size))
	1155
	1156	grid_particles(k,j,i)%particles(1:number_of_particles) = tmp_particles_s(1:number_of_particles)
	1157	grid_particles(k,j,i)%particles(number_of_particles+1:new_size) = zero_particle
	1158
	1159	ELSE
	1160
	1161	ALLOCATE(tmp_particles_d(number_of_particles))
	1162	tmp_particles_d(1:number_of_particles) = grid_particles(k,j,i)%particles(1:number_of_particles)
	1163
	1164	DEALLOCATE(grid_particles(k,j,i)%particles)
	1165	ALLOCATE(grid_particles(k,j,i)%particles(new_size))
	1166
	1167	grid_particles(k,j,i)%particles(1:number_of_particles) = tmp_particles_d(1:number_of_particles)
	1168	grid_particles(k,j,i)%particles(number_of_particles+1:new_size) = zero_particle
	1169
	1170	DEALLOCATE(tmp_particles_d)
	1171
	1172	ENDIF
	1173
	1174	ENDIF
	1175	ENDDO
	1176	ENDDO
	1177	ENDDO
	1178
	1179	END SUBROUTINE dealloc_particles_array
	1180
	1181
[1359]	1182	END MODULE lpm_exchange_horiz_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |