Skip to content

Commit 71437fd

Browse files
danieljvickersdanieljvickers
committed
Modified the cuboid levelset to use local coords
1 parent ab29c8a commit 71437fd

File tree

1 file changed

+41
-73
lines changed

1 file changed

+41
-73
lines changed

‎src/common/m_compute_levelset.fpp‎

Lines changed: 41 additions & 73 deletions
Original file line numberDiff line numberDiff line change
@@ -287,9 +287,9 @@ contains
287287
(xy_local(2) > bottom_left(2) .and. xy_local(2) < top_right(2))) then
288288

289289
side_dists(1) = bottom_left(1) - xy_local(1)
290-
side_dists(2) = xy_local(1) - top_right(1)
290+
side_dists(2) = top_right(1) - xy_local(1)
291291
side_dists(3) = bottom_left(2) - xy_local(2)
292-
side_dists(4) = xy_local(2) - top_right(2)
292+
side_dists(4) = top_right(2) - xy_local(2)
293293
min_dist = initial_distance_buffer
294294
idx = 1
295295

@@ -312,8 +312,8 @@ contains
312312
abs(side_dists(idx))
313313
end if
314314
! convert the normal vector back into the global coordinate system
315-
levelset_norm%sf(i, j, 0, ib_patch_id, 1) = &
316-
matmul(rotation, levelset_norm%sf(i, j, 0, ib_patch_id, 1))
315+
levelset_norm%sf(i, j, 0, ib_patch_id, :) = &
316+
matmul(rotation, levelset_norm%sf(i, j, 0, ib_patch_id, :))
317317
else
318318
levelset_norm%sf(i, j, 0, ib_patch_id, :) = 0._wp
319319
end if
@@ -338,7 +338,7 @@ contains
338338
real(wp), dimension(1:3) :: xyz_local !< x and y coordinates in local IB frame
339339
real(wp), dimension(1:3, 1:3) :: rotation, inverse_rotation
340340

341-
integer :: i, j, k !< Loop index variables
341+
integer :: i, j, k, l, idx !< Loop index variables
342342

343343
length_x = patch_ib(ib_patch_id)%length_x
344344
length_y = patch_ib(ib_patch_id)%length_y
@@ -364,81 +364,49 @@ contains
364364
do j = 0, n
365365
do k = 0, p
366366

