• VTF
  • FSI
  • AMROC
  • SFC
  • Motion
  • STLIB
  • Main Page

    • src/generic/Generic_CellWallDerivatives.f90

    ! ----- THE X-CELL WALLS


! --- VALUES AT CELL WALLS
MODULE Generic_XCellWallAv

  ! Returns xcwfave - the value of f evaluated at 
  !                    cell walls.  
  !               xcwfave(i) is f evaluated at the interface
  !                    between the i, and i-1 cells.
  !                               



  INTERFACE XCellWallAv
     MODULE PROCEDURE OneDXCellWallAv
     MODULE PROCEDURE OneDXCellWallAvIdx
     MODULE PROCEDURE TwoDXCellWallAv
     MODULE PROCEDURE TwoDXCellWallAvIdx
     MODULE PROCEDURE ThreeDXCellWallAv
     MODULE PROCEDURE ThreeDXCellWallAvIdx
  END INTERFACE

CONTAINS


  SUBROUTINE OneDXCellWallAv(f,xcwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi)
    
    integer i
    ! initialize the cell wall values to be zero
  
    xcwfave(ixlo) = f(ixlo)
    
    ! use second-order interpolation
    do i = ixlo, ixhi-1
       xcwfave(i+1) = (f(i+1) + f(i))*0.5d0
    end do

    RETURN

  END SUBROUTINE OneDXCellWallAv

  SUBROUTINE OneDXCellWallAvIdx(f,xcwfave,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    integer n, id

    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi)
    
    integer i
    ! initialize the cell wall values to be zero
  
    xcwfave(ixlo) = f(id,ixlo)
    
    ! use second-order interpolation
    do i = ixlo, ixhi-1
       xcwfave(i+1) = (f(id,i+1) + f(id,i))*0.5d0
    end do

    RETURN

  END SUBROUTINE OneDXCellWallAvIdx

  SUBROUTINE TwoDXCellWallAv(f,xcwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi,iylo:iyhi)
  
    INTEGER :: i,j
  
    ! initialize the cell wall values to be zero
    do j=iylo,iyhi
       xcwfave(ixlo,j) = f(ixlo,j)
    end do
    
    ! use second-order interpolation
    
    DO j=iylo,iyhi
       do i = ixlo,ixhi-1
          xcwfave(i+1,j) = (f(i+1,j)+f(i,j))*0.5d0
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDXCellWallAv

    SUBROUTINE TwoDXCellWallAvIdx(f,xcwfave,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
  
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi,iylo:iyhi)
  
    INTEGER :: i,j
  
    ! initialize the cell wall values to be zero
    do j=iylo,iyhi
       xcwfave(ixlo,j) = f(id,ixlo,j)
    end do
    
    ! use second-order interpolation
    
    DO j=iylo,iyhi
       do i = ixlo,ixhi-1
          xcwfave(i+1,j) = (f(id,i+1,j)+f(id,i,j))*0.5d0
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDXCellWallAvIdx

  SUBROUTINE ThreeDXCellWallAv(f,xcwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallAv')
  
    ! initialize the cell wall values to be zero
    
    do k = izlo,izhi
       do j = iylo,iyhi
      
          xcwfave(ixlo,j,k) = f(ixlo,j,k)
       end do
    end do
    ! use second-order interpolation
    DO k=izlo,izhi
       DO j=iylo,iyhi
          do i= ixlo,ixhi-1
             xcwfave(i+1,j,k) = (f(i+1,j,k) + f(i,j,k))*0.5d0
          end do
       END DO
    END DO
    
    call cleslog_log_exit('ThreeDXCellWallAv')

    RETURN

  END SUBROUTINE ThreeDXCellWallAv

  SUBROUTINE ThreeDXCellWallAvIdx(f,xcwfave, n, id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallAvIdx')
  
    ! initialize the cell wall values to be zero
    
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwfave(ixlo,j,k) = f(id, ixlo,j,k)
       end do
    end do
    ! use second-order interpolation
    DO k=izlo,izhi
       DO j=iylo,iyhi
          do i= ixlo,ixhi-1
             xcwfave(i+1,j,k) = (f(id,i+1,j,k) + f(id,i,j,k))*0.5d0
          end do
       END DO
    END DO
    
    call cleslog_log_exit('ThreeDXCellWallAvIdx')

    RETURN

  END SUBROUTINE ThreeDXCellWallAvIdx

END MODULE Generic_XCellWallAv

! DERIVATIVES ON CELL WALLS

MODULE Generic_XCellWallDx
  ! Returns xcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               xcwdf(i) is df/dx evaluated at the interface
  !                 between the i, and i-1 cells.
  !                               



  INTERFACE XCellWallDx
     MODULE PROCEDURE OneDXCellWallDx
     MODULE PROCEDURE OneDXCellWallDxIdx
     MODULE PROCEDURE TwoDXCellWallDx
     MODULE PROCEDURE TwoDXCellWallDxIdx
     MODULE PROCEDURE ThreeDXCellWallDx
     MODULE PROCEDURE ThreeDXCellWallDxIdx
  END INTERFACE

CONTAINS

  SUBROUTINE OneDXCellWallDx(f,xcwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi)
    DOUBLE PRECISION :: cdx
  
    integer i

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
  
    xcwdf(ixlo) = 0.0d0
  
    ! use second-order center difference
    do i = ixlo,ixhi -1
       xcwdf(i+1) = (f(i+1) - f(i))*cdx
    end do
    
    RETURN
    
    
  END SUBROUTINE OneDXCellWallDx

  SUBROUTINE OneDXCellWallDxIdx(f,xcwdf,n,id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi)
    DOUBLE PRECISION :: cdx
  
    integer i

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
  
    xcwdf(ixlo) = 0.0d0
  
    ! use second-order center difference
    do i = ixlo,ixhi -1
       xcwdf(i+1) = (f(id,i+1) - f(id,i))*cdx
    end do
    
    RETURN
    
    
  END SUBROUTINE OneDXCellWallDxIdx


  SUBROUTINE TwoDXCellWallDx(f,xcwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
    
    do j = iylo,iyhi
       xcwdf(ixlo,j) = 0.0d0
    end do

    ! use second-order center difference
    DO j=iylo,iyhi
       do i = ixlo, ixhi-1
          xcwdf(i+1,j) = (f(i+1,j) - f(i,j))*cdx
       end do
    END DO
    
    RETURN
    
  END SUBROUTINE TwoDXCellWallDx

  SUBROUTINE TwoDXCellWallDxIdx(f,xcwdf,n,id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
    
    do j = iylo,iyhi
       xcwdf(ixlo,j) = 0.0d0
    end do

    ! use second-order center difference
    DO j=iylo,iyhi
       do i = ixlo, ixhi-1
          xcwdf(i+1,j) = (f(id,i+1,j) - f(id,i,j))*cdx
       end do
    END DO
    
    RETURN
    
  END SUBROUTINE TwoDXCellWallDxIdx

  SUBROUTINE ThreeDXCellWallDx(f,xcwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallDx')

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    DO k = izlo, izhi
       DO j = iylo,iyhi
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(i+1,j,k) - f(i,j,k))*cdx
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDXCellWallDx')

    RETURN
    
  END SUBROUTINE ThreeDXCellWallDx

  SUBROUTINE ThreeDXCellWallDxIdx(f,xcwdf, n, id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallDxIdx')

    cdx = 1.0d0/dx

    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    DO k = izlo, izhi
       DO j = iylo,iyhi
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(id,i+1,j,k) - f(id,i,j,k))*cdx
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDXCellWallDxIdx')

    RETURN
    
  END SUBROUTINE ThreeDXCellWallDxIdx

END MODULE Generic_XCellWallDx


MODULE Generic_XCellWallDy

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns xcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               xcwdf(i,j) is df/dy evaluated at the interface
  !                 between the (i,j), and (i-1,j) cells.
  !                               



  INTERFACE XCellWallDy
     MODULE PROCEDURE TwoDXCellWallDy
     MODULE PROCEDURE TwoDXCellWallDyIdx
     MODULE PROCEDURE ThreeDXCellWallDy
     MODULE PROCEDURE ThreeDXCellWallDyIdx
  END INTERFACE

CONTAINS



  SUBROUTINE TwoDXCellWallDy(f,xcwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j
  
    cdy = 0.25d0/dy
    
    ! initialize the cell wall derivative to be zero
    do j = iylo,iyhi
       xcwdf(ixlo,j) = 0.0d0
    end do
    do i = ixlo,ixhi
       xcwdf(i,iylo) = 0.0d0
       xcwdf(i,iyhi) = 0.0d0
    end do

    ! use second-order center difference
    DO j = iylo+1,iyhi-1
       do i = ixlo, ixhi-1 
          xcwdf(i+1,j) = (f(i+1,j+1)+f(i,j+1)-f(i+1,j-1)-f(i,j-1))*cdy
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDXCellWallDy

  SUBROUTINE TwoDXCellWallDyIdx(f,xcwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j
  
    cdy = 0.25d0/dy
    
    ! initialize the cell wall derivative to be zero
    do j = iylo,iyhi
       xcwdf(ixlo,j) = 0.0d0
    end do
    do i = ixlo,ixhi
       xcwdf(i,iylo) = 0.0d0
       xcwdf(i,iyhi) = 0.0d0
    end do

    ! use second-order center difference
    DO j = iylo+1,iyhi-1
       do i = ixlo, ixhi-1 
          xcwdf(i+1,j) = (f(id,i+1,j+1)+f(id,i,j+1)-f(id,i+1,j-1)-f(id,i,j-1))*cdy
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDXCellWallDyIdx

  SUBROUTINE ThreeDXCellWallDy(f,xcwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER ::i, j,k

    call cleslog_log_enter('ThreeDXCellWallDy')
  
    cdy = 0.25d0/dy
    
    ! initialize the cell wall derivative to be zero
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
       do i= ixlo,ixhi
          xcwdf(i,iylo,k) = 0.0d0
          xcwdf(i,iyhi,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k= izlo, izhi
       DO j = iylo+1,iyhi-1
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(i+1,j+1,k)+f(i,j+1,k)&
                  &-f(i+1,j-1,k)-f(i,j-1,k))*cdy
          end do
       end do
    end do

    call cleslog_log_exit('ThreeDXCellWallDy')

    RETURN

  END SUBROUTINE ThreeDXCellWallDy

  SUBROUTINE ThreeDXCellWallDyIdx(f,xcwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER ::i,j,k

    call cleslog_log_enter('ThreeDXCellWallDyIdx')
  
    cdy = 0.25d0/dy
    
    ! initialize the cell wall derivative to be zero
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
       do i= ixlo,ixhi
          xcwdf(i,iylo,k) = 0.0d0
          xcwdf(i,iyhi,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k= izlo, izhi
       DO j = iylo+1,iyhi-1
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(id,i+1,j+1,k)+f(id,i,j+1,k)&
                  &-f(id,i+1,j-1,k)-f(id,i,j-1,k))*cdy
          end do
       end do
    end do

    call cleslog_log_exit('ThreeDXCellWallDyIdx')

    RETURN

  END SUBROUTINE ThreeDXCellWallDyIdx

END MODULE Generic_XCellWallDy



MODULE Generic_XCellWallDz

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns xcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               xcwdf(i,j,k) is df/dz evaluated at the interface
  !                 between the (i,j,k), and (i-1,j,k) cells.
  !                               



  INTERFACE XCellWallDz
     MODULE PROCEDURE ThreeDXCellWallDz
     MODULE PROCEDURE ThreeDXCellWallDzIdx
  END INTERFACE

CONTAINS


  SUBROUTINE ThreeDXCellWallDz(f,xcwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallDz')
  
    cdz = 0.25d0/dz
    
    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
    end do
    do j = iylo,iyhi
       do i = ixlo,ixhi
          xcwdf(i,j,izlo) = 0.0d0
          xcwdf(i,j,izhi) = 0.0d0
       end do
    end do

    ! use second-order center difference
    do k= izlo+1, izhi-1
       DO j = iylo,iyhi
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(i+1,j,k+1)+f(i,j,k+1) &
                  &-f(i+1,j,k-1)-f(i,j,k-1))*cdz
          end do
       END DO
    end do
    
    call cleslog_log_exit('ThreeDXCellWallDz')

    RETURN

  END SUBROUTINE ThreeDXCellWallDz

  SUBROUTINE ThreeDXCellWallDzIdx(f,xcwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: xcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDXCellWallDzIdx')
  
    cdz = 0.25d0/dz
    
    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do j = iylo,iyhi
          xcwdf(ixlo,j,k) = 0.0d0
       end do
    end do
    do j = iylo,iyhi
       do i = ixlo,ixhi
          xcwdf(i,j,izlo) = 0.0d0
          xcwdf(i,j,izhi) = 0.0d0
       end do
    end do

    ! use second-order center difference
    do k= izlo+1, izhi-1
       DO j = iylo,iyhi
          do i = ixlo,ixhi-1
             xcwdf(i+1,j,k) = (f(id,i+1,j,k+1)+f(id,i,j,k+1) &
                  &-f(id,i+1,j,k-1)-f(id,i,j,k-1))*cdz
          end do
       END DO
    end do
    
    call cleslog_log_exit('ThreeDXCellWallDzIdx')

    RETURN

  END SUBROUTINE ThreeDXCellWallDzIdx

END MODULE Generic_XCellWallDz

! ----- THE Y-CELL WALLS


MODULE Generic_YCellWallAv

  ! Returns ycwfave - the value of f evaluated at 
  !                    cell walls.  
  !               ycwfave(i,j,k) is f evaluated at the interface
  !                    between the (i,j,k), and (i,j-1,k) cells.
  !                               


  
  INTERFACE YCellWallAv
     MODULE PROCEDURE TwoDYCellWallAv
     MODULE PROCEDURE TwoDYCellWallAvIdx
     MODULE PROCEDURE ThreeDYCellWallAv
     MODULE PROCEDURE ThreeDYCellWallAvIdx
  END INTERFACE

CONTAINS
  

  SUBROUTINE TwoDYCellWallAv(f,ycwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwfave(ixlo:ixhi,iylo:iyhi)
  
    INTEGER :: i,j
  
    ! initialize the cell wall values to be zero
    
    do i = ixlo, ixhi
       ycwfave(i,iylo) = f(i,iylo)
    end do

    ! use second-order interpolation
    DO j = iylo,iyhi-1
       do i = ixlo, ixhi
          ycwfave(i,j+1) = (f(i,j+1) + f(i,j))*0.5d0
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDYCellWallAv

  SUBROUTINE TwoDYCellWallAvIdx(f,ycwfave,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwfave(ixlo:ixhi,iylo:iyhi)
  
    INTEGER :: i,j
  
    ! initialize the cell wall values to be zero
    
    do i = ixlo, ixhi
       ycwfave(i,iylo) = f(id,i,iylo)
    end do

    ! use second-order interpolation
    DO j = iylo,iyhi-1
       do i = ixlo, ixhi
          ycwfave(i,j+1) = (f(id,i,j+1) + f(id,i,j))*0.5d0
       end do
    END DO
    
    RETURN

  END SUBROUTINE TwoDYCellWallAvIdx

  SUBROUTINE ThreeDYCellWallAv(f,ycwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallAv')
  
    ! initialize the cell wall values to be zero
  
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwfave(i,iylo,k) = f(i,iylo,k)
       end do
    end do
    
    ! use second-order interpolation
    DO k=izlo,izhi
       DO j=iylo,iyhi-1
          do i = ixlo,ixhi
             ycwfave(i,j+1,k) = (f(i,j+1,k) + f(i,j,k))*0.5d0
          end do
       END DO
       
    END DO

    call cleslog_log_exit('ThreeDYCellWallAv')

    RETURN

  END SUBROUTINE ThreeDYCellWallAv

  SUBROUTINE ThreeDYCellWallAvIdx(f,ycwfave,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallAvIdx')
  
    ! initialize the cell wall values to be zero
  
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwfave(i,iylo,k) = f(id,i,iylo,k)
       end do
    end do
    
    ! use second-order interpolation
    DO k=izlo,izhi
       DO j=iylo,iyhi-1
          do i = ixlo,ixhi
             ycwfave(i,j+1,k) = (f(id,i,j+1,k) + f(id,i,j,k))*0.5d0
          end do
       END DO
       
    END DO

    call cleslog_log_exit('ThreeDYCellWallAvIdx')

    RETURN

  END SUBROUTINE ThreeDYCellWallAvIdx

END MODULE Generic_YCellWallAv

! DERIVATIVES ON CELL WALLS

MODULE Generic_YCellWallDx
  ! Returns ycwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               ycwdf(i) is df/dx evaluated at the interface
  !                 between the j, and j-1 cells.
  !                               



  INTERFACE YCellWallDx
     MODULE PROCEDURE TwoDYCellWallDx
     MODULE PROCEDURE TwoDYCellWallDxIdx
     MODULE PROCEDURE ThreeDYCellWallDx
     MODULE PROCEDURE ThreeDYCellWallDxIdx
  END INTERFACE

CONTAINS

  SUBROUTINE TwoDYCellWallDx(f,ycwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
    
    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j

    cdx = 0.25d0/dx
    
    ! initialize the cell wall derivative to be zero
    
    do i = ixlo,ixhi
       ycwdf(i,iylo) = 0.0d0
    end do
    do j = iylo,iyhi
       ycwdf(ixlo,j) = 0.0d0
       ycwdf(ixhi,j) = 0.0d0
    end do
    ! use second-order center difference
    
    DO j = iylo,iyhi-1
       do i = ixlo+1, ixhi-1
          ycwdf(i,j+1) = (f(i+1,j+1)+f(i+1,j) &
               &-f(i-1,j+1)-f(i-1,j))*cdx
       enddo
    end do
    RETURN
    
  END SUBROUTINE TwoDYCellWallDx

    SUBROUTINE TwoDYCellWallDxIdx(f,ycwdf,n,id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
    
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j

    cdx = 0.25d0/dx
    
    ! initialize the cell wall derivative to be zero
    
    do i = ixlo,ixhi
       ycwdf(i,iylo) = 0.0d0
    end do
    do j = iylo,iyhi
       ycwdf(ixlo,j) = 0.0d0
       ycwdf(ixhi,j) = 0.0d0
    end do
    ! use second-order center difference
    
    DO j = iylo,iyhi-1
       do i = ixlo+1, ixhi-1
          ycwdf(i,j+1) = (f(id,i+1,j+1)+f(id,i+1,j) &
               &-f(id,i-1,j+1)-f(id,i-1,j))*cdx
       enddo
    end do
    RETURN
    
  END SUBROUTINE TwoDYCellWallDxIdx

  SUBROUTINE ThreeDYCellWallDx(f,ycwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms

    ! ---- Shared Procedures
    
    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDx')

    cdx = 0.25d0/dx

    ! initialize the cell wall derivative to be zero
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do

       do j = iylo,iyhi
          ycwdf(ixlo,j,k) = 0.0d0
          ycwdf(ixhi,j,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k=izlo, izhi
       DO j = iylo,iyhi-1
          do i = ixlo+1,ixhi-1
             ycwdf(i,j+1,k) = (f(i+1,j+1,k)+f(i+1,j,k) &
                  &-f(i-1,j+1,k)-f(i-1,j,k))*cdx
          end do
       enddo
    enddo
    
    call cleslog_log_exit('ThreeDYCellWallDx')

    RETURN
    
  END SUBROUTINE ThreeDYCellWallDx

    SUBROUTINE ThreeDYCellWallDxIdx(f,ycwdf,n,id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms

    ! ---- Shared Procedures
    
    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDxIdx')

    cdx = 0.25d0/dx

    ! initialize the cell wall derivative to be zero
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do

       do j = iylo,iyhi
          ycwdf(ixlo,j,k) = 0.0d0
          ycwdf(ixhi,j,k) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k=izlo, izhi
       DO j = iylo,iyhi-1
          do i = ixlo+1,ixhi-1
             ycwdf(i,j+1,k) = (f(id,i+1,j+1,k)+f(id,i+1,j,k) &
                  &-f(id,i-1,j+1,k)-f(id,i-1,j,k))*cdx
          end do
       enddo
    enddo
    
    call cleslog_log_exit('ThreeDYCellWallDxIdx')

    RETURN
    
  END SUBROUTINE ThreeDYCellWallDxIdx

END MODULE Generic_YCellWallDx


MODULE Generic_YCellWallDy

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns ycwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               ycwdf(i,j) is df/dy evaluated at the interface
  !                 between the (i,j,k), and (i,j-1,k) cells.
  !                               



  INTERFACE YCellWallDy
     MODULE PROCEDURE TwoDYCellWallDy
     MODULE PROCEDURE TwoDYCellWallDyIdx
     MODULE PROCEDURE ThreeDYCellWallDy
     MODULE PROCEDURE ThreeDYCellWallDyIdx
  END INTERFACE

CONTAINS

  SUBROUTINE TwoDYCellWallDy(f,ycwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi)
   
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j
  
    cdy = 1.0d0/dy
    
    ! initialize the cell wall derivative to be zero
    
    do i = ixlo,ixhi
       ycwdf(i,iylo) = 0.0d0
    end do
    ! use second-order center difference
    ! to compute the derivative at the wall between j,j-1
    DO j = iylo,iyhi-1
       do i = ixlo,ixhi
          ycwdf(i,j+1) =(f(i,j+1)-f(i,j))*cdy
       end do
    END DO

    RETURN

  END SUBROUTINE TwoDYCellWallDy

  SUBROUTINE TwoDYCellWallDyIdx(f,ycwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi)
   
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j
  
    cdy = 1.0d0/dy
    
    ! initialize the cell wall derivative to be zero
    
    do i = ixlo,ixhi
       ycwdf(i,iylo) = 0.0d0
    end do
    ! use second-order center difference
    ! to compute the derivative at the wall between j,j-1
    DO j = iylo,iyhi-1
       do i = ixlo,ixhi
          ycwdf(i,j+1) =(f(id,i,j+1)-f(id,i,j))*cdy
       end do
    END DO

    RETURN

  END SUBROUTINE TwoDYCellWallDyIdx

  SUBROUTINE ThreeDYCellWallDy(f,ycwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDy')
  
    cdy = 1.0d0/dy
    
    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do
    end do
    
    ! use second-order center difference
    ! to compute the derivative at the cell wall between j,j-1
    DO k = izlo,izhi
       DO j = iylo,iyhi-1
          do i = ixlo,ixhi
             ycwdf(i,j+1,k) =(f(i,j+1,k)-f(i,j,k))*cdy
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDYCellWallDy')
    
    RETURN
    
  END SUBROUTINE ThreeDYCellWallDy

  SUBROUTINE ThreeDYCellWallDyIdx(f,ycwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDyIdx')
  
    cdy = 1.0d0/dy
    
    ! initialize the cell wall derivative to be zero
  
    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do
    end do
    
    ! use second-order center difference
    ! to compute the derivative at the cell wall between j,j-1
    DO k = izlo,izhi
       DO j = iylo,iyhi-1
          do i = ixlo,ixhi
             ycwdf(i,j+1,k) =(f(id,i,j+1,k)-f(id,i,j,k))*cdy
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDYCellWallDyIdx')
    
    RETURN
    
  END SUBROUTINE ThreeDYCellWallDyIdx

END MODULE Generic_YCellWallDy

MODULE Generic_YCellWallDz

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns xcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               ycwdf(i,j,k) is df/dz evaluated at the interface
  !                 between the (i,j,k), and (i,j-1,k) cells.
  !                               



  INTERFACE YCellWallDz
     MODULE PROCEDURE ThreeDYCellWallDz
     MODULE PROCEDURE ThreeDYCellWallDzIdx
  END INTERFACE

CONTAINS


  SUBROUTINE ThreeDYCellWallDz(f,ycwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
    
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDz')
  
    cdz = 0.25d0/dz

    ! initialize the cell wall derivative to be zero

    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do
    end do

    do j = iylo,iyhi
       do i = ixlo,ixhi
          ycwdf(i,j,izlo) = 0.0d0
          ycwdf(i,j,izhi) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k=izlo+1, izhi-1
       DO j = iylo,iyhi-1
          do i = ixlo, ixhi
             ycwdf(i,j+1,k) = (f(i,j+1,k+1)+f(i,j,k+1) &
                  &-f(i,j+1,k-1)-f(i,j,k-1))*cdz
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDYCellWallDz')
    
    RETURN

  END SUBROUTINE ThreeDYCellWallDz

    SUBROUTINE ThreeDYCellWallDzIdx(f,ycwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures
    
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: ycwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDYCellWallDzIdx')
  
    cdz = 0.25d0/dz

    ! initialize the cell wall derivative to be zero

    do k = izlo,izhi
       do i = ixlo,ixhi
          ycwdf(i,iylo,k) = 0.0d0
       end do
    end do

    do j = iylo,iyhi
       do i = ixlo,ixhi
          ycwdf(i,j,izlo) = 0.0d0
          ycwdf(i,j,izhi) = 0.0d0
       end do
    end do
    ! use second-order center difference
    
    do k=izlo+1, izhi-1
       DO j = iylo,iyhi-1
          do i = ixlo, ixhi
             ycwdf(i,j+1,k) = (f(id,i,j+1,k+1)+f(id,i,j,k+1) &
                  &-f(id,i,j+1,k-1)-f(id,i,j,k-1))*cdz
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDYCellWallDzIdx')
    
    RETURN

  END SUBROUTINE ThreeDYCellWallDzIdx

END MODULE Generic_YCellWallDz

! ----- THE Z-CELL WALLS


MODULE Generic_ZCellWallAv

  ! Returns zcwfave - the value of f evaluated at 
  !                    cell walls.  
  !               zcwfave(i,j,k) is f evaluated at the interface
  !                    between the (i,j,k), and (i,j,k-1) cells.
  !                               


  
  INTERFACE ZCellWallAv
     MODULE PROCEDURE ThreeDZCellWallAv
     MODULE PROCEDURE ThreeDZCellWallAvIdx
  END INTERFACE

CONTAINS
  
  SUBROUTINE ThreeDZCellWallAv(f,zcwfave)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallAv')
  
    ! initialize the cell wall values to be zero
    do j = iylo,iyhi
       do i = ixlo,ixhi
          zcwfave(i,j,izlo) = f(i,j,izlo)
       end do
    end do
    
    ! use second-order interpolation
    DO k=izlo,izhi-1
       DO j=iylo,iyhi
          do i=ixlo,ixhi
             zcwfave(i,j,k+1) = (f(i,j,k+1) + f(i,j,k))*0.5d0
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDZCellWallAv')

    RETURN

  END SUBROUTINE ThreeDZCellWallAv

  SUBROUTINE ThreeDZCellWallAvIdx(f,zcwfave,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwfave(ixlo:ixhi,iylo:iyhi,izlo:izhi)
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallAvIdx')
  
    ! initialize the cell wall values to be zero
    do j = iylo,iyhi
       do i = ixlo,ixhi
          zcwfave(i,j,izlo) = f(id,i,j,izlo)
       end do
    end do
    
    ! use second-order interpolation
    DO k=izlo,izhi-1
       DO j=iylo,iyhi
          do i=ixlo,ixhi
             zcwfave(i,j,k+1) = (f(id,i,j,k+1) + f(id,i,j,k))*0.5d0
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDZCellWallAvIdx')

    RETURN

  END SUBROUTINE ThreeDZCellWallAvIdx

END MODULE Generic_ZCellWallAv


! DERIVATIVES ON CELL WALLS

MODULE Generic_ZCellWallDx
  ! Returns zcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               zcwdf(i) is df/dx evaluated at the interface
  !                 between the k, and k-1 cells.
  !                               



  INTERFACE ZCellWallDx
     MODULE PROCEDURE ThreeDZCellWallDx
     MODULE PROCEDURE ThreeDZCellWallDxIdx
  END INTERFACE

CONTAINS

  SUBROUTINE ThreeDZCellWallDx(f,zcwdf)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures

    IMPLICIT NONE

    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDx')

    cdx = 0.25d0/dx

    ! initialize the cell wall derivative to be zero
  
    do j = iylo, iyhi
       do i = ixlo,ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
   
    do k = izlo,izhi
       do j = iylo,iyhi
          zcwdf(ixlo,j,k) = 0.0d0
          zcwdf(ixhi,j,k) = 0.0d0
       end do
    end do
    
    ! use second-order center difference
    
    do k=izlo, izhi-1
       DO j = iylo,iyhi
          do i = ixlo+1,ixhi-1
             zcwdf(i,j,k+1) = (f(i+1,j,k+1)+f(i+1,j,k)&
                  &-f(i-1,j,k+1)-f(i-1,j,k))*cdx
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDZCellWallDx')
    
    RETURN
    
  END SUBROUTINE ThreeDZCellWallDx

  SUBROUTINE ThreeDZCellWallDxIdx(f,zcwdf,n,id)


    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
    ! ---- Shared Procedures

    IMPLICIT NONE

    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdx
  
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDxIdx')

    cdx = 0.25d0/dx

    ! initialize the cell wall derivative to be zero
  
    do j = iylo, iyhi
       do i = ixlo,ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
   
    do k = izlo,izhi
       do j = iylo,iyhi
          zcwdf(ixlo,j,k) = 0.0d0
          zcwdf(ixhi,j,k) = 0.0d0
       end do
    end do
    
    ! use second-order center difference
    
    do k=izlo, izhi-1
       DO j = iylo,iyhi
          do i = ixlo+1,ixhi-1
             zcwdf(i,j,k+1) = (f(id,i+1,j,k+1)+f(id,i+1,j,k)&
                  &-f(id,i-1,j,k+1)-f(id,i-1,j,k))*cdx
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDZCellWallDxIdx')
    
    RETURN
    
  END SUBROUTINE ThreeDZCellWallDxIdx

END MODULE Generic_ZCellWallDx


MODULE Generic_ZCellWallDy

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns ycwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               zcwdf(i,j,k) is df/dy evaluated at the interface
  !                 between the (i,j,k), and (i,j,k-1) cells.
  !                               



  INTERFACE ZCellWallDy
     MODULE PROCEDURE ThreeDZCellWallDy
     MODULE PROCEDURE ThreeDZCellWallDyIdx
  END INTERFACE

CONTAINS



  SUBROUTINE ThreeDZCellWallDy(f,zcwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms

    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDy')
  
    cdy = 0.25d0/dy
    
    do k = izlo,izhi
       do i = ixlo,ixhi
          zcwdf(i,iylo,k) = 0.0d0
          zcwdf(i,iyhi,k) = 0.0d0
       end do
    end do
    do j = iylo,iyhi
       do i = ixlo, ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
 
    ! use second-order center difference
    
    do k=izlo, izhi-1
       DO j = iylo+1,iyhi-1
          do i = ixlo,ixhi
             zcwdf(i,j,k+1) = (f(i,j+1,k+1)+f(i,j+1,k) &
                  &-f(i,j-1,k+1)-f(i,j-1,k))*cdy
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDZCellWallDy')

    RETURN
    
  END SUBROUTINE ThreeDZCellWallDy

  SUBROUTINE ThreeDZCellWallDyIdx(f,zcwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms

    ! ---- Shared Procedures
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION :: cdy
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDyIdx')
  
    cdy = 0.25d0/dy
    
    do k = izlo,izhi
       do i = ixlo,ixhi
          zcwdf(i,iylo,k) = 0.0d0
          zcwdf(i,iyhi,k) = 0.0d0
       end do
    end do
    do j = iylo,iyhi
       do i = ixlo, ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
 
    ! use second-order center difference
    
    do k=izlo, izhi-1
       DO j = iylo+1,iyhi-1
          do i = ixlo,ixhi
             zcwdf(i,j,k+1) = (f(id,i,j+1,k+1)+f(id,i,j+1,k) &
                  &-f(id,i,j-1,k+1)-f(id,i,j-1,k))*cdy
          end do
       enddo
    enddo

    call cleslog_log_exit('ThreeDZCellWallDyIdx')

    RETURN
    
  END SUBROUTINE ThreeDZCellWallDyIdx

END MODULE Generic_ZCellWallDy

MODULE Generic_ZCellWallDz

  ! ---     | j+1 |
  ! ---   __|_____|__
  ! ---     |  i  |
  ! ---  i-1|  j  |i+1
  ! ---   __|_____|__
  ! ---     |     |
  ! ---     | j-1 |

  ! Returns zcwdf - the derivative of f evaluated at 
  !                 cell walls.  
  !               zcwdf(i,j,k) is df/dz evaluated at the interface
  !                 between the (i,j,k), and (i,j,k-1) cells.
  !                               



  INTERFACE ZCellWallDz
     MODULE PROCEDURE ThreeDZCellWallDz
     MODULE PROCEDURE ThreeDZCellWallDzIdx
  END INTERFACE

CONTAINS


  SUBROUTINE ThreeDZCellWallDz(f,zcwdf)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    DOUBLE PRECISION, INTENT(IN) :: f(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDz')
  
    cdz = 1.0d0/dz
    
    ! initialize the cell wall derivative to be zero
    do j = iylo,iyhi
       do i = ixlo, ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
    ! use second-order center difference
    DO k = izlo,izhi-1
       DO j = iylo,iyhi
          do i = ixlo,ixhi
             zcwdf(i,j,k+1) =(f(i,j,k+1)-f(i,j,k))*cdz
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDZCellWallDz')

    RETURN
    
  END SUBROUTINE ThreeDZCellWallDz

  SUBROUTINE ThreeDZCellWallDzIdx(f,zcwdf,n,id)

    ! ----  Shared Variables
    USE mesh
    USE array_bounds
    USE method_parms
    
  
    IMPLICIT NONE
    
    INTEGER :: n,id
    DOUBLE PRECISION, INTENT(IN) :: f(n,ixlo:ixhi,iylo:iyhi,izlo:izhi)
    DOUBLE PRECISION, INTENT(OUT):: zcwdf(ixlo:ixhi,iylo:iyhi,izlo:izhi)
    
    DOUBLE PRECISION :: cdz
    
    INTEGER :: i,j,k

    call cleslog_log_enter('ThreeDZCellWallDzIdx')
  
    cdz = 1.0d0/dz
    
    ! initialize the cell wall derivative to be zero
    do j = iylo,iyhi
       do i = ixlo, ixhi
          zcwdf(i,j,izlo) = 0.0d0
       end do
    end do
    ! use second-order center difference
    DO k = izlo,izhi-1
       DO j = iylo,iyhi
          do i = ixlo,ixhi
             zcwdf(i,j,k+1) =(f(id,i,j,k+1)-f(id,i,j,k))*cdz
          end do
       END DO
    END DO

    call cleslog_log_exit('ThreeDZCellWallDzIdx')

    RETURN
    
  END SUBROUTINE ThreeDZCellWallDzIdx

END MODULE Generic_ZCellWallDz

<