diff --git a/src/simulation/m_rhs.fpp b/src/simulation/m_rhs.fpp
index 7e74b33ed..66c0fe934 100644
--- a/src/simulation/m_rhs.fpp
+++ b/src/simulation/m_rhs.fpp
@@ -599,11 +599,12 @@ contains
 
         if (.not. igr) then
             ! Association/Population of Working Variables
-            $:GPU_PARALLEL_LOOP(private='[i, j, k, l]', collapse=4)
-            do i = 1, sys_size
-                do l = idwbuff(3)%beg, idwbuff(3)%end
-                    do k = idwbuff(2)%beg, idwbuff(2)%end
-                        do j = idwbuff(1)%beg, idwbuff(1)%end
+            $:GPU_PARALLEL_LOOP(private='[i, j, k, l]', collapse=3)
+            do l = idwbuff(3)%beg, idwbuff(3)%end
+                do k = idwbuff(2)%beg, idwbuff(2)%end
+                    do j = idwbuff(1)%beg, idwbuff(1)%end
+                        $:GPU_LOOP(parallelism='[seq]')
+                        do i = 1, sys_size
                             q_cons_qp%vf(i)%sf(j, k, l) = q_cons_vf(i)%sf(j, k, l)
                         end do
                     end do
@@ -1124,12 +1125,13 @@ contains
             end if
 
             if (hypo_nc_mode /= hypo_nc_mode_dual_pass) then
-                $:GPU_PARALLEL_LOOP(collapse=4,private='[j, k_loop, l_loop, q_loop, inv_ds, flux_face1, flux_face2]')
-                do j = 1, sys_size
-                    do q_loop = 0, p
-                        do l_loop = 0, n
-                            do k_loop = 0, m
-                                inv_ds = 1._wp/dx(k_loop)
+                $:GPU_PARALLEL_LOOP(collapse=3,private='[j, k_loop, l_loop, q_loop, inv_ds, flux_face1, flux_face2]')
+                do q_loop = 0, p
+                    do l_loop = 0, n
+                        do k_loop = 0, m
+                            inv_ds = 1._wp/dx(k_loop)
+                            $:GPU_LOOP(parallelism='[seq]')
+                            do j = 1, sys_size
                                 flux_face1 = flux_n(1)%vf(j)%sf(k_loop - 1, l_loop, q_loop)
                                 flux_face2 = flux_n(1)%vf(j)%sf(k_loop, l_loop, q_loop)
                                 rhs_vf(j)%sf(k_loop, l_loop, q_loop) = inv_ds*(flux_face1 - flux_face2)
@@ -1200,12 +1202,13 @@ contains
             end if
 
             if (hypo_nc_mode /= hypo_nc_mode_dual_pass) then
-                $:GPU_PARALLEL_LOOP(collapse=4,private='[j, k, l, q, inv_ds, flux_face1, flux_face2]')
-                do j = 1, sys_size
-                    do l = 0, p
-                        do k = 0, n
-                            do q = 0, m
-                                inv_ds = 1._wp/dy(k)
+                $:GPU_PARALLEL_LOOP(collapse=3,private='[j, k, l, q, inv_ds, flux_face1, flux_face2]')
+                do l = 0, p
+                    do k = 0, n
+                        do q = 0, m
+                            inv_ds = 1._wp/dy(k)
+                            $:GPU_LOOP(parallelism='[seq]')
+                            do j = 1, sys_size
                                 flux_face1 = flux_n(2)%vf(j)%sf(q, k - 1, l)
                                 flux_face2 = flux_n(2)%vf(j)%sf(q, k, l)
                                 rhs_vf(j)%sf(q, k, l) = rhs_vf(j)%sf(q, k, l) + inv_ds*(flux_face1 - flux_face2)
@@ -1362,12 +1365,13 @@ contains
                 $:END_GPU_PARALLEL_LOOP()
             else  ! Cartesian Coordinates
                 if (hypo_nc_mode /= hypo_nc_mode_dual_pass) then
