Add new parameters derivatives_zero_y and derivatives_zero_z in order to

ecplictly setting any derivatives in the y- and/or z-direction to zero
when simulating 1- or 2-D problems with a 3D code. This avoids issues with
errors at roundoff level in the y- and z-directions. Should be used with care
and together with the parameter check_grid_sizes to avoid complaints about
too small grids in the y- and z-directions.

Has only been tested in a few cases. If problems are caused by this commit,
please report them immediately.


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

1 files changed, 114 insertions, 110 deletions

diff --git a/src/Dissipation_4_2_alt.F90 b/src/Dissipation_4_2_alt.F90
index c0a68ba..884f6b3 100644
--- a/src/Dissipation_4_2_alt.F90
+++ b/src/Dissipation_4_2_alt.F90
@@ -82,119 +82,123 @@ subroutine dissipation_4_2_alt (var, ni, nj, nk, bb, gsize, delta, epsilon, rhs)
                                 ( var(il-3:ir-3,:,:) + &
                                   var(il+3:ir+3,:,:) ) ) * idel
 
-  call set_coeff ( a )
-
-  idel = epsilon / ( 64 * delta(2) )
-
-  if ( bb(3) == 0 ) then
-    jl = 1 + gsize(2)
-  else
-    rhs(:,1,:) = rhs(:,1,:) + &
-                 ( a(1,1) * var(:,1,:) + a(2,1) * var(:,2,:) + &
-                   a(3,1) * var(:,3,:) + a(4,1) * var(:,4,:) ) * idel
-    rhs(:,2,:) = rhs(:,2,:) + &
-                 ( a(1,2) * var(:,1,:) + a(2,2) * var(:,2,:) + &
-                   a(3,2) * var(:,3,:) + a(4,2) * var(:,4,:) + &
-                   a(5,2) * var(:,5,:) ) * idel
-    rhs(:,3,:) = rhs(:,3,:) + &
-                 ( a(1,3) * var(:,1,:) + a(2,3) * var(:,2,:) + &
-                   a(3,3) * var(:,3,:) + a(4,3) * var(:,4,:) + &
-                   a(5,3) * var(:,5,:) + a(6,3) * var(:,6,:) ) * idel
-    rhs(:,4,:) = rhs(:,4,:) + &
-                 ( a(1,4) * var(:,1,:) + a(2,4) * var(:,2,:) + &
-                   a(3,4) * var(:,3,:) + a(4,4) * var(:,4,:) + &
-                   a(5,4) * var(:,5,:) + a(6,4) * var(:,6,:) + &
-                   a(7,4) * var(:,7,:) ) * idel
-
-    jl = 5
-  end if
-  if ( bb(4) == 0 ) then
-    jr = nj - gsize(2)
-  else
-    rhs(:,nj-3,:) = rhs(:,nj-3,:) + &
-                 ( a(1,4) * var(:,nj,:) + a(2,4) * var(:,nj-1,:) + &
-                   a(3,4) * var(:,nj-2,:) + a(4,4) * var(:,nj-3,:) + &
-                   a(5,4) * var(:,nj-4,:) + a(6,4) * var(:,nj-5,:) + &
-                   a(7,4) * var(:,nj-6,:) ) * idel
-    rhs(:,nj-2,:) = rhs(:,nj-2,:) + &
-                 ( a(1,3) * var(:,nj,:) + a(2,3) * var(:,nj-1,:) + &
-                   a(3,3) * var(:,nj-2,:) + a(4,3) * var(:,nj-3,:) + &
-                   a(5,3) * var(:,nj-4,:) + a(6,3) * var(:,nj-5,:) ) * idel
-    rhs(:,nj-1,:) = rhs(:,nj-1,:) + &
-                 ( a(1,2) * var(:,nj,:) + a(2,2) * var(:,nj-1,:) + &
-                   a(3,2) * var(:,nj-2,:) + a(4,2) * var(:,nj-3,:) + &
-                   a(5,2) * var(:,nj-4,:) ) * idel
-    rhs(:,nj,:)   = rhs(:,nj,:) + &
-                 ( a(1,1) * var(:,nj,:) + a(2,1) * var(:,nj-1,:) + &
-                   a(3,1) * var(:,nj-2,:) + a(4,1) * var(:,nj-3,:) ) * idel
-
-    jr = nj - 4
+  if ( zero_derivs_y == 0 ) then
+    call set_coeff ( a )
+  
+    idel = epsilon / ( 64 * delta(2) )
+  
+    if ( bb(3) == 0 ) then
+      jl = 1 + gsize(2)
+    else
+      rhs(:,1,:) = rhs(:,1,:) + &
+                   ( a(1,1) * var(:,1,:) + a(2,1) * var(:,2,:) + &
+                     a(3,1) * var(:,3,:) + a(4,1) * var(:,4,:) ) * idel
+      rhs(:,2,:) = rhs(:,2,:) + &
+                   ( a(1,2) * var(:,1,:) + a(2,2) * var(:,2,:) + &
+                     a(3,2) * var(:,3,:) + a(4,2) * var(:,4,:) + &
+                     a(5,2) * var(:,5,:) ) * idel
+      rhs(:,3,:) = rhs(:,3,:) + &
+                   ( a(1,3) * var(:,1,:) + a(2,3) * var(:,2,:) + &
+                     a(3,3) * var(:,3,:) + a(4,3) * var(:,4,:) + &
+                     a(5,3) * var(:,5,:) + a(6,3) * var(:,6,:) ) * idel
+      rhs(:,4,:) = rhs(:,4,:) + &
+                   ( a(1,4) * var(:,1,:) + a(2,4) * var(:,2,:) + &
+                     a(3,4) * var(:,3,:) + a(4,4) * var(:,4,:) + &
+                     a(5,4) * var(:,5,:) + a(6,4) * var(:,6,:) + &
+                     a(7,4) * var(:,7,:) ) * idel
+  
+      jl = 5
+    end if
+    if ( bb(4) == 0 ) then
+      jr = nj - gsize(2)
+    else
+      rhs(:,nj-3,:) = rhs(:,nj-3,:) + &
+                   ( a(1,4) * var(:,nj,:) + a(2,4) * var(:,nj-1,:) + &
+                     a(3,4) * var(:,nj-2,:) + a(4,4) * var(:,nj-3,:) + &
+                     a(5,4) * var(:,nj-4,:) + a(6,4) * var(:,nj-5,:) + &
+                     a(7,4) * var(:,nj-6,:) ) * idel
+      rhs(:,nj-2,:) = rhs(:,nj-2,:) + &
+                   ( a(1,3) * var(:,nj,:) + a(2,3) * var(:,nj-1,:) + &
+                     a(3,3) * var(:,nj-2,:) + a(4,3) * var(:,nj-3,:) + &
+                     a(5,3) * var(:,nj-4,:) + a(6,3) * var(:,nj-5,:) ) * idel
+      rhs(:,nj-1,:) = rhs(:,nj-1,:) + &
+                   ( a(1,2) * var(:,nj,:) + a(2,2) * var(:,nj-1,:) + &
+                     a(3,2) * var(:,nj-2,:) + a(4,2) * var(:,nj-3,:) + &
+                     a(5,2) * var(:,nj-4,:) ) * idel
+      rhs(:,nj,:)   = rhs(:,nj,:) + &
+                   ( a(1,1) * var(:,nj,:) + a(2,1) * var(:,nj-1,:) + &
+                     a(3,1) * var(:,nj-2,:) + a(4,1) * var(:,nj-3,:) ) * idel
+  
+      jr = nj - 4
+    end if
+    rhs(:,jl:jr,:) = rhs(:,jl:jr,:) + &
+                     ( -20.0_wp * var(:,jl:jr,:) + &
+                        15.0_wp * ( var(:,jl-1:jr-1,:) + &
+                                    var(:,jl+1:jr+1,:) ) - &
+                         6.0_wp * ( var(:,jl-2:jr-2,:) + &
+                                    var(:,jl+2:jr+2,:) ) + &
+                                  ( var(:,jl-3:jr-3,:) + &
+                                    var(:,jl+3:jr+3,:) ) ) * idel
   end if
