!---------------------------------------------------------------------
!     Copyright (C) GFD Dennou Club, 2004. All rights reserved.
!---------------------------------------------------------------------
                                                                 !=begin
!= Module Differentiate_Center4
!
!   * Developer: KITAMORI Taichi (kitamo@ep.sci.hokudai.ac.jp)
!   * Version: $Id: differentiate_center4.f90,v 1.2 2005/04/22 15:01:52 sugiyama Exp $ 
!   * Tag Name: $Name:  $
!   * Change History: 
!
!== Overview 
!
!微分演算を行うための関数群をまとめたパッケージ型モジュール
!水平 Arakawa-C, 鉛直 Lorentz グリッドでの微分計算
!
!== Error Handling
!
!== Known Bugs
!
!== Note
!
!本モジュールで定義されている関数は, モジュール differentiate_center2 で
!定義された関数に同じである. 
!
!== Future Plans
!
                                                                 !=end
module differentiate_center4
  
  interface ss_dx_fs
     function ss_dx_Dim1CenteredOrder4_fs(var)     
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax)
       real(8)             :: ss_dx_Dim1CenteredOrder4_fs(DimXMin:DimXMax)
     end function ss_dx_Dim1CenteredOrder4_fs

     function ss_dx_Dim2CenteredOrder4_fs(var) 
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: ss_dx_Dim2CenteredOrder4_fs(DimXMin:DimXMax, DimZMin:DimZMax)
     end function ss_dx_Dim2CenteredOrder4_fs
  end interface
  
  interface ss_dz_sf
     function ss_dz_Dim1CenteredOrder4_sf(var)     
       use gridset
       real(8), intent(in) :: var(DimZMin:DimZMax)
       real(8)             :: ss_dz_Dim1CenteredOrder4_sf(DimZMin:DimZMax)
     end function ss_dz_Dim1CenteredOrder4_sf
     
     function ss_dz_Dim2CenteredOrder4_sf(var) 
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: ss_dz_Dim2CenteredOrder4_sf(DimXMin:DimXMax, DimZMin:DimZMax)
     end function ss_dz_Dim2CenteredOrder4_sf
  end interface

  interface fs_dx_ss
     function fs_dx_Dim2CenteredOrder4_ss(var) 
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: fs_dx_Dim2CenteredOrder4_ss(DimXMin:DimXMax, DimZMin:DimZMax)
     end function fs_dx_Dim2CenteredOrder4_ss
  end interface

  interface fs_dz_ff
     function fs_dz_Dim2CenteredOrder4_ff(var)
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: fs_dz_Dim2CenteredOrder4_ff(DimXMin:DimXMax, DimZMin:DimZMax)
     end function fs_dz_Dim2CenteredOrder4_ff
  end interface

  interface sf_dz_ss
     function sf_dz_Dim2CenteredOrder4_ss(var) 
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: sf_dz_Dim2CenteredOrder4_ss(DimXMin:DimXMax, DimZMin:DimZMax)
     end function sf_dz_Dim2CenteredOrder4_ss
  end interface

  interface sf_dx_ff
     function sf_dx_Dim2CenteredOrder4_ff(var) 
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: sf_dx_Dim2CenteredOrder4_ff(DimXMin:DimXMax, DimZMin:DimZMax)
     end function sf_dx_Dim2CenteredOrder4_ff
  end interface

  interface ff_dx_sf
     function ff_dx_Dim2CenteredOrder4_sf(var)
       use gridset 
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: ff_dx_Dim2CenteredOrder4_sf(DimXMin:DimXMax, DimZMin:DimZMax)
     end function ff_dx_Dim2CenteredOrder4_sf
  end interface

  interface ff_dz_fs
     function ff_dz_Dim2CenteredOrder4_fs(var)
       use gridset
       real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
       real(8)             :: ff_dz_Dim2CenteredOrder4_fs(DimXMin:DimXMax, DimZMin:DimZMax)
     end function ff_dz_Dim2CenteredOrder4_fs
    end interface

end module differentiate_center4


function ss_dx_Dim1CenteredOrder4_fs(var)     
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax)
  real(8)             :: ss_dx_Dim1CenteredOrder4_fs(DimXMin:DimXMax)
  
  ss_dx_Dim1CenteredOrder4_fs = 0.0d0
  ss_dx_Dim1CenteredOrder4_fs(DimXMin+2:DimXMax-1) =       &
       &  (                       &
       &     var(DimXMin+2:DimXMax-1)     &
       &   - var(DimXMin+1:DimXMax-2)     &
       &  ) * 9.0d0 / ( 8.0d0 *DelX )     &
       &- (                             &
       &     var(DimXMin+3:DimXMax)     &
       &   - var(DimXMin:DimXMax-3)    &
       &  ) / (24.0d0 * DelX)
