Skip to content

Commit 04a5bf0

Browse files
author
Daniel Vickers
committed
Some loop optimizations and debug multi-particle periodicity
1 parent 4380b0a commit 04a5bf0

1 file changed

Lines changed: 32 additions & 9 deletions

File tree

src/simulation/m_ib_patches.fpp

Lines changed: 32 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -66,6 +66,9 @@ contains
6666
! find the indices to the left and right of the IB in i, j, k
6767
call get_bounding_indices(center, patch_id, il, ir, jl, jr, kl, kr)
6868

69+
! skip patches whose bounding box does not overlap this rank's domain
70+
if (ir < il .or. jr < jl .or. kr < kl) cycle
71+
6972
$:GPU_PARALLEL_LOOP(private='[i, j, k, xyz_local, length, radius, airfoil_id, eta]', &
7073
& copyin='[patch_id, encoded_patch_id, center]', collapse=3)
7174
do k = kl, kr
@@ -80,6 +83,7 @@ contains
8083
if (patch_ib(patch_id)%geometry == 8) then
8184
! sphere geometry
8285
radius = patch_ib(patch_id)%radius
86+
8387
if (f_is_inside_sphere(xyz_local(1), xyz_local(2), xyz_local(3), &
8488
& radius)) ib_markers%sf(i, j, k) = encoded_patch_id
8589
else if (patch_ib(patch_id)%geometry == 9) then
@@ -134,6 +138,9 @@ contains
134138
! find the indices to the left and right of the IB in i, j, k
135139
call get_bounding_indices(center, patch_id, il, ir, jl, jr, kl, kr)
136140