367-
x = x_cc(i)
368-
y = y_cc(j)
369-
z = z_cc(k)
370-
371-
if ((x > Left .and. x < Right) .or. &
372-
(y > Bottom .and. y < Top) .or. &
373-
(z > Back .and. z < Front)) then
374-
375-
side_dists(1) = Left - x
376-
side_dists(2) = x - Right
377-
side_dists(3) = Bottom - y
378-
side_dists(4) = y - Top
379-
side_dists(5) = Back - z
380-
side_dists(6) = z - Front
381-
382-
min_dist = minval(abs(side_dists))
383-
384-
if (f_approx_equal(min_dist, abs(side_dists(1)))) then
385-
levelset%sf(i, j, k, ib_patch_id) = side_dists(1)
386-
if (f_approx_equal(side_dists(1), 0._wp)) then
387-
levelset_norm%sf(i, j, k, ib_patch_id, 1) = 0._wp
388-
else
389-
levelset_norm%sf(i, j, k, ib_patch_id, 1) = side_dists(1)/ &
390-
abs(side_dists(1))
391-
end if
392-
393-
else if (f_approx_equal(min_dist, abs(side_dists(2)))) then
394-
levelset%sf(i, j, k, ib_patch_id) = side_dists(2)
395-
if (f_approx_equal(side_dists(2), 0._wp)) then
396-
levelset_norm%sf(i, j, k, ib_patch_id, 1) = 0._wp
397-
else
398-
levelset_norm%sf(i, j, k, ib_patch_id, 1) = -side_dists(2)/ &
399-
abs(side_dists(2))
400-
end if
401-
402-
else if (f_approx_equal(min_dist, abs(side_dists(3)))) then
403-
levelset%sf(i, j, k, ib_patch_id) = side_dists(3)
404-
if (f_approx_equal(side_dists(3), 0._wp)) then
405-
levelset_norm%sf(i, j, k, ib_patch_id, 2) = 0._wp
406-
else
407-
levelset_norm%sf(i, j, k, ib_patch_id, 2) = side_dists(3)/ &
408-
abs(side_dists(3))
409-
end if
410-
411-
else if (f_approx_equal(min_dist, abs(side_dists(4)))) then
412-
levelset%sf(i, j, k, ib_patch_id) = side_dists(4)
413-
if (f_approx_equal(side_dists(4), 0._wp)) then
414-
levelset_norm%sf(i, j, k, ib_patch_id, 2) = 0._wp
415-
else
416-
levelset_norm%sf(i, j, k, ib_patch_id, 2) = -side_dists(4)/ &
417-
abs(side_dists(4))
418-
end if
367+
xyz_local = [x_cc(i) - x_centroid, cart_y - y_centroid, cart_z - z_centroid] ! get coordinate frame centered on IB
368+
xyz_local = matmul(inverse_rotation, xyz_local) ! rotate the frame into the IB's coordiantes
419369
420-
else if (f_approx_equal(min_dist, abs(side_dists(5)))) then
421-
levelset%sf(i, j, k, ib_patch_id) = side_dists(5)
422-
if (f_approx_equal(side_dists(5), 0._wp)) then
423-
levelset_norm%sf(i, j, k, ib_patch_id, 3) = 0._wp
424-
else
425-
levelset_norm%sf(i, j, k, ib_patch_id, 3) = side_dists(5)/ &
426-
abs(side_dists(5))
370+
if ((xyz_local(1) > Left .and. xyz_local(1) < Right) .or. &
371+
(xyz_local(2) > Bottom .and. xyz_local(2) < Top) .or. &
372+
(xyz_local(3) > Back .and. xyz_local(3) < Front)) then
373+
374+
side_dists(1) = Left - xyz_local(1)
375+
side_dists(2) = Right - xyz_local(1)
376+
side_dists(3) = Bottom - xyz_local(2)
377+
side_dists(4) = Top - xyz_local(2)
378+
side_dists(5) = Back - xyz_local(3)
379+
side_dists(6) = Front - xyz_local(3)
380+
381+
! get the minimal_distance
382+
min_dist = side_dists(1)
383+
idx = 1
384+
do l = 2, 6
385+
if (abs(side_dists(l)) < abs(min_dist)) then
386+
idx = l
387+
min_dist = side_dists(idx)
427388
end if
389+
end do
428390
429-
else if (f_approx_equal(min_dist, abs(side_dists(6)))) then
430-
levelset%sf(i, j, k, ib_patch_id) = side_dists(6)
431-
if (f_approx_equal(side_dists(6), 0._wp)) then
432-
levelset_norm%sf(i, j, k, ib_patch_id, 3) = 0._wp
391+
levelset%sf(i, j, k, ib_patch_id) = min_dist
392+
if (f_approx_equal(min_dist, 0._wp)) then
393+
levelset_norm%sf(i, j, k, ib_patch_id, :) = 0._wp
394+
else
395+
if (idx == 1 .or. idx == 2) then
396+
levelset_norm%sf(i, j, k, ib_patch_id, 1) = side_dists(idx)/ &
397+
abs(side_dists(idx))
398+
else if (idx == 3 .or. idx == 4) then
399+
levelset_norm%sf(i, j, k, ib_patch_id, 2) = side_dists(idx)/ &
400+
abs(side_dists(idx))
433401
else
434-
levelset_norm%sf(i, j, k, ib_patch_id, 3) = -side_dists(6)/ &
435-
abs(side_dists(6))
402+
levelset_norm%sf(i, j, k, ib_patch_id, 2) = side_dists(idx)/ &
403+
abs(side_dists(idx))
436404
end if
437-
405+
! convert the normal vector back into the global coordinate system
406+
levelset_norm%sf(i, j, k, ib_patch_id, :) = &
407+
matmul(rotation, levelset_norm%sf(i, j, k, ib_patch_id, :))
438408
end if
439-
440409
end if
441-
442410
end do
443411
end do
444412
end do

0 commit comments

Comments
 (0)