-                    $:GPU_PARALLEL_LOOP(collapse=4,private='[j, k, l, q, inv_ds, flux_face1, flux_face2]')
-                    do j = 1, sys_size
-                        do k = 0, p
-                            do q = 0, n
-                                do l = 0, m
-                                    inv_ds = 1._wp/dz(k)
+                    $:GPU_PARALLEL_LOOP(collapse=3,private='[j, k, l, q, inv_ds, flux_face1, flux_face2]')
+                    do k = 0, p
+                        do q = 0, n
+                            do l = 0, m
+                                inv_ds = 1._wp/dz(k)
+                                $:GPU_LOOP(parallelism='[seq]')
+                                do j = 1, sys_size
                                     flux_face1 = flux_n(3)%vf(j)%sf(l, q, k - 1)
                                     flux_face2 = flux_n(3)%vf(j)%sf(l, q, k)
                                     rhs_vf(j)%sf(l, q, k) = rhs_vf(j)%sf(l, q, k) + inv_ds*(flux_face1 - flux_face2)
diff --git a/src/simulation/m_riemann_state.fpp b/src/simulation/m_riemann_state.fpp
index aeb34c8a9..1c589d4f7 100644
--- a/src/simulation/m_riemann_state.fpp
+++ b/src/simulation/m_riemann_state.fpp
@@ -1345,11 +1345,12 @@ contains
                 ! Reshaping Outputted Data in y-direction
 
                 if (norm_dir == 2) then
-                    $:GPU_PARALLEL_LOOP(collapse=4)
-                    do i = 1, sys_size
-                        do l = is3%beg, is3%end
-                            do j = is1%beg, is1%end
-                                do k = is2%beg, is2%end
+                    $:GPU_PARALLEL_LOOP(collapse=3)
+                    do l = is3%beg, is3%end
+                        do j = is1%beg, is1%end
+                            do k = is2%beg, is2%end
+                                $:GPU_LOOP(parallelism='[seq]')
+                                do i = 1, sys_size
                                     flux_vf(i)%sf(k, j, l) = flux_${INFIX}$rsx_vf(k, j, l, i)
                                 end do
                             end do
@@ -1401,11 +1402,12 @@ contains
                     #:endif
                     ! Reshaping Outputted Data in z-direction
                 else if (norm_dir == 3) then
-                    $:GPU_PARALLEL_LOOP(collapse=4)
-                    do i = 1, sys_size
-                        do j = is1%beg, is1%end
-                            do k = is2%beg, is2%end
-                                do l = is3%beg, is3%end
+                    $:GPU_PARALLEL_LOOP(collapse=3)
+                    do j = is1%beg, is1%end
+                        do k = is2%beg, is2%end
+                            do l = is3%beg, is3%end
+                                $:GPU_LOOP(parallelism='[seq]')
+                                do i = 1, sys_size
                                     flux_vf(i)%sf(l, k, j) = flux_${INFIX}$rsx_vf(l, k, j, i)
                                 end do
                             end do
@@ -1455,11 +1457,12 @@ contains
                         end if
                     #:endif
                 else if (norm_dir == 1) then
-                    $:GPU_PARALLEL_LOOP(collapse=4)
-                    do i = 1, sys_size
-                        do l = is3%beg, is3%end
-                            do k = is2%beg, is2%end
-                                do j = is1%beg, is1%end
+                    $:GPU_PARALLEL_LOOP(collapse=3)
+                    do l = is3%beg, is3%end
+                        do k = is2%beg, is2%end
+                            do j = is1%beg, is1%end
+                                $:GPU_LOOP(parallelism='[seq]')
+                                do i = 1, sys_size
                                     flux_vf(i)%sf(j, k, l) = flux_${INFIX}$rsx_vf(j, k, l, i)
                                 end do
                             end do
diff --git a/src/simulation/m_weno.fpp b/src/simulation/m_weno.fpp
index 89cc45f77..70e838379 100644
--- a/src/simulation/m_weno.fpp
+++ b/src/simulation/m_weno.fpp
@@ -868,11 +868,12 @@ contains
         type(scalar_field), dimension(1:), intent(in) :: v_vf
         integer                                       :: i, j, k, l, n_vars
 
-        $:GPU_PARALLEL_LOOP(collapse=4)
-        do i = 1, v_size
-            do l = idwbuff(3)%beg, idwbuff(3)%end
-                do k = idwbuff(2)%beg, idwbuff(2)%end
-                    do j = idwbuff(1)%beg, idwbuff(1)%end
+        $:GPU_PARALLEL_LOOP(collapse=3)
+        do l = idwbuff(3)%beg, idwbuff(3)%end
+            do k = idwbuff(2)%beg, idwbuff(2)%end
+                do j = idwbuff(1)%beg, idwbuff(1)%end
+                    $:GPU_LOOP(parallelism='[seq]')
+                    do i = 1, v_size
                         v_rs_weno(j, k, l, i) = v_vf(i)%sf(j, k, l)
                     end do
                 end do