Skip to content

Commit 84c1577

Browse files
wilfonbawilfonba
committed
kahan sum in void fraction BCs
1 parent 521d61c commit 84c1577

3 files changed

Lines changed: 100 additions & 53 deletions

File tree

‎src/common/m_boundary_common.fpp‎

Lines changed: 73 additions & 41 deletions
Original file line numberDiff line numberDiff line change
@@ -1194,26 +1194,27 @@ contains
11941194

11951195
end subroutine s_qbmm_extrapolation
11961196

1197-
impure subroutine s_populate_beta_buffers(q_beta, bc_type, nvar)
1197+
impure subroutine s_populate_beta_buffers(q_beta, bc_type, nvar, kcomp)
11981198

11991199
type(scalar_field), dimension(:), intent(inout) :: q_beta
12001200
type(integer_field), dimension(1:num_dims, 1:2), intent(in) :: bc_type
12011201
integer, intent(in) :: nvar
1202+
type(scalar_field), dimension(:), intent(inout) :: kcomp
12021203

12031204
integer :: k, l
12041205

12051206
!< x-direction
12061207
if (bc_x%beg >= 0) then
1207-
call s_mpi_reduce_beta_variables_buffers(q_beta, 1, -1, nvar)
1208+
call s_mpi_reduce_beta_variables_buffers(q_beta, 1, -1, nvar, kcomp)
12081209
else
12091210
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
12101211
do l = beta_bc_bounds(3)%beg, beta_bc_bounds(3)%end
12111212
do k = beta_bc_bounds(2)%beg, beta_bc_bounds(2)%end
12121213
select case (bc_x%beg)
12131214
case (BC_PERIODIC)
1214-
call s_beta_periodic(q_beta, 1, -1, k, l, nvar)
1215+
call s_beta_periodic(q_beta, 1, -1, k, l, nvar, kcomp)
12151216
case (BC_REFLECTIVE)
1216-
call s_beta_reflective(q_beta, 1, -1, k, l, nvar)
1217+
call s_beta_reflective(q_beta, 1, -1, k, l, nvar, kcomp)
12171218
case default
12181219
end select
12191220
end do
@@ -1222,16 +1223,16 @@ contains
12221223
end if
12231224

12241225
if (bc_x%end >= 0) then
1225-
call s_mpi_reduce_beta_variables_buffers(q_beta, 1, 1, nvar)
1226+
call s_mpi_reduce_beta_variables_buffers(q_beta, 1, 1, nvar, kcomp)
12261227
else
12271228
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
12281229
do l = beta_bc_bounds(3)%beg, beta_bc_bounds(3)%end
12291230
do k = beta_bc_bounds(2)%beg, beta_bc_bounds(2)%end
12301231
select case (bc_x%end)
12311232
case (BC_PERIODIC)
1232-
call s_beta_periodic(q_beta, 1, 1, k, l, nvar)
1233+
call s_beta_periodic(q_beta, 1, 1, k, l, nvar, kcomp)
12331234
case (BC_REFLECTIVE)
1234-
call s_beta_reflective(q_beta, 1, 1, k, l, nvar)
1235+
call s_beta_reflective(q_beta, 1, 1, k, l, nvar, kcomp)
12351236
case default
12361237
end select
12371238
end do
@@ -1241,16 +1242,16 @@ contains
12411242

12421243
!< y-direction
12431244
if (bc_y%beg >= 0) then
1244-
call s_mpi_reduce_beta_variables_buffers(q_beta, 2, -1, nvar)
1245+
call s_mpi_reduce_beta_variables_buffers(q_beta, 2, -1, nvar, kcomp)
12451246
else
12461247
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
12471248
do l = beta_bc_bounds(3)%beg, beta_bc_bounds(3)%end
12481249
do k = beta_bc_bounds(1)%beg, beta_bc_bounds(1)%end
12491250
select case (bc_y%beg)
12501251
case (BC_PERIODIC)
1251-
call s_beta_periodic(q_beta, 2, -1, k, l, nvar)
1252+
call s_beta_periodic(q_beta, 2, -1, k, l, nvar, kcomp)
12521253
case (BC_REFLECTIVE)
1253-
call s_beta_reflective(q_beta, 2, -1, k, l, nvar)
1254+
call s_beta_reflective(q_beta, 2, -1, k, l, nvar, kcomp)
12541255
case default
12551256
end select
12561257
end do
@@ -1259,16 +1260,16 @@ contains
12591260
end if
12601261

12611262
if (bc_y%end >= 0) then
1262-
call s_mpi_reduce_beta_variables_buffers(q_beta, 2, 1, nvar)
1263+
call s_mpi_reduce_beta_variables_buffers(q_beta, 2, 1, nvar, kcomp)
12631264
else
12641265
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
12651266
do l = beta_bc_bounds(3)%beg, beta_bc_bounds(3)%end
12661267
do k = beta_bc_bounds(1)%beg, beta_bc_bounds(1)%end
12671268
select case (bc_y%end)
12681269
case (BC_PERIODIC)
1269-
call s_beta_periodic(q_beta, 2, 1, k, l, nvar)
1270+
call s_beta_periodic(q_beta, 2, 1, k, l, nvar, kcomp)
12701271
case (BC_REFLECTIVE)
1271-
call s_beta_reflective(q_beta, 2, 1, k, l, nvar)
1272+
call s_beta_reflective(q_beta, 2, 1, k, l, nvar, kcomp)
12721273
case default
12731274
end select
12741275
end do
@@ -1281,16 +1282,16 @@ contains
12811282
#:if not MFC_CASE_OPTIMIZATION or num_dims > 2
12821283
!< z-direction
12831284
if (bc_z%beg >= 0) then
1284-
call s_mpi_reduce_beta_variables_buffers(q_beta, 3, -1, nvar)
1285+
call s_mpi_reduce_beta_variables_buffers(q_beta, 3, -1, nvar, kcomp)
12851286
else
12861287
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
12871288
do l = beta_bc_bounds(2)%beg, beta_bc_bounds(2)%end
12881289
do k = beta_bc_bounds(1)%beg, beta_bc_bounds(1)%end
12891290
select case (bc_type(3, 1)%sf(k, l, 0))
12901291
case (BC_PERIODIC)
1291-
call s_beta_periodic(q_beta, 3, -1, k, l, nvar)
1292+
call s_beta_periodic(q_beta, 3, -1, k, l, nvar, kcomp)
12921293
case (BC_REFLECTIVE)
1293-
call s_beta_reflective(q_beta, 3, -1, k, l, nvar)
1294+
call s_beta_reflective(q_beta, 3, -1, k, l, nvar, kcomp)
12941295
case default
12951296
end select
12961297
end do
@@ -1299,16 +1300,16 @@ contains
12991300
end if
13001301

13011302
if (bc_z%end >= 0) then
1302-
call s_mpi_reduce_beta_variables_buffers(q_beta, 3, 1, nvar)
1303+
call s_mpi_reduce_beta_variables_buffers(q_beta, 3, 1, nvar, kcomp)
13031304
else
13041305
$:GPU_PARALLEL_LOOP(private='[l,k]', collapse=2)
13051306
do l = beta_bc_bounds(2)%beg, beta_bc_bounds(2)%end
13061307
do k = beta_bc_bounds(1)%beg, beta_bc_bounds(1)%end
13071308
select case (bc_type(3, 2)%sf(k, l, 0))
13081309
case (BC_PERIODIC)
1309-
call s_beta_periodic(q_beta, 3, 1, k, l, nvar)
1310+
call s_beta_periodic(q_beta, 3, 1, k, l, nvar, kcomp)
13101311
case (BC_REFLECTIVE)
1311-
call s_beta_reflective(q_beta, 3, 1, k, l, nvar)
1312+
call s_beta_reflective(q_beta, 3, 1, k, l, nvar, kcomp)
13121313
case default
13131314
end select
13141315
end do
@@ -1319,58 +1320,69 @@ contains
13191320

13201321
end subroutine s_populate_beta_buffers
13211322

1322-
subroutine s_beta_periodic(q_beta, bc_dir, bc_loc, k, l, nvar)
1323+
subroutine s_beta_periodic(q_beta, bc_dir, bc_loc, k, l, nvar, kcomp)
13231324
$:GPU_ROUTINE(function_name='s_beta_periodic', &
13241325
& parallelism='[seq]', cray_inline=True)
13251326
type(scalar_field), dimension(num_dims + 1), intent(inout) :: q_beta
13261327
integer, intent(in) :: bc_dir, bc_loc
13271328
integer, intent(in) :: k, l
13281329
integer, intent(in) :: nvar
1330+
type(scalar_field), dimension(:), intent(inout) :: kcomp
13291331

13301332
integer :: j, i
1333+
real(wp) :: y_kahan, t_kahan
13311334