end function ss_dx_Dim1CenteredOrder4_fs


function ss_dz_Dim1CenteredOrder4_sf(var)     
  use gridset
  implicit none

  real(8), intent(in) :: var(DimZMin:DimZMax)
  real(8)             :: ss_dz_Dim1CenteredOrder4_sf(DimZMin:DimZMax)
  
  ss_dz_Dim1CenteredOrder4_sf = 0.0d0
  ss_dz_Dim1CenteredOrder4_sf(DimZMin+2:DimZMax-1) =        &
       &  (                       &
       &     var(DimZMin+2:DimZMax-1)     &
       &   - var(DimZMin+1:DimZMax-2)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelZ ) &
       &- (                              &
       &     var(DimZMin+3:DimZMax)    &
       &   - var(DimZMin:DimZMax-3)    &
       &  ) / ( 24.0d0 * DelZ ) 
end function ss_dz_Dim1CenteredOrder4_sf


function ss_dx_Dim2CenteredOrder4_fs(var) 
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: ss_dx_Dim2CenteredOrder4_fs(DimXMin:DimXMax, DimZMin:DimZMax)
  
  ss_dx_Dim2CenteredOrder4_fs = 0.0d0 
  ss_dx_Dim2CenteredOrder4_fs(DimXMin+2:DimXMax-1,DimZMin:DimZMax) =       &
       &  (                                 &
       &     var(DimXMin+2:DimXMax-1,DimZMin:DimZMax)     &
       &   - var(DimXMin+1:DimXMax-2,DimZMin:DimZMax)     &
       &  ) * 9.0d0 / (8.0d0 * DelX )                                &
       &- (                                               &
       &     var(DimXMin+3:DimXMax,DimZMin:DimZMax)     &
       &   - var(DimXMin:DimXMax-3,DimZMin:DimZMax)     &
       &  ) / ( 24.0d0 * DelX )
end function ss_dx_Dim2CenteredOrder4_fs


function ss_dz_Dim2CenteredOrder4_sf(var) 
  use gridset
  implicit none
  
  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: ss_dz_Dim2CenteredOrder4_sf(DimXMin:DimXMax, DimZMin:DimZMax)
  
  ss_dz_Dim2CenteredOrder4_sf = 0.0d0 
  ss_dz_Dim2CenteredOrder4_sf(DimXMin:DimXMax,DimZMin+2:DimZMax-1) =        &
       &  (                                 &
       &     var(DimXMin:DimXMax,DimZMin+2:DimZMax-1)     &
       &   - var(DimXMin:DimXMax,DimZMin+1:DimZMax-2)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelZ )      &
       &- (                                 &
       &     var(DimXMin:DimXMax,DimZMin+3:DimZMax)       &
       &   - var(DimXMin:DimXMax,DimZMin:DimZMax-3)       &
       &  ) / ( 24.0d0 * DelZ )
end function ss_dz_Dim2CenteredOrder4_sf



!!!
!!!速度 u の位置での微分演算
!!!
function fs_dx_Dim2CenteredOrder4_ss(var) 
  use gridset
  implicit none
  
  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: fs_dx_Dim2CenteredOrder4_ss(DimXMin:DimXMax, DimZMin:DimZMax)

  fs_dx_Dim2CenteredOrder4_ss = 0.0d0
  fs_dx_Dim2CenteredOrder4_ss(DimXMin+1:DimXMax-2,DimZMin:DimZMax) =        &
       &  (                                 &
       &     var(DimXMin+2:DimXMax-1,DimZMin:DimZMax)     &
       &   - var(DimXMin+1:DimXMax-2,DimZMin:DimZMax)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelX )      &
       &- (                                 &
       &     var(DimXMin+3:DimXMax,DimZMin:DimZMax)     &
       &   - var(DimXMin:DimXMax-3,DimZMin:DimZMax)     &
       &  ) / ( 24.0d0 * DelX)
end function fs_dx_Dim2CenteredOrder4_ss


