Add non-sbp upwinding operators. These have the same order as the centered

operators but does not satisfy any SBP property.


git-svn-id: https://svn.cct.lsu.edu/repos/numrel/LSUThorns/SummationByParts/trunk@114 f69c4107-0314-4c4f-9ad4-17e986b73f4a

7 files changed, 750 insertions, 17 deletions

diff --git a/param.ccl b/param.ccl
index 2ffd8d2..eba0e40 100644
--- a/param.ccl
+++ b/param.ccl
@@ -17,6 +17,10 @@ BOOLEAN sbp_2nd_deriv "Should the 2nd derivative operator be SBP" STEERABLE=alwa
 {
 } "yes"
   
+BOOLEAN sbp_upwind_deriv "Should the upwind derivative operator be SBP" STEERABLE=always
+{
+} "yes"
+  
 INT order "Order of accuracy" STEERABLE=always
 {
   2:8:2 :: ""
diff --git a/src/All_Coeffs_mod.F90 b/src/All_Coeffs_mod.F90
index 539e46e..1886497 100644
--- a/src/All_Coeffs_mod.F90
+++ b/src/All_Coeffs_mod.F90
@@ -59,6 +59,41 @@ module All_Coeffs_mod
 
     end subroutine coeffs_up_2_1
 
+    subroutine coeffs_up_2 ( a1, q1, a2, q2 )
+
+      CCTK_REAL, dimension(-2:2), intent(OUT) :: a1, a2
+      CCTK_REAL, dimension(3,2), intent(OUT) :: q1, q2
+
+      a1(-2) = 0.5_wp
+      a1(-1) = -2.0_wp
+      a1(0) = 1.5_wp
+      a1(1) = zero
+      a1(2) = zero
+
+      q1(1,1) = zero
+      q1(2,1) = zero
+      q1(3,1) = zero
+
+      q1(1,2) = -1.0_wp
+      q1(2,2) = 1.0_wp
+      q1(3,2) = zero
+
+      a2(-2) = zero
+      a2(-1) = zero
+      a2(0) = -1.5_wp
+      a2(1) = 2.0_wp
+      a2(2) = -0.5_wp
+
+      q2(1,1) = zero
+      q2(2,1) = zero
+      q2(3,1) = zero
+
+      q2(1,2) = zero
+      q2(2,2) = -1.0_wp
+      q2(3,2) = 1.0_wp
+
+    end subroutine coeffs_up_2
+
     subroutine coeffs_2_2_1 ( a, q )
 
       CCTK_REAL, dimension(2), intent(OUT) :: a
@@ -201,6 +236,59 @@ module All_Coeffs_mod
 
     end subroutine coeffs_up_4_2
 
+    subroutine coeffs_up_4 ( a1, q1, a2, q2 )
+
+      CCTK_REAL, dimension(-3:3), intent(OUT) :: a1, a2
+      CCTK_REAL, dimension(4,3), intent(OUT) :: q1, q2
+
+      a1(-3) = -0.08333333333333333333333333333333333333333_wp
+      a1(-2) = 0.5_wp
+      a1(-1) = -1.5_wp
+      a1(0) = 0.8333333333333333333333333333333333333333_wp
+      a1(1) = 0.25_wp
+      a1(2) = zero
+      a1(3) = zero
+
+      q1(1,1) = zero
+      q1(2,1) = zero
+      q1(3,1) = zero
+      q1(4,1) = zero
+
+      q1(1,2) = -1.0_wp
+      q1(2,2) = 1.0_wp
+      q1(3,2) = zero
+      q1(4,2) = zero
+
+      q1(1,3) = 0.5_wp
+      q1(2,3) = -2.0_wp
+      q1(3,3) = 1.5_wp
+      q1(4,3) = zero
+
+      a2(-3) = zero
+      a2(-2) = zero
+      a2(-1) = -0.25_wp
+      a2(0) = -0.8333333333333333333333333333333333333333_wp
+      a2(1) = 1.5_wp
+      a2(2) = -0.5_wp
+      a2(3) = 0.08333333333333333333333333333333333333333_wp
+
+      q2(1,1) = zero
+      q2(2,1) = zero
+      q2(3,1) = zero
+      q2(4,1) = zero
+
+      q2(1,2) = zero
+      q2(2,2) = zero
+      q2(3,2) = -1.0_wp
+      q2(4,2) = 1.0_wp
+
+      q2(1,3) = zero
+      q2(2,3) = -1.5_wp
+      q2(3,3) = 2.0_wp
+      q2(4,3) = -0.5_wp
+
+    end subroutine coeffs_up_4
+
     subroutine coeffs_2_4_2 ( a, q )
 
       CCTK_REAL, dimension(3), intent(OUT) :: a
@@ -581,6 +669,81 @@ module All_Coeffs_mod
 
     end subroutine coeffs_up_6_3_opt
 
+    subroutine coeffs_up_6 ( a1, q1, a2, q2 )
+
+      CCTK_REAL, dimension(-4:4), intent(OUT) :: a1, a2
+      CCTK_REAL, dimension(5,4), intent(OUT) :: q1, q2
+
+      a1(-4) = 0.01666666666666666666666666666666666666667_wp
+      a1(-3) = -0.1333333333333333333333333333333333333333_wp
+      a1(-2) = 0.5_wp
+      a1(-1) = -1.333333333333333333333333333333333333333_wp
+      a1(0) = 0.5833333333333333333333333333333333333333_wp
+      a1(1) = 0.4_wp
+      a1(2) = -0.03333333333333333333333333333333333333333_wp
+      a1(3) = zero
+      a1(4) = zero
+
+      q1(1,1) = zero
+      q1(2,1) = zero
+      q1(3,1) = zero
+      q1(4,1) = zero
+      q1(5,1) = zero
+
+      q1(1,2) = -1.0_wp
+      q1(2,2) = 1.0_wp
+      q1(3,2) = zero
+      q1(4,2) = zero
+      q1(5,2) = zero
+
+      q1(1,3) = 0.5_wp
+      q1(2,3) = -2.0_wp
+      q1(3,3) = 1.5_wp
+      q1(4,3) = zero
+      q1(5,3) = zero
+
+      q1(1,4) = -0.08333333333333333333333333333333333333333_wp
+      q1(2,4) = 0.5_wp
+      q1(3,4) = -1.5_wp
+      q1(4,4) = 0.8333333333333333333333333333333333333333_wp
+      q1(5,4) = 0.25_wp
+
+      a2(-4) = zero
+      a2(-3) = zero
+      a2(-2) = 0.03333333333333333333333333333333333333333_wp
+      a2(-1) = -0.4_wp
+      a2(0) = -0.5833333333333333333333333333333333333333_wp
+      a2(1) = 1.333333333333333333333333333333333333333_wp
+      a2(2) = -0.5_wp
+      a2(3) = 0.1333333333333333333333333333333333333333_wp
+      a2(4) = -0.01666666666666666666666666666666666666667_wp
+
+      q2(1,1) = zero
+      q2(2,1) = zero
+      q2(3,1) = zero
+      q2(4,1) = zero
+      q2(5,1) = zero
+
+      q2(1,2) = zero
+      q2(2,2) = zero
+      q2(3,2) = zero
+      q2(4,2) = -1.0_wp
+      q2(5,2) = 1.0_wp
+
+      q2(1,3) = zero
+      q2(2,3) = zero
+      q2(3,3) = -1.5_wp
+      q2(4,3) = 2.0_wp
+      q2(5,3) = -0.5_wp
+
+      q2(1,4) = -0.25_wp 
+      q2(2,4) = -0.8333333333333333333333333333333333333333_wp
+      q2(3,4) = 1.5_wp
+      q2(4,4) = -0.5_wp 
+      q2(5,4) = 0.08333333333333333333333333333333333333333_wp
+
+    end subroutine coeffs_up_6
+
     subroutine coeffs_2_6_3 ( a, q )
 
       CCTK_REAL, dimension(4), intent(OUT) :: a
@@ -1201,6 +1364,116 @@ module All_Coeffs_mod
 
     end subroutine coeffs_up_8_4_opt
 
+    subroutine coeffs_up_8 ( a1, q1, a2, q2 )
+
+      CCTK_REAL, dimension(-5:5), intent(OUT) :: a1, a2
+      CCTK_REAL, dimension(7,5), intent(OUT) :: q1, q2
+
+      a1(-5) = -0.003571428571428571428571428571428571428571_wp
+      a1(-4) = 0.03571428571428571428571428571428571428571_wp
+      a1(-3) = -0.1666666666666666666666666666666666666667_wp
+      a1(-2) = 0.5_wp
+      a1(-1) = -1.25_wp
+      a1(0) = 0.45_wp
+      a1(1) = 0.5_wp
+      a1(2) = -0.07142857142857142857142857142857142857143_wp
+      a1(3) = 0.005952380952380952380952380952380952380952_wp
+      a1(4) = zero
+      a1(5) = zero
+
+      q1(2,1) = zero
+      q1(3,1) = zero
+      q1(4,1) = zero
+      q1(5,1) = zero
+      q1(6,1) = zero
+      q1(7,1) = zero
+
+      q1(1,2) = -1.0_wp
+      q1(2,2) = 1.0_wp
+      q1(3,2) = zero
+      q1(4,2) = zero
+      q1(5,2) = zero
+      q1(6,2) = zero
+      q1(7,2) = zero
+
+      q1(1,3) = 0.5_wp
+      q1(2,3) = -2.0_wp
+      q1(3,3) = 1.5_wp
+      q1(4,3) = zero
+      q1(5,3) = zero
+      q1(6,3) = zero
+      q1(7,3) = zero
+
+      q1(1,4) = -0.08333333333333333333333333333333333333333_wp
+      q1(2,4) = 0.5_wp
+      q1(3,4) = -1.5_wp
+      q1(4,4) = 0.8333333333333333333333333333333333333333_wp
+      q1(5,4) = 0.25_wp
+      q1(6,4) = zero
+      q1(7,4) = zero
+
+      q1(1,5) = 0.01666666666666666666666666666666666666667_wp
+      q1(2,5) = -0.1333333333333333333333333333333333333333_wp
+      q1(3,5) = 0.5_wp
+      q1(4,5) = -1.333333333333333333333333333333333333333_wp
+      q1(5,5) = 0.5833333333333333333333333333333333333333_wp
+      q1(6,5) = 0.4_wp
+      q1(7,5) = -0.03333333333333333333333333333333333333333_wp
+
+      a2(-5) = zero
+      a2(-4) = zero
+      a2(-3) = -0.005952380952380952380952380952380952380952_wp
+      a2(-2) = 0.07142857142857142857142857142857142857143_wp
+      a2(-1) = -0.5_wp
+      a2(0) = -0.45_wp
+      a2(1) = 1.25_wp
+      a2(2) = -0.5_wp
+      a2(3) = 0.1666666666666666666666666666666666666667_wp
+      a2(4) = -0.03571428571428571428571428571428571428571_wp
+      a2(5) = 0.003571428571428571428571428571428571428571_wp
+
+      q2(1,1) = zero
+      q2(2,1) = zero
+      q2(3,1) = zero
+      q2(4,1) = zero
+      q2(5,1) = zero
+      q2(6,1) = zero
+      q2(7,1) = zero
+
+      q2(1,2) = zero
+      q2(2,2) = zero
+      q2(3,2) = zero
+      q2(4,2) = zero
+      q2(5,2) = zero
+      q2(6,2) = -1.0_wp
+      q2(7,2) = 1.0_wp
+
+      q2(1,3) = zero
+      q2(2,3) = zero
+      q2(3,3) = zero
+      q2(4,3) = zero
+      q2(5,3) = -1.5_wp
+      q2(6,3) = 2.0_wp
+      q2(7,3) = -0.5_wp
+
+      q2(1,4) = zero 
+      q2(2,4) = zero
+      q2(3,4) = -0.25_wp 
+      q2(4,4) = -0.8333333333333333333333333333333333333333_wp
+      q2(5,4) = 1.5_wp
+      q2(6,4) = -0.5_wp 
+      q2(7,4) = 0.08333333333333333333333333333333333333333_wp
+
+      q2(1,5) = 0.03333333333333333333333333333333333333333_wp
+      q2(2,5) = -0.4_wp
+      q2(3,5) = -0.5833333333333333333333333333333333333333_wp
+      q2(4,5) = 1.333333333333333333333333333333333333333_wp
+      q2(5,5) = -0.5_wp
+      q2(6,5) = 0.1333333333333333333333333333333333333333_wp
+      q2(7,5) = -0.01666666666666666666666666666666666666667_wp
+
+    end subroutine coeffs_up_8
+
     subroutine coeffs_2_8_4 ( a, q )
 
       CCTK_REAL, dimension(5), intent(OUT) :: a
diff --git a/src/Coefficients_2_1.F90 b/src/Coefficients_2_1.F90
index 99f21c7..3f89b83 100644
--- a/src/Coefficients_2_1.F90
+++ b/src/Coefficients_2_1.F90
@@ -136,3 +136,85 @@ subroutine set_coeff_up_2_1 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
     end if
   end do
 end subroutine set_coeff_up_2_1
+
+subroutine set_coeff_up_2 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
+
+  use All_Coeffs_mod
+
+  implicit none
+
+  DECLARE_CCTK_FUNCTIONS
+  DECLARE_CCTK_PARAMETERS
+
+  CCTK_INT, intent(IN) :: nsize, dir, loc_order
+  CCTK_INT, dimension(2), intent(IN) :: bb
+  CCTK_INT, intent(IN) :: gsize
+  CCTK_INT, dimension(nsize), intent(OUT) :: imin, imax
+  CCTK_REAL, dimension(nsize,nsize), intent(OUT) :: dd
+
+  CCTK_REAL, dimension(-2:2), save :: a1, a2
+  CCTK_REAL, dimension(3,2), save :: q1, q2
+
+  CCTK_INT :: i, il, ir
+
+  logical, save :: first = .true.
+
+  if ( .not. ( dir == -1 .or. dir == 1 ) ) then
+    call CCTK_WARN(0, 'Internal error: set_coeff_up_2 called with invalid dir value')
+  end if
+
+  if ( first ) then
+    call coeffs_up_2 ( a1, q1, a2, q2 )
+    first = .false.
+  end if
+
+  dd = zero
+  imin = 0
+  imax = -1
+
+  if ( bb(1) == 0 ) then
+    il = 1 + gsize
+  else
+    if ( dir == -1 ) then
+      dd(1,1) = q1(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:2,2) = q1(1:2,2);
+      imin(2) = 1; imax(2) = 2
+    else
+      dd(1,1) = q2(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(2:4,2) = a2(0:2);
+      imin(2) = 2; imax(2) = 4
+    end if
+    il = 3
+  end if
+  if ( bb(2) == 0 ) then
+    ir = nsize - gsize
+  else
+    if ( dir == -1 ) then
+      dd(nsize,nsize) = q1(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-3:nsize-1,nsize-1) = a1(-2:0)
+      imin(nsize-1) = nsize-3; imax(nsize-1) = nsize-1
+    else
+      dd(nsize,nsize) = q2(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-1:nsize,nsize-1) = q2(2:3,2)
+      imin(nsize-1) = nsize-1; imax(nsize-1) = nsize
+    end if
+    ir = nsize - 2
+  end if
+  do i = il, ir
+    if ( dir == -1 ) then
+      dd(i-2:i,i) = a1(-2:0)
+      imin(i) = i-2; imax(i) = i
+    else
+      dd(i:i+2,i) = a2(0:2)
+      imin(i) = i; imax(i) = i+2
+    end if
+  end do
+end subroutine set_coeff_up_2
diff --git a/src/Coefficients_4_2.F90 b/src/Coefficients_4_2.F90
index 099c5b6..d2b1b07 100644
--- a/src/Coefficients_4_2.F90
+++ b/src/Coefficients_4_2.F90
@@ -178,3 +178,97 @@ subroutine set_coeff_up_4_2 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
     end if
   end do
 end subroutine set_coeff_up_4_2
+
+subroutine set_coeff_up_4 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
+
+  use All_Coeffs_mod
+
+  implicit none
+
+  DECLARE_CCTK_FUNCTIONS
+  DECLARE_CCTK_PARAMETERS
+
+  CCTK_INT, intent(IN) :: nsize, dir, loc_order
+  CCTK_INT, dimension(2), intent(IN) :: bb
+  CCTK_INT, intent(IN) :: gsize
+  CCTK_INT, dimension(nsize), intent(OUT) :: imin, imax
+  CCTK_REAL, dimension(nsize,nsize), intent(OUT) :: dd
+
+  CCTK_REAL, dimension(-3:3), save :: a1, a2
+  CCTK_REAL, dimension(4,3), save :: q1, q2
+
+  CCTK_INT :: i, il, ir
+
+  logical, save :: first = .true.
+
+  if ( .not. ( dir == -1 .or. dir == 1 ) ) then
+    call CCTK_WARN(0, 'Internal error: set_coeff_up_2 called with invalid dir value')
+  end if
+
+  if ( first ) then
+    call coeffs_up_4 ( a1, q1, a2, q2 )
+    first = .false.
+  end if
+
+  dd = zero
+  imin = 0
+  imax = -1
+
+  if ( bb(1) == 0 ) then
+    il = 1 + gsize
+  else
+    if ( dir == -1 ) then
+      dd(1,1) = q1(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:2,2) = q1(1:2,2);
+      imin(2) = 1; imax(2) = 2
+
+      dd(1:3,3) = q1(1:3,3);
+      imin(3) = 1; imax(3) = 3
+    else
+      dd(1,1) = q2(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:5,2) = a2(-1:3);
+      imin(2) = 1; imax(2) = 5
+
+      dd(2:6,3) = a2(-1:3);
+      imin(3) = 2; imax(3) = 6
+    end if
+    il = 4
+  end if
+  if ( bb(2) == 0 ) then
+    ir = nsize - gsize
+  else
+    if ( dir == -1 ) then
+      dd(nsize,nsize) = q1(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-4:nsize,nsize-1) = a1(-3:1)
+      imin(nsize-1) = nsize-4; imax(nsize-1) = nsize
+
+      dd(nsize-5:nsize-1,nsize-2) = a1(-3:1)
+      imin(nsize-2) = nsize-5; imax(nsize-2) = nsize-1
+    else
+      dd(nsize,nsize) = q2(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-1:nsize,nsize-1) = q2(3:4,2)
+      imin(nsize-1) = nsize-1; imax(nsize-1) = nsize
+
+      dd(nsize-2:nsize,nsize-2) = q2(2:4,3)
+      imin(nsize-2) = nsize-2; imax(nsize-2) = nsize
+    end if
+    ir = nsize - 3
+  end if
+  do i = il, ir
+    if ( dir == -1 ) then
+      dd(i-3:i+1,i) = a1(-3:1)
+      imin(i) = i-3; imax(i) = i+1
+    else
+      dd(i-1:i+3,i) = a2(-1:3)
+      imin(i) = i-1; imax(i) = i+3
+    end if
+  end do
+end subroutine set_coeff_up_4
diff --git a/src/Coefficients_6_3_min_err_coeff.F90 b/src/Coefficients_6_3_min_err_coeff.F90
index c0eef1e..04d4e4f 100644
--- a/src/Coefficients_6_3_min_err_coeff.F90
+++ b/src/Coefficients_6_3_min_err_coeff.F90
@@ -215,3 +215,109 @@ subroutine set_coeff_up_6_3_opt ( nsize, dir, loc_order, bb, gsize, imin, imax,
     end if
   end do
 end subroutine set_coeff_up_6_3_opt
+
+subroutine set_coeff_up_6 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
+
+  use All_Coeffs_mod
+
+  implicit none
+
+  DECLARE_CCTK_FUNCTIONS
+  DECLARE_CCTK_PARAMETERS
+
+  CCTK_INT, intent(IN) :: nsize, dir, loc_order
+  CCTK_INT, dimension(2), intent(IN) :: bb
+  CCTK_INT, intent(IN) :: gsize
+  CCTK_INT, dimension(nsize), intent(OUT) :: imin, imax
+  CCTK_REAL, dimension(nsize,nsize), intent(OUT) :: dd
+
+  CCTK_REAL, dimension(-4:4), save :: a1, a2
+  CCTK_REAL, dimension(5,4), save :: q1, q2
+
+  CCTK_INT :: i, il, ir
+
+  logical, save :: first = .true.
+
+  if ( .not. ( dir == -1 .or. dir == 1 ) ) then
+    call CCTK_WARN(0, 'Internal error: set_coeff_up_2 called with invalid dir value')
+  end if
+
+  if ( first ) then
+    call coeffs_up_6 ( a1, q1, a2, q2 )
+    first = .false.
+  end if
+
+  dd = zero
+  imin = 0
+  imax = -1
+
+  if ( bb(1) == 0 ) then
+    il = 1 + gsize
+  else
+    if ( dir == -1 ) then
+      dd(1,1) = q1(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:2,2) = q1(1:2,2);
+      imin(2) = 1; imax(2) = 2
+
+      dd(1:3,3) = q1(1:3,3);
+      imin(3) = 1; imax(3) = 3
+
+      dd(1:5,4) = q1(1:5,4);
+      imin(4) = 1; imax(4) = 5
+    else
+      dd(1,1) = q2(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:5,2) = q2(1:5,4);
+      imin(2) = 1; imax(2) = 5
+
+      dd(1:7,3) = a2(-2:4);
+      imin(3) = 1; imax(3) = 7
+
+      dd(2:8,4) = a2(-2:4);
+      imin(4) = 2; imax(4) = 8
+    end if
+    il = 5
+  end if
+  if ( bb(2) == 0 ) then
+    ir = nsize - gsize
+  else
+    if ( dir == -1 ) then
+      dd(nsize,nsize) = q1(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-4:nsize,nsize-1) = q1(1:5,4)
+      imin(nsize-1) = nsize-4; imax(nsize-1) = nsize
+
+      dd(nsize-6:nsize,nsize-2) = a1(-4:2)
+      imin(nsize-2) = nsize-6; imax(nsize-2) = nsize
+
+      dd(nsize-7:nsize-1,nsize-3) = a1(-4:2)
+      imin(nsize-3) = nsize-7; imax(nsize-3) = nsize-1
+    else
+      dd(nsize,nsize) = q2(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-1:nsize,nsize-1) = q2(4:5,2)
+      imin(nsize-1) = nsize-1; imax(nsize-1) = nsize
+
+      dd(nsize-2:nsize,nsize-2) = q2(3:5,3)
+      imin(nsize-2) = nsize-2; imax(nsize-2) = nsize
+
+      dd(nsize-4:nsize,nsize-3) = q2(1:5,4)
+      imin(nsize-3) = nsize-4; imax(nsize-3) = nsize
+    end if
+    ir = nsize - 4
+  end if
+  do i = il, ir
+    if ( dir == -1 ) then
+      dd(i-4:i+2,i) = a1(-4:2)
+      imin(i) = i-4; imax(i) = i+2
+    else
+      dd(i-2:i+4,i) = a2(-2:4)
+      imin(i) = i-2; imax(i) = i+4
+    end if
+  end do
+end subroutine set_coeff_up_6
diff --git a/src/Coefficients_8_4_min_err_coeff.F90 b/src/Coefficients_8_4_min_err_coeff.F90
index 5bae502..8f0c5bf 100644
--- a/src/Coefficients_8_4_min_err_coeff.F90
+++ b/src/Coefficients_8_4_min_err_coeff.F90
@@ -251,3 +251,121 @@ subroutine set_coeff_up_8_4_opt ( nsize, dir, loc_order, bb, gsize, imin, imax,
     end if
   end do
 end subroutine set_coeff_up_8_4_opt
+
+subroutine set_coeff_up_8 ( nsize, dir, loc_order, bb, gsize, imin, imax, dd )
+
+  use All_Coeffs_mod
+
+  implicit none
+
+  DECLARE_CCTK_FUNCTIONS
+  DECLARE_CCTK_PARAMETERS
+
+  CCTK_INT, intent(IN) :: nsize, dir, loc_order
+  CCTK_INT, dimension(2), intent(IN) :: bb
+  CCTK_INT, intent(IN) :: gsize
+  CCTK_INT, dimension(nsize), intent(OUT) :: imin, imax
+  CCTK_REAL, dimension(nsize,nsize), intent(OUT) :: dd
+
+  CCTK_REAL, dimension(-5:5), save :: a1, a2
+  CCTK_REAL, dimension(7,5), save :: q1, q2
+
+  CCTK_INT :: i, il, ir
+
+  logical, save :: first = .true.
+
+  if ( .not. ( dir == -1 .or. dir == 1 ) ) then
+    call CCTK_WARN(0, 'Internal error: set_coeff_up_8 called with invalid dir value')
+  end if
+
+  if ( first ) then
+    call coeffs_up_8 ( a1, q1, a2, q2 )
+    first = .false.
+  end if
+
+  dd = zero
+  imin = 0
+  imax = -1
+
+  if ( bb(1) == 0 ) then
+    il = 1 + gsize
+  else
+    if ( dir == -1 ) then
+      dd(1,1) = q1(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:2,2) = q1(1:2,2);
+      imin(2) = 1; imax(2) = 2
+
+      dd(1:3,3) = q1(1:3,3);
+      imin(3) = 1; imax(3) = 3
+
+      dd(1:5,4) = q1(1:5,4);
+      imin(4) = 1; imax(4) = 5
+
+      dd(1:7,5) = q1(1:7,5);
+      imin(5) = 1; imax(5) = 7
+    else
+      dd(1,1) = q2(1,1);
+      imin(1) = 1; imax(1) = 1
+
+      dd(1:5,2) = q2(3:7,4);
+      imin(2) = 1; imax(2) = 5
+
+      dd(1:7,3) = q2(1:7,5);
+      imin(3) = 1; imax(3) = 7
+
+      dd(1:9,4) = a2(-3:5);
+      imin(4) = 1; imax(4) = 9
+
+      dd(2:10,5) = a2(-3:5);
+      imin(5) = 2; imax(5) = 10
+    end if
+    il = 6
+  end if
+  if ( bb(2) == 0 ) then
+    ir = nsize - gsize
+  else
+    if ( dir == -1 ) then
+      dd(nsize,nsize) = q1(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-4:nsize,nsize-1) = q1(1:5,4)
+      imin(nsize-1) = nsize-4; imax(nsize-1) = nsize
+
+      dd(nsize-6:nsize,nsize-2) = q1(1:7,5)
+      imin(nsize-2) = nsize-6; imax(nsize-2) = nsize
+
+      dd(nsize-8:nsize,nsize-3) = a1(-5:3)
+      imin(nsize-3) = nsize-8; imax(nsize-3) = nsize
+
+      dd(nsize-9:nsize-1,nsize-4) = a1(-5:3)
+      imin(nsize-4) = nsize-9; imax(nsize-4) = nsize-1
+    else
+      dd(nsize,nsize) = q2(1,1)
+      imin(nsize) = nsize; imax(nsize) = nsize
+
+      dd(nsize-1:nsize,nsize-1) = q2(6:7,2)
+      imin(nsize-1) = nsize-1; imax(nsize-1) = nsize
+
+      dd(nsize-2:nsize,nsize-2) = q2(5:7,3)
+      imin(nsize-2) = nsize-2; imax(nsize-2) = nsize
+
+      dd(nsize-4:nsize,nsize-3) = q2(3:7,4)
+      imin(nsize-3) = nsize-4; imax(nsize-3) = nsize
+
+      dd(nsize-6:nsize,nsize-4) = q2(1:7,5)
+      imin(nsize-4) = nsize-6; imax(nsize-4) = nsize
+    end if
+    ir = nsize - 5
+  end if
+  do i = il, ir
+    if ( dir == -1 ) then
+      dd(i-5:i+3,i) = a1(-5:3)
+      imin(i) = i-5; imax(i) = i+3
+    else
+      dd(i-3:i+5,i) = a2(-3:5)
+      imin(i) = i-3; imax(i) = i+5
+    end if
+  end do
+end subroutine set_coeff_up_8
diff --git a/src/get_up_coeffs.c b/src/get_up_coeffs.c
index 3a78f6a..e75ab9a 100644
--- a/src/get_up_coeffs.c
+++ b/src/get_up_coeffs.c
@@ -60,6 +60,39 @@ void DiffUpCoeff ( const CCTK_POINTER_TO_CONST cctkGH_, const CCTK_INT dir,
                                                CCTK_INT *imax,
                                                CCTK_REAL *q);
 
+  void CCTK_FCALL CCTK_FNAME(set_coeff_up_2)(const CCTK_INT *nsize,
+                                             const CCTK_INT *up_dir,
+                                             const CCTK_INT *loc_order,
+                                             const CCTK_INT *bb,
+                                             const CCTK_INT *gsize,
+                                             CCTK_INT *imin,
+                                             CCTK_INT *imax,
+                                             CCTK_REAL *q);
+  void CCTK_FCALL CCTK_FNAME(set_coeff_up_4)(const CCTK_INT *nsize,
+                                             const CCTK_INT *up_dir,
+                                             const CCTK_INT *loc_order,
+                                             const CCTK_INT *bb,
+                                             const CCTK_INT *gsize,
+                                             CCTK_INT *imin,
+                                             CCTK_INT *imax,
+                                             CCTK_REAL *q);
+  void CCTK_FCALL CCTK_FNAME(set_coeff_up_6)(const CCTK_INT *nsize,
+                                             const CCTK_INT *up_dir,
+                                             const CCTK_INT *loc_order,
+                                             const CCTK_INT *bb,
+                                             const CCTK_INT *gsize,
+                                             CCTK_INT *imin,
+                                             CCTK_INT *imax,
+                                             CCTK_REAL *q);
+  void CCTK_FCALL CCTK_FNAME(set_coeff_up_8)(const CCTK_INT *nsize,
+                                             const CCTK_INT *up_dir,
+                                             const CCTK_INT *loc_order,
+                                             const CCTK_INT *bb,
+                                             const CCTK_INT *gsize,
+                                             CCTK_INT *imin,
+                                             CCTK_INT *imax,
+                                             CCTK_REAL *q);
+
   ni = cctk_lsh[0]; nj = cctk_lsh[1]; nk = cctk_lsh[2];
 
   if ( table_handle >=0 ) {
@@ -99,35 +132,58 @@ void DiffUpCoeff ( const CCTK_POINTER_TO_CONST cctkGH_, const CCTK_INT dir,
   }
 
   if ( CCTK_Equals(norm_type,"Diagonal") ) {
-    if ( CCTK_Equals(operator_type,"Minimal Bandwidth") ) {
-      switch(loc_order) {
-      case 2: {
-        CCTK_FNAME(set_coeff_up_2_1)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
-        break;
-      }
-      case 4: {
-        CCTK_FNAME(set_coeff_up_4_2)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
-        break;
-      }
-      default:
-        CCTK_WARN (0, "Unknown upwinding stencil specified");
+    if ( sbp_upwind_deriv ) {
+      if ( CCTK_Equals(operator_type,"Minimal Bandwidth") ) {
+        switch(loc_order) {
+        case 2: {
+          CCTK_FNAME(set_coeff_up_2_1)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        case 4: {
+          CCTK_FNAME(set_coeff_up_4_2)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        default:
+          CCTK_WARN (0, "Unknown upwinding stencil specified");
+        }
+      } else {
+        switch(loc_order) {
+        case 2: {
+          CCTK_FNAME(set_coeff_up_2_1)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        case 4: {
+          CCTK_FNAME(set_coeff_up_4_2)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        case 6: {
+          CCTK_FNAME(set_coeff_up_6_3_opt)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        case 8: {
+          CCTK_FNAME(set_coeff_up_8_4_opt)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+          break;
+        }
+        default:
+          CCTK_WARN (0, "Unknown upwinding stencil specified");
+        }
       }
-    } else {
+   } else {
       switch(loc_order) {
       case 2: {
-        CCTK_FNAME(set_coeff_up_2_1)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+        CCTK_FNAME(set_coeff_up_2)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
         break;
       }
       case 4: {
-        CCTK_FNAME(set_coeff_up_4_2)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+        CCTK_FNAME(set_coeff_up_4)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
         break;
       }
       case 6: {
-        CCTK_FNAME(set_coeff_up_6_3_opt)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+        CCTK_FNAME(set_coeff_up_6)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
         break;
       }
       case 8: {
-        CCTK_FNAME(set_coeff_up_8_4_opt)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
+        CCTK_FNAME(set_coeff_up_8)(&nsize,&up,&loc_order,bb,&gsize,imin,imax,q);
         break;
       }
       default: