anchor.diff:        theta_sp = -40._r8/300._r8 
         theta_eq =  30._r8/300._r8 
         theta_np = -20._r8/300._r8 
  -      if(testcase==7)then
  +      if(testcase==4)then
           u00 = 80.0
           lat0 = pi/7.0
           lat1 = pi/2.0 - lat0
         else
  -        if(testcase==8)then ! Jet in Southern Hemisphere
  +        if(testcase==5)then ! Jet in Southern Hemisphere
              u00 = 80.0
              lat0 = -5.d0*deg2rad
              lat1 = -45.d0*deg2rad
           end if
         end if
  -	
  +
         en = exp(-4/(lat1 - lat0)**2)
         umen = u00/en
         totvol = 0.0D0
  
positive.code:           call plot_scalarfield(bt, mesh)
          end if
       
          if(testcase==2 .or. testcase==3)then
            h_error%name=trim(swmname)//"_h_error_t"//trim(adjustl(trim(atime)))
            call plot_scalarfield(h_error, mesh)
  
  
negative.code:     real(r8),allocatable  :: p3(:)      
      real(r8),allocatable  :: vn(:)      
  
      !Reconstruction vector
      real (r8) :: urecon(1:3)
  
      allocate (p1(3),p2(3),p3(3),vn(3))
      FEPS = 1.0D-8
      dot = 0.0D0
  