141+
! skip patches whose bounding box does not overlap this rank's domain
142+
if (ir < il .or. jr < jl) cycle
143+
137144
$:GPU_PARALLEL_LOOP(private='[i, j, xyz_local, airfoil_id, eta, length]', copyin='[patch_id, &
138145
& encoded_patch_id, center]', collapse=2)
139146
do j = jl, jr
@@ -199,8 +206,8 @@ contains
199206
do xp = xp_lower, xp_upper
200207
do yp = yp_lower, yp_upper
201208
do zp = zp_lower, zp_upper
202-
$:GPU_PARALLEL_LOOP(private='[xp, yp, zp, i, il, ir, j, jl, jr, k, kl, kr, xyz_local, length, radius, &
203-
& patch_id, airfoil_id, model_id, encoded_patch_id, center, eta]')
209+
$:GPU_PARALLEL_LOOP(private='[i, il, ir, j, jl, jr, k, kl, kr, xyz_local, length, radius, patch_id, &
210+
& airfoil_id, model_id, encoded_patch_id, center, eta]', copyin='[xp, yp, zp]')
204211
do patch_id = 1, num_ibs
205212
center(1) = patch_ib(patch_id)%x_centroid + real(xp, wp)*(x_domain%end - x_domain%beg)
206213
center(2) = patch_ib(patch_id)%y_centroid + real(yp, wp)*(y_domain%end - y_domain%beg)
@@ -224,6 +231,7 @@ contains
224231
if (patch_ib(patch_id)%geometry == 8) then
225232
! sphere geometry
226233
radius = patch_ib(patch_id)%radius
234+
227235
if (f_is_inside_sphere(xyz_local(1), xyz_local(2), xyz_local(3), &
228236
& radius)) ib_markers%sf(i, j, k) = encoded_patch_id
229237
else if (patch_ib(patch_id)%geometry == 9) then
@@ -267,8 +275,8 @@ contains
267275
268276
do xp = xp_lower, xp_upper
269277
do yp = yp_lower, yp_upper
270-
$:GPU_PARALLEL_LOOP(private='[xp, yp, i, il, ir, j, jl, jr, xyz_local, length, radius, patch_id, airfoil_id, &
271-
& model_id, encoded_patch_id, center, eta]')
278+
$:GPU_PARALLEL_LOOP(private='[i, il, ir, j, jl, jr, xyz_local, length, radius, patch_id, airfoil_id, &
279+
& model_id, encoded_patch_id, center, eta]', copyin='[xp, yp]')
272280
do patch_id = 1, num_ibs
273281
center(1) = patch_ib(patch_id)%x_centroid + real(xp, wp)*(x_domain%end - x_domain%beg)
274282
center(2) = patch_ib(patch_id)%y_centroid + real(yp, wp)*(y_domain%end - y_domain%beg)
@@ -293,31 +301,31 @@ contains
293301
! circular geometries
294302
radius = patch_ib(patch_id)%radius
295303
if (f_is_inside_cylinder(xyz_local(1), xyz_local(2), 0._wp, radius, 0._wp)) ib_markers%sf(i, &
296-
& j, k) = encoded_patch_id
304+
& j, 0) = encoded_patch_id
297305
else if (patch_ib(patch_id)%geometry == 3) then
298306
! rectangular geometries
299307
length = [patch_ib(patch_id)%length_x, patch_ib(patch_id)%length_y, 0._wp]
300308
if (f_is_inside_cuboid(xyz_local(1), xyz_local(2), xyz_local(3), length)) ib_markers%sf(i, j, &
301-
& k) = encoded_patch_id
309+
& 0) = encoded_patch_id
302310
else if (patch_ib(patch_id)%geometry == 4) then
303311
! 2D airfoil geometry
304312
airfoil_id = patch_ib(patch_id)%airfoil_id
305313
xyz_local = xyz_local - patch_ib(patch_id)%centroid_offset
306314
if (f_is_inside_airfoil(xyz_local(1), xyz_local(2), 0._wp, 0._wp, &
307-
& airfoil_id)) ib_markers%sf(i, j, k) = encoded_patch_id
315+
& airfoil_id)) ib_markers%sf(i, j, 0) = encoded_patch_id
308316
else if (patch_ib(patch_id)%geometry == 5) then
309317
! STL model geometry
310318
xyz_local = xyz_local - patch_ib(patch_id)%centroid_offset
311319
model_id = patch_ib(patch_id)%model_id
312320
eta = f_model_is_inside(gpu_ntrs(model_id), model_id, xyz_local)
313321
if (eta > stl_models(model_id)%model_threshold) then
314-
ib_markers%sf(i, j, k) = encoded_patch_id
322+
ib_markers%sf(i, j, 0) = encoded_patch_id
315323
end if
316324
else if (patch_ib(patch_id)%geometry == 6) then
317325
! ellipse geometry
318326
length = [patch_ib(patch_id)%length_x, patch_ib(patch_id)%length_y, 0._wp]
319327
if (f_is_inside_ellipse(xyz_local(1), xyz_local(2), length)) ib_markers%sf(i, j, &
320-
& k) = encoded_patch_id
328+
& 0) = encoded_patch_id
321329
end if
322330
end do
323331
end do
@@ -460,6 +468,7 @@ contains
460468
real(wp), dimension(3) :: bbox_min, bbox_max, local_corner, world_corner
461469
real(wp), dimension(2) :: lx, ly, lz
462470
integer :: cx, cy, cz
471+
logical :: outside_domain
463472

464473
if (patch_ib(patch_id)%geometry == 2 .or. patch_ib(patch_id)%geometry == 8) then
465474
! circle and sphere geometries
@@ -535,6 +544,20 @@ contains
535544
end do
536545
end if
537546

547+
! completely skip patches whose bounding box does not overlap this rank's domain
548+
outside_domain = bbox_min(1) > x_cc(m + gp_layers + 1) .or. bbox_max(1) < x_cc(-gp_layers - 1) .or. bbox_min(2) > y_cc(n &
549+
& + gp_layers + 1) .or. bbox_max(2) < y_cc(-gp_layers - 1)
550+
if (num_dims == 3) then
551+
outside_domain = outside_domain .or. bbox_min(3) > z_cc(p + gp_layers + 1) .or. bbox_max(3) < z_cc(-gp_layers - 1)
552+
end if
553+
554+
if (outside_domain) then
555+
il = 1; ir = 0
556+
jl = 1; jr = 0
557+
kl = 1; kr = 0
558+
return
559+
end if
560+
538561
il = -gp_layers - 1
539562
jl = -gp_layers - 1
540563
kl = -gp_layers - 1

0 commit comments

Comments
 (0)