From 420c339e2faee868526e35583d234db75b34ad14 Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Fri, 21 Aug 2026 11:40:35 -0600 Subject: [PATCH 1/5] Add timers --- src/core_atmosphere/dynamics/mpas_atm_time_integration.F | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index 227fbde862..7fb44bef01 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -3972,7 +3972,7 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart end if !$OMP BARRIER - + call mpas_timer_start('atm_advance_acoustic_step_3976') !$acc parallel default(present) !$acc loop gang worker private(ts,rs) do iCell=cellSolveStart,cellSolveEnd ! loop over all owned cells to solve @@ -4104,7 +4104,7 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart end do ! end of loop over cells !$acc end parallel - + call mpas_timer_stop('atm_advance_acoustic_step_3976') end subroutine atm_advance_acoustic_step_work From c8a070433b47b3a38eb2bf2a54e2f6c239a3268a Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 25 Aug 2026 16:58:50 -0600 Subject: [PATCH 2/5] Opt 1 - Split up parallel region into two separate ones This commit optimizes the GPU performance of the relevant kernel in atm_advance_acoustic_step_work by splitting up a single parallel region into two separate regions. A single larger parallel region, especially with the presence of worker-level private variables ts and rs, results in higher register usage per thread, which in turn reduces the occupancy (number of threads executing at any given time on the GPU relative to max capacity) on the GPU. Splitting up into two parallel regionsm, and removing private variables ts and rs from the second region, reduces the register usage and improves the occupancy and performance in this instance. It is also required to save rs and ts to rho_pp and rtheta_pp respectively, in order for the second loop to retain the respective values. This optimization degrades performance on CPUs, and is addressed in the subsequent commit, Co-authored-by: Pranay Reddy Kommera --- .../dynamics/mpas_atm_time_integration.F | 18 ++++++++++++++++-- 1 file changed, 16 insertions(+), 2 deletions(-) diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index 7fb44bef01..f7b30d2480 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -4046,8 +4046,22 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart + cofwt(k-1,iCell)*(etp(k-1)*ts(k-1)+etm(k-1)*rtheta_pp(k-1,iCell))) end do + !$acc loop vector + do k=1,nVertLevels + rho_pp(k,iCell) = rs(k) + rtheta_pp(k,iCell) = ts(k) + end do + + end if + enddo + !$acc end parallel + ! tridiagonal solve sweeping up and then down the column + !$acc parallel default(present) + !$acc loop gang worker + do iCell=cellSolveStart,cellSolveEnd ! loop over all owned cells to solve + if(specZoneMaskCell(iCell) == 0.0) then ! not specified zone, compute... !MGD VECTOR DEPENDENCE !$acc loop seq do k=2,nVertLevels @@ -4084,9 +4098,9 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart !DIR$ IVDEP !$acc loop vector do k=1,nVertLevels - rho_pp(k,iCell) = rs(k) - dts*cofrz(k) *( ewp(k+1)*rw_p(k+1,iCell) & + rho_pp(k,iCell) = rho_pp(k,iCell) - dts*cofrz(k) *( ewp(k+1)*rw_p(k+1,iCell) & -ewp(k )*rw_p(k ,iCell)) - rtheta_pp(k,iCell) = ts(k) - dts*rdzw(k)*( ewp(k+1)*coftz(k+1,iCell)*rw_p(k+1,iCell) & + rtheta_pp(k,iCell) = rtheta_pp(k,iCell) - dts*rdzw(k)*( ewp(k+1)*coftz(k+1,iCell)*rw_p(k+1,iCell) & -ewp(k )*coftz(k ,iCell)*rw_p(k ,iCell)) end do From 167ed11bca5118a1992b059344baa379a217d8bf Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Thu, 27 Aug 2026 17:05:33 -0600 Subject: [PATCH 3/5] Introduce preprocessor directive to selectively split loops for GPU builds The previous commit split up a parallel region into two as it is more optimal on GPUs, however this degrades performance on CPUs. This commit introduces a preprocessor directive that conditionally splits up the loops only in the case of an OPENACC GPU build, in order to preserve the current performance on CPUs. Co-authored-by: Michael Duda --- src/core_atmosphere/dynamics/mpas_atm_time_integration.F | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index f7b30d2480..b7dcbdf239 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -4052,16 +4052,19 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart rtheta_pp(k,iCell) = ts(k) end do +#ifdef MPAS_OPENACC end if enddo !$acc end parallel - ! tridiagonal solve sweeping up and then down the column + !$acc parallel default(present) !$acc loop gang worker do iCell=cellSolveStart,cellSolveEnd ! loop over all owned cells to solve if(specZoneMaskCell(iCell) == 0.0) then ! not specified zone, compute... +#endif + ! tridiagonal solve sweeping up and then down the column !MGD VECTOR DEPENDENCE !$acc loop seq do k=2,nVertLevels From ea3604b07094c679021cdb648749c439411f84c9 Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Fri, 28 Aug 2026 13:28:22 -0600 Subject: [PATCH 4/5] Opt 2 + introduce macros The second optimization involves swapping the sequential and vector loops in order to increase the DRAM bandwidth on GPUs. This optimization, however, results in performance degradation on CPUS. Hence, two preprocessor macros are introduced to try to preserve the same loop order on CPUs but interchange the loops on GPUs. Co-authored-by: Pranay Reddy Kommera --- .../dynamics/mpas_atm_time_integration.F | 22 ++++++++++++++----- 1 file changed, 17 insertions(+), 5 deletions(-) diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index b7dcbdf239..6de9c386b3 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -9,9 +9,13 @@ #ifdef MPAS_OPENACC #define MPAS_ACC_TIMER_START(X) call mpas_timer_start(X) #define MPAS_ACC_TIMER_STOP(X) call mpas_timer_stop(X) +#define NEDGESONCELL_INNER_LOOP() do i=1,nEdgesOnCell(iCell) +#define NEDGESONCELL_OUTER_LOOP() #else #define MPAS_ACC_TIMER_START(X) #define MPAS_ACC_TIMER_STOP(X) +#define NEDGESONCELL_OUTER_LOOP() do i=1,nEdgesOnCell(iCell) +#define NEDGESONCELL_INNER_LOOP() #endif module atm_time_integration @@ -3997,14 +4001,22 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart rs(k) = 0.0 end do - !$acc loop seq - do i=1,nEdgesOnCell(iCell) - iEdge = edgesOnCell(i,iCell) - cell1 = cellsOnEdge(1,iEdge) - cell2 = cellsOnEdge(2,iEdge) + ! When building this code for CPUs, the below macro will expand to "do i=1,nEdgesOnCell(iCell)", + ! preserving the nested loop order as originally written, and optimized for CPU execution. + ! However, when building for GPUs (presently with OpenACC), the order of the nested loops are + ! effectively interchanged to yield better GPU performance, without affecting correctness. In + ! the latter case, NEDGESONCELL_OUTER_LOOP expands to an empty string and NEDGESONCELL_INNER_LOOP + ! now expands to "do i=1,nEdgesOnCell(iCell)". + NEDGESONCELL_OUTER_LOOP() !DIR$ IVDEP !$acc loop vector do k=1,nVertLevels + !$acc loop seq + NEDGESONCELL_INNER_LOOP() + + iEdge = edgesOnCell(i,iCell) + cell1 = cellsOnEdge(1,iEdge) + cell2 = cellsOnEdge(2,iEdge) flux = edgesOnCell_sign(i,iCell)*dts*dvEdge(iEdge)*ru_p(k,iEdge) * invAreaCell(iCell) rs(k) = rs(k)-flux ts(k) = ts(k)-flux*0.5*(theta_m(k,cell2)+theta_m(k,cell1)) From a242478f8ea5fc3846225f3fcc66cbc932b76fda Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 8 Sep 2026 10:44:48 -0600 Subject: [PATCH 5/5] Opt 3 - Fuse two vector loops This commit fuses two vertical level vector loops in order to optimize GPU performance. It does not noticeably degrade CPU performance. Co-authored-by: Pranay Reddy Kommera --- src/core_atmosphere/dynamics/mpas_atm_time_integration.F | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index 6de9c386b3..9e88e115c7 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -4099,12 +4099,8 @@ subroutine atm_advance_acoustic_step_work(nCells, nEdges, nCellsSolve, cellStart *(fzm(k)*rho_zz(k,iCell)+fzp(k)*rho_zz(k-1,iCell)) & *w(k,iCell) )/(1.0+dts*dss(k,iCell)) & - (rw_save(k ,iCell) - rw(k ,iCell)) - end do - + ! accumulate (rho*omega)' for use later in scalar transport -!DIR$ IVDEP - !$acc loop vector - do k=2,nVertLevels wwAvg(k,iCell) = wwAvg(k,iCell) + ewp(k)*rw_p(k,iCell) end do