function fs_dz_Dim2CenteredOrder4_ff(var)
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: fs_dz_Dim2CenteredOrder4_ff(DimXMin:DimXMax, DimZMin:DimZMax)

   fs_dz_Dim2CenteredOrder4_ff = 0.0d0
   fs_dz_Dim2CenteredOrder4_ff(DimXMin:DimXMax, DimZMin+2:DimZMax-1) =  &
        &  (  var(DimXMin:DimXMax,DimZMin+2:DimZMax-1)     &
        &   - var(DimXMin:DimXMax,DimZMin+1:DimZMax-2)     &
        &  ) * 9.0d0 / ( 8.0d0 * DelX )      &
        &- (                                 &
        &     var(DimXMin:DimXMax,DimZMin+3:DimZMax)     &
        &   - var(DimXMin:DimXMax,DimZMin:DimZMax-3)     &
        &  ) / ( 24.0d0 * DelX)

end function fs_dz_Dim2CenteredOrder4_ff
  
!!!
!!!速度 w の位置での微分操作
!!!
function sf_dx_Dim2CenteredOrder4_ff(var) 
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: sf_dx_Dim2CenteredOrder4_ff(DimXMin:DimXMax, DimZMin:DimZMax)

  sf_dx_Dim2CenteredOrder4_ff = 0.0d0
  sf_dx_Dim2CenteredOrder4_ff(DimXMin+2:DimXMax-1, DimZMin:DimZMax) =  &
       &  (                                 &
       &     var(DimXMin+2:DimXMax-1,DimZMin:DimZMax)     &
       &   - var(DimXMin+1:DimXMax-2,DimZMin:DimZMax)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelX )      &
       &- (                                 &
       &     var(DimXMin+3:DimXMax,DimZMin:DimZMax)     &
       &   - var(DimXMin:DimXMax-3,DimZMin:DimZMax)     &
       &  ) / ( 24.0d0 * DelX)

end function sf_dx_Dim2CenteredOrder4_ff

function sf_dz_Dim2CenteredOrder4_ss(var) 
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: sf_dz_Dim2CenteredOrder4_ss(DimXMin:DimXMax, DimZMin:DimZMax)

  sf_dz_Dim2CenteredOrder4_ss = 0.0d0 
  sf_dz_Dim2CenteredOrder4_ss(DimXMin:DimXMax,DimZMin+1:DimZMax-2) =        &
       &  (                                 &
       &     var(DimXMin:DimXMax,DimZMin+2:DimZMax-1)     &
       &   - var(DimXMin:DimXMax,DimZMin+1:DimZMax-2)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelZ )      &
       &- (                                 &
       &     var(DimXMin:DimXMax,DimZMin+3:DimZMax)     &
       &   - var(DimXMin:DimXMax,DimZMin:DimZMax-3)     &
       &  ) / ( 24.0d0 * DelZ )

end function sf_dz_Dim2CenteredOrder4_ss

!!!
!!!格子点の角の位置での微分操作
!!!
function ff_dx_Dim2CenteredOrder4_sf(var)
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: ff_dx_Dim2CenteredOrder4_sf(DimXMin:DimXMax, DimZMin:DimZMax)

  ff_dx_Dim2CenteredOrder4_sf = 0.0d0
  ff_dx_Dim2CenteredOrder4_sf(DimXMin+1:DimXMax-2, DimZMin:DimZMax) =   &
       &  (                                 &
       &     var(DimXMin+2:DimXMax-1,DimZMin:DimZMax)     &
       &   - var(DimXMin+1:DimXMax-2,DimZMin:DimZMax)     &
       &  ) * 9.0d0 / ( 8.0d0 * DelX )      &
       &- (                                 &
       &     var(DimXMin+3:DimXMax,DimZMin:DimZMax)     &
       &   - var(DimXMin:DimXMax-3,DimZMin:DimZMax)     &
       &  ) / ( 24.0d0 * DelX)

end function ff_dx_Dim2CenteredOrder4_sf


function ff_dz_Dim2CenteredOrder4_fs(var)
  use gridset
  implicit none

  real(8), intent(in) :: var(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)             :: ff_dz_Dim2CenteredOrder4_fs(DimXMin:DimXMax, DimZMin:DimZMax)

   ff_dz_Dim2CenteredOrder4_fs = 0.0d0
   ff_dz_Dim2CenteredOrder4_fs(DimXMin:DimXMax, DimZMin+1:DimZMax-2) =  &
        &  (  var(DimXMin:DimXMax,DimZMin+2:DimZMax-1)     &
        &   - var(DimXMin:DimXMax,DimZMin+1:DimZMax-2)     &
        &  ) * 9.0d0 / ( 8.0d0 * DelX )      &
        &- (                                 &
        &     var(DimXMin:DimXMax,DimZMin+3:DimZMax)     &
        &   - var(DimXMin:DimXMax,DimZMin:DimZMax-3)     &
        &  ) / ( 24.0d0 * DelX)


end function ff_dz_Dim2CenteredOrder4_fs