-  rhs(:,jl:jr,:) = rhs(:,jl:jr,:) + &
-                   ( -20.0_wp * var(:,jl:jr,:) + &
-                      15.0_wp * ( var(:,jl-1:jr-1,:) + &
-                                  var(:,jl+1:jr+1,:) ) - &
-                       6.0_wp * ( var(:,jl-2:jr-2,:) + &
-                                  var(:,jl+2:jr+2,:) ) + &
-                                ( var(:,jl-3:jr-3,:) + &
-                                  var(:,jl+3:jr+3,:) ) ) * idel
 
-  call set_coeff ( a )
-
-  idel = epsilon / ( 64 * delta(3) )
-
-  if ( bb(5) == 0 ) then
-    kl = 1 + gsize(3)
-  else
-    rhs(:,:,1) = rhs(:,:,1) + &
-                 ( a(1,1) * var(:,:,1) + a(2,1) * var(:,:,2) + &
-                   a(3,1) * var(:,:,3) + a(4,1) * var(:,:,4) ) * idel
-    rhs(:,:,2) = rhs(:,:,2) + &
-                 ( a(1,2) * var(:,:,1) + a(2,2) * var(:,:,2) + &
-                   a(3,2) * var(:,:,3) + a(4,2) * var(:,:,4) + &
-                   a(5,2) * var(:,:,5) ) * idel
-    rhs(:,:,3) = rhs(:,:,3) + &
-                 ( a(1,3) * var(:,:,1) + a(2,3) * var(:,:,2) + &
-                   a(3,3) * var(:,:,3) + a(4,3) * var(:,:,4) + &
-                   a(5,3) * var(:,:,5) + a(6,3) * var(:,:,6) ) * idel
-    rhs(:,:,4) = rhs(:,:,4) + &
-                 ( a(1,4) * var(:,:,1) + a(2,4) * var(:,:,2) + &
-                   a(3,4) * var(:,:,3) + a(4,4) * var(:,:,4) + &
-                   a(5,4) * var(:,:,5) + a(6,4) * var(:,:,6) + &
-                   a(7,4) * var(:,:,7) ) * idel
-
-    kl = 5
-  end if
-  if ( bb(6) == 0 ) then
-    kr = nk - gsize(3)
-  else
-    rhs(:,:,nk-3) = rhs(:,:,nk-3) + &
-                 ( a(1,4) * var(:,:,nk) + a(2,4) * var(:,:,nk-1) + &
-                   a(3,4) * var(:,:,nk-2) + a(4,4) * var(:,:,nk-3) + &
-                   a(5,4) * var(:,:,nk-4) + a(6,4) * var(:,:,nk-5) + &
-                   a(7,4) * var(:,:,nk-6) ) * idel
-    rhs(:,:,nk-2) = rhs(:,:,nk-2) + &
-                 ( a(1,3) * var(:,:,nk) + a(2,3) * var(:,:,nk-1) + &
-                   a(3,3) * var(:,:,nk-2) + a(4,3) * var(:,:,nk-3) + &
-                   a(5,3) * var(:,:,nk-4) + a(6,3) * var(:,:,nk-5) ) * idel
-    rhs(:,:,nk-1) = rhs(:,:,nk-1) + &
-                 ( a(1,2) * var(:,:,nk) + a(2,2) * var(:,:,nk-1) + &
-                   a(3,2) * var(:,:,nk-2) + a(4,2) * var(:,:,nk-3) + &
-                   a(5,2) * var(:,:,nk-4) ) * idel
-    rhs(:,:,nk)   = rhs(:,:,nk) + &
-                 ( a(1,1) * var(:,:,nk) + a(2,1) * var(:,:,nk-1) + &
-                   a(3,1) * var(:,:,nk-2) + a(4,1) * var(:,:,nk-3) ) * idel
-
-    kr = nk - 4
+  if ( zero_derivs_z == 0 ) then
+    call set_coeff ( a )
+  
+    idel = epsilon / ( 64 * delta(3) )
+  
+    if ( bb(5) == 0 ) then
+      kl = 1 + gsize(3)
+    else
+      rhs(:,:,1) = rhs(:,:,1) + &
+                   ( a(1,1) * var(:,:,1) + a(2,1) * var(:,:,2) + &
+                     a(3,1) * var(:,:,3) + a(4,1) * var(:,:,4) ) * idel
+      rhs(:,:,2) = rhs(:,:,2) + &
+                   ( a(1,2) * var(:,:,1) + a(2,2) * var(:,:,2) + &
+                     a(3,2) * var(:,:,3) + a(4,2) * var(:,:,4) + &
+                     a(5,2) * var(:,:,5) ) * idel
+      rhs(:,:,3) = rhs(:,:,3) + &
+                   ( a(1,3) * var(:,:,1) + a(2,3) * var(:,:,2) + &
+                     a(3,3) * var(:,:,3) + a(4,3) * var(:,:,4) + &
+                     a(5,3) * var(:,:,5) + a(6,3) * var(:,:,6) ) * idel
+      rhs(:,:,4) = rhs(:,:,4) + &
+                   ( a(1,4) * var(:,:,1) + a(2,4) * var(:,:,2) + &
+                     a(3,4) * var(:,:,3) + a(4,4) * var(:,:,4) + &
+                     a(5,4) * var(:,:,5) + a(6,4) * var(:,:,6) + &
+                     a(7,4) * var(:,:,7) ) * idel
+  
+      kl = 5
+    end if
+    if ( bb(6) == 0 ) then
+      kr = nk - gsize(3)
+    else
+      rhs(:,:,nk-3) = rhs(:,:,nk-3) + &
+                   ( a(1,4) * var(:,:,nk) + a(2,4) * var(:,:,nk-1) + &
+                     a(3,4) * var(:,:,nk-2) + a(4,4) * var(:,:,nk-3) + &
+                     a(5,4) * var(:,:,nk-4) + a(6,4) * var(:,:,nk-5) + &
+                     a(7,4) * var(:,:,nk-6) ) * idel
+      rhs(:,:,nk-2) = rhs(:,:,nk-2) + &
+                   ( a(1,3) * var(:,:,nk) + a(2,3) * var(:,:,nk-1) + &
+                     a(3,3) * var(:,:,nk-2) + a(4,3) * var(:,:,nk-3) + &
+                     a(5,3) * var(:,:,nk-4) + a(6,3) * var(:,:,nk-5) ) * idel
+      rhs(:,:,nk-1) = rhs(:,:,nk-1) + &
+                   ( a(1,2) * var(:,:,nk) + a(2,2) * var(:,:,nk-1) + &
+                     a(3,2) * var(:,:,nk-2) + a(4,2) * var(:,:,nk-3) + &
+                     a(5,2) * var(:,:,nk-4) ) * idel
+      rhs(:,:,nk)   = rhs(:,:,nk) + &
+                   ( a(1,1) * var(:,:,nk) + a(2,1) * var(:,:,nk-1) + &
+                     a(3,1) * var(:,:,nk-2) + a(4,1) * var(:,:,nk-3) ) * idel
+  
+      kr = nk - 4
+    end if
+    rhs(:,:,kl:kr) = rhs(:,:,kl:kr) + &
+                     ( -20.0_wp * var(:,:,kl:kr) + &
+                        15.0_wp * ( var(:,:,kl-1:kr-1) + &
+                                    var(:,:,kl+1:kr+1) ) - &
+                         6.0_wp * ( var(:,:,kl-2:kr-2) + &
+                                    var(:,:,kl+2:kr+2) ) + &
+                                  ( var(:,:,kl-3:kr-3) + &
+                                    var(:,:,kl+3:kr+3) ) ) * idel
   end if
-  rhs(:,:,kl:kr) = rhs(:,:,kl:kr) + &
-                   ( -20.0_wp * var(:,:,kl:kr) + &
-                      15.0_wp * ( var(:,:,kl-1:kr-1) + &
-                                  var(:,:,kl+1:kr+1) ) - &
-                       6.0_wp * ( var(:,:,kl-2:kr-2) + &
-                                  var(:,:,kl+2:kr+2) ) + &
-                                ( var(:,:,kl-3:kr-3) + &
-                                  var(:,:,kl+3:kr+3) ) ) * idel
 
 contains