13321335
if (bc_dir == 1) then !< x-direction
13331336
if (bc_loc == -1) then !bc_x%beg
13341337
do i = 1, nvar
13351338
do j = -mapCells - 1, mapCells
1336-
q_beta(beta_vars(i))%sf(j, k, l) = q_beta(beta_vars(i))%sf(j, k, l) + &
1337-
q_beta(beta_vars(i))%sf(m + j + 1, k, l)
1339+
y_kahan = q_beta(beta_vars(i))%sf(m + j + 1, k, l) - kcomp(beta_vars(i))%sf(j, k, l)
1340+
t_kahan = q_beta(beta_vars(i))%sf(j, k, l) + y_kahan
1341+
kcomp(beta_vars(i))%sf(j, k, l) = (t_kahan - q_beta(beta_vars(i))%sf(j, k, l)) - y_kahan
1342+
q_beta(beta_vars(i))%sf(j, k, l) = t_kahan
13381343
end do
13391344
end do
13401345
else !< bc_x%end
13411346
do i = 1, nvar
13421347
do j = -mapcells, mapcells + 1
13431348
q_beta(beta_vars(i))%sf(m + j, k, l) = q_beta(beta_vars(i))%sf(j - 1, k, l)
1349+
kcomp(beta_vars(i))%sf(m + j, k, l) = kcomp(beta_vars(i))%sf(j - 1, k, l)
13441350
end do
13451351
end do
13461352
end if
13471353
elseif (bc_dir == 2) then !< y-direction
13481354
if (bc_loc == -1) then !< bc_y%beg
13491355
do i = 1, nvar
13501356
do j = -mapcells - 1, mapcells
1351-
q_beta(beta_vars(i))%sf(k, j, l) = q_beta(beta_vars(i))%sf(k, j, l) + &
1352-
q_beta(beta_vars(i))%sf(k, n + j + 1, l)
1357+
y_kahan = q_beta(beta_vars(i))%sf(k, n + j + 1, l) - kcomp(beta_vars(i))%sf(k, j, l)
1358+
t_kahan = q_beta(beta_vars(i))%sf(k, j, l) + y_kahan
1359+
kcomp(beta_vars(i))%sf(k, j, l) = (t_kahan - q_beta(beta_vars(i))%sf(k, j, l)) - y_kahan
1360+
q_beta(beta_vars(i))%sf(k, j, l) = t_kahan
13531361
end do
13541362
end do
13551363
else !< bc_y%end
13561364
do i = 1, nvar
13571365
do j = -mapcells, mapcells + 1
13581366
q_beta(beta_vars(i))%sf(k, n + j, l) = q_beta(beta_vars(i))%sf(k, j - 1, l)
1367+
kcomp(beta_vars(i))%sf(k, n + j, l) = kcomp(beta_vars(i))%sf(k, j - 1, l)
13591368
end do
13601369
end do
13611370
end if
13621371
elseif (bc_dir == 3) then !< z-direction
13631372
if (bc_loc == -1) then !< bc_z%beg
13641373
do i = 1, nvar
13651374
do j = -mapcells - 1, mapcells
1366-
q_beta(beta_vars(i))%sf(k, l, j) = q_beta(beta_vars(i))%sf(k, l, j) + &
1367-
q_beta(beta_vars(i))%sf(k, l, p + j + 1)
1375+
y_kahan = q_beta(beta_vars(i))%sf(k, l, p + j + 1) - kcomp(beta_vars(i))%sf(k, l, j)
1376+
t_kahan = q_beta(beta_vars(i))%sf(k, l, j) + y_kahan
1377+
kcomp(beta_vars(i))%sf(k, l, j) = (t_kahan - q_beta(beta_vars(i))%sf(k, l, j)) - y_kahan
1378+
q_beta(beta_vars(i))%sf(k, l, j) = t_kahan
13681379
end do
13691380
end do
13701381
else !< bc_z%end
13711382
do i = 1, nvar
13721383
do j = -mapcells, mapcells + 1
13731384
q_beta(beta_vars(i))%sf(k, l, p + j) = q_beta(beta_vars(i))%sf(k, l, j - 1)
1385+
kcomp(beta_vars(i))%sf(k, l, p + j) = kcomp(beta_vars(i))%sf(k, l, j - 1)
13741386
end do
13751387
end do
13761388
end if
@@ -1436,83 +1448,103 @@ contains
14361448

14371449
end subroutine s_beta_extrapolation
14381450

1439-
subroutine s_beta_reflective(q_beta, bc_dir, bc_loc, k, l, nvar)
1451+
subroutine s_beta_reflective(q_beta, bc_dir, bc_loc, k, l, nvar, kcomp)
14401452
$:GPU_ROUTINE(function_name='s_beta_reflective', &
14411453
& parallelism='[seq]', cray_inline=True)
14421454
type(scalar_field), dimension(num_dims + 1), intent(inout) :: q_beta
14431455
integer, intent(in) :: bc_dir, bc_loc
14441456
integer, intent(in) :: k, l
14451457
integer, intent(in) :: nvar
1458+
type(scalar_field), dimension(:), intent(inout) :: kcomp
14461459

14471460
integer :: j, i
1461+
real(wp) :: y_kahan, t_kahan
14481462

14491463
! Reflective BC for void fraction:
1450-
! 1) Fold ghost-cell contributions back onto their mirror interior cells
1451-
! 2) Set ghost cells = mirror of (now-folded) interior values
1464+
! 1) Fold ghost-cell contributions back onto their mirror interior cells (Kahan ADD)
1465+
! 2) Set ghost cells = mirror of (now-folded) interior values (REPLACE + copy comp)
14521466

14531467
if (bc_dir == 1) then !< x-direction
14541468
if (bc_loc == -1) then !< bc_x%beg
14551469
do i = 1, nvar
14561470
do j = 1, mapCells + 1
1457-
q_beta(beta_vars(i))%sf(j - 1, k, l) = q_beta(beta_vars(i))%sf(j - 1, k, l) + &
1458-
q_beta(beta_vars(i))%sf(-j, k, l)
1471+
y_kahan = q_beta(beta_vars(i))%sf(-j, k, l) - kcomp(beta_vars(i))%sf(j - 1, k, l)
1472+
t_kahan = q_beta(beta_vars(i))%sf(j - 1, k, l) + y_kahan
1473+
kcomp(beta_vars(i))%sf(j - 1, k, l) = (t_kahan - q_beta(beta_vars(i))%sf(j - 1, k, l)) - y_kahan
1474+
q_beta(beta_vars(i))%sf(j - 1, k, l) = t_kahan
14591475
end do
14601476
do j = 1, mapCells + 1
14611477
q_beta(beta_vars(i))%sf(-j, k, l) = q_beta(beta_vars(i))%sf(j - 1, k, l)
1478+
kcomp(beta_vars(i))%sf(-j, k, l) = kcomp(beta_vars(i))%sf(j - 1, k, l)
14621479
end do
14631480
end do
14641481
else !< bc_x%end
14651482
do i = 1, nvar
14661483
do j = 1, mapCells + 1
1467-
q_beta(beta_vars(i))%sf(m - (j - 1), k, l) = q_beta(beta_vars(i))%sf(m - (j - 1), k, l) + &
1468-
q_beta(beta_vars(i))%sf(m + j, k, l)
1484+
y_kahan = q_beta(beta_vars(i))%sf(m + j, k, l) - kcomp(beta_vars(i))%sf(m - (j - 1), k, l)
1485+
t_kahan = q_beta(beta_vars(i))%sf(m - (j - 1), k, l) + y_kahan
1486+
kcomp(beta_vars(i))%sf(m - (j - 1), k, l) = (t_kahan - q_beta(beta_vars(i))%sf(m - (j - 1), k, l)) - y_kahan
1487+
q_beta(beta_vars(i))%sf(m - (j - 1), k, l) = t_kahan
14691488
end do
14701489
do j = 1, mapCells + 1
14711490
q_beta(beta_vars(i))%sf(m + j, k, l) = q_beta(beta_vars(i))%sf(m - (j - 1), k, l)
1491+
kcomp(beta_vars(i))%sf(m + j, k, l) = kcomp(beta_vars(i))%sf(m - (j - 1), k, l)
14721492
end do
14731493
end do
14741494
end if
14751495
elseif (bc_dir == 2) then !< y-direction
14761496
if (bc_loc == -1) then !< bc_y%beg
14771497
do i = 1, nvar
14781498
do j = 1, mapCells + 1
1479-
q_beta(beta_vars(i))%sf(k, j - 1, l) = q_beta(beta_vars(i))%sf(k, j - 1, l) + &
1480-
q_beta(beta_vars(i))%sf(k, -j, l)
1499+
y_kahan = q_beta(beta_vars(i))%sf(k, -j, l) - kcomp(beta_vars(i))%sf(k, j - 1, l)
1500+
t_kahan = q_beta(beta_vars(i))%sf(k, j - 1, l) + y_kahan
1501+
kcomp(beta_vars(i))%sf(k, j - 1, l) = (t_kahan - q_beta(beta_vars(i))%sf(k, j - 1, l)) - y_kahan
1502+
q_beta(beta_vars(i))%sf(k, j - 1, l) = t_kahan
14811503
end do
14821504
do j = 1, mapCells + 1
14831505
q_beta(beta_vars(i))%sf(k, -j, l) = q_beta(beta_vars(i))%sf(k, j - 1, l)
1506+
kcomp(beta_vars(i))%sf(k, -j, l) = kcomp(beta_vars(i))%sf(k, j - 1, l)
14841507
end do
14851508
end do
14861509
else !< bc_y%end
14871510
do i = 1, nvar
14881511
do j = 1, mapCells + 1
1489-
q_beta(beta_vars(i))%sf(k, n - (j - 1), l) = q_beta(beta_vars(i))%sf(k, n - (j - 1), l) + &
1490-
q_beta(beta_vars(i))%sf(k, n + j, l)
1512+
y_kahan = q_beta(beta_vars(i))%sf(k, n + j, l) - kcomp(beta_vars(i))%sf(k, n - (j - 1), l)
1513+
t_kahan = q_beta(beta_vars(i))%sf(k, n - (j - 1), l) + y_kahan
1514+
kcomp(beta_vars(i))%sf(k, n - (j - 1), l) = (t_kahan - q_beta(beta_vars(i))%sf(k, n - (j - 1), l)) - y_kahan
1515+
q_beta(beta_vars(i))%sf(k, n - (j - 1), l) = t_kahan
14911516
end do
14921517
do j = 1, mapCells + 1
14931518
q_beta(beta_vars(i))%sf(k, n + j, l) = q_beta(beta_vars(i))%sf(k, n - (j - 1), l)
1519+
kcomp(beta_vars(i))%sf(k, n + j, l) = kcomp(beta_vars(i))%sf(k, n - (j - 1), l)
14941520
end do
14951521
end do
14961522
end if
14971523
elseif (bc_dir == 3) then !< z-direction
14981524
if (bc_loc == -1) then !< bc_z%beg
14991525
do i = 1, nvar
15001526
do j = 1, mapCells + 1
1501-
q_beta(beta_vars(i))%sf(k, l, j - 1) = q_beta(beta_vars(i))%sf(k, l, j - 1) + &
1502-
q_beta(beta_vars(i))%sf(k, l, -j)
1527+
y_kahan = q_beta(beta_vars(i))%sf(k, l, -j) - kcomp(beta_vars(i))%sf(k, l, j - 1)
1528+
t_kahan = q_beta(beta_vars(i))%sf(k, l, j - 1) + y_kahan
1529+
kcomp(beta_vars(i))%sf(k, l, j - 1) = (t_kahan - q_beta(beta_vars(i))%sf(k, l, j - 1)) - y_kahan
1530+
q_beta(beta_vars(i))%sf(k, l, j - 1) = t_kahan
15031531
end do
15041532
do j = 1, mapCells + 1
15051533
q_beta(beta_vars(i))%sf(k, l, -j) = q_beta(beta_vars(i))%sf(k, l, j - 1)
1534+
kcomp(beta_vars(i))%sf(k, l, -j) = kcomp(beta_vars(i))%sf(k, l, j - 1)
15061535
end do
15071536
end do
15081537
else !< bc_z%end
15091538
do i = 1, nvar
15101539
do j = 1, mapCells + 1
1511-
q_beta(beta_vars(i))%sf(k, l, p - (j - 1)) = q_beta(beta_vars(i))%sf(k, l, p - (j - 1)) + &
1512-
q_beta(beta_vars(i))%sf(k, l, p + j)
1540+
y_kahan = q_beta(beta_vars(i))%sf(k, l, p + j) - kcomp(beta_vars(i))%sf(k, l, p - (j - 1))
1541+
t_kahan = q_beta(beta_vars(i))%sf(k, l, p - (j - 1)) + y_kahan
1542+
kcomp(beta_vars(i))%sf(k, l, p - (j - 1)) = (t_kahan - q_beta(beta_vars(i))%sf(k, l, p - (j - 1))) - y_kahan
1543+
q_beta(beta_vars(i))%sf(k, l, p - (j - 1)) = t_kahan
15131544
end do
15141545
do j = 1, mapCells + 1
15151546
q_beta(beta_vars(i))%sf(k, l, p + j) = q_beta(beta_vars(i))%sf(k, l, p - (j - 1))
1547+
kcomp(beta_vars(i))%sf(k, l, p + j) = kcomp(beta_vars(i))%sf(k, l, p - (j - 1))
15161548
end do
15171549
end do
15181550
end if

0 commit comments

Comments
 (0)