
      MODULE types
      
!    symbolic names for 4, 2, and 1 byte integers
	
	INTEGER, PARAMETER :: I4B = SELECTED_INT_KIND(9)
	INTEGER, PARAMETER :: I2B = SELECTED_INT_KIND(4)
	INTEGER, PARAMETER :: I1B = SELECTED_INT_KIND(2)

! and for single and double precision reals
!	INTEGER, PARAMETER :: SP = KIND(1.0)
!	INTEGER, PARAMETER :: DP = KIND(1.0D0)
	INTEGER, PARAMETER :: SP = KIND(1.0D0)

      END MODULE types


	
module hpfprocessors
!HPF$ PROCESSORS PROC(4)
end module hpfprocessors


module diff_interface

INTERFACE
subroutine diff_x(nx,ny,nz,h,f,df)
use types
use hpfprocessors
integer(I4B)                :: nx,ny,nz
real(SP) h
real(SP), dimension(:,:,:), intent(IN)  :: f
real(SP), dimension(:,:,:), intent(OUT) :: df

!HPF$ DISTRIBUTE *(*,*,BLOCK) ONTO *PROC :: f
!HPF$ ALIGN WITH *f                      :: df
end subroutine diff_x

subroutine diff_y(nx,ny,nz,h,f,df)
use types
use hpfprocessors
integer(I4B)                :: nx,ny,nz
real(SP) h
real(SP), dimension(:,:,:), intent(IN)  :: f
real(SP), dimension(:,:,:), intent(OUT) :: df

!HPF$ DISTRIBUTE *(*,*,BLOCK) ONTO *PROC :: f
!HPF$ ALIGN WITH *f                      :: df
end subroutine diff_y

subroutine diff_z(nx,ny,nz,h,f,df)
use types
use hpfprocessors
integer(I4B)                :: nx,ny,nz
real(SP) h
real(SP), dimension(:,:,:), intent(IN)  :: f
real(SP), dimension(:,:,:), intent(OUT) :: df

!HPF$ DISTRIBUTE *(*,*,BLOCK) ONTO *PROC :: f
!HPF$ ALIGN WITH *f                      :: df
end subroutine diff_z

END INTERFACE
end module diff_interface







subroutine connect5(g11,g12,g13,g22,g23,g33, &
                   gup11,gup12,gup13,gup22,gup23,gup33, &
                   c111,c112,c113,c122,c123,c133, & 
                   c211,c212,c213,c222,c223,c233, & 
                   c311,c312,c313,c322,c323,c333, &
                   nx,ny,nz,h)

use types
use hpfprocessors
use diff_interface

real(SP), dimension(:,:,:), intent(IN) :: g11,g12,g13,g22,g23,g33
real(SP), dimension(:,:,:), intent(IN) :: gup11,gup12,gup13,gup22,gup23,gup33
real(SP), dimension(:,:,:), intent(INOUT):: c111,c112,c113,c122,c123,c133
real(SP), dimension(:,:,:), intent(INOUT):: c211,c212,c213,c222,c223,c233
real(SP), dimension(:,:,:), intent(INOUT):: c311,c312,c313,c322,c323,c333

!HPF$ DISTRIBUTE g11 *(*,*,BLOCK) ONTO *PROC
!HPF$ ALIGN WITH *g11 :: g12,g13,g22,g23,g33
!HPF$ ALIGN WITH *g11 :: gup11,gup12,gup13,gup22,gup23,gup33
!HPF$ ALIGN WITH *g11 :: c111,c112,c113,c122,c123,c133
!HPF$ ALIGN WITH *g11 :: c211,c212,c213,c222,c223,c233
!HPF$ ALIGN WITH *g11 :: c311,c312,c313,c322,c323,c333

integer(I4B)  nx, ny, nz
real(SP) h

real(SP), dimension(size(g11,1),size(g11,2),size(g11,3)) :: & 
     t11,t12,t13,t21,t22,t23,t31,t32,t33

!HPF$ ALIGN WITH  g11 ::  t11,t12,t13,t21,t22,t23,t31,t32,t33

      call diff_y(nx,ny,nz,h,g11,t11)
      call diff_y(nx,ny,nz,h,g12,t12)
      call diff_y(nx,ny,nz,h,g13,t13)

      call diff_y(nx,ny,nz,h,g22,t21)
      call diff_y(nx,ny,nz,h,g23,t22)
      call diff_y(nx,ny,nz,h,g33,t23)

! --------- connect8 -------------------------------------------------

      c111 = c111  - 0.5d0*gup12 *t11
      c112 = c112  - 0.5d0*gup12 *t12 
      c113 = c113  - 0.5d0*gup12 *t13 
      c122 = c122  - 0.5d0*gup12 *t21 
      c123 = c123  - 0.5d0*gup12 *t22 
      c133 = c133  - 0.5d0*gup12 *t23 

      c211 = c211  - 0.5d0*gup22 *t11 
      c212 = c212  - 0.5d0*gup22 *t12 
      c213 = c213  - 0.5d0*gup22 *t13 
      c222 = c222  - 0.5d0*gup22 *t21 
      c223 = c223  - 0.5d0*gup22 *t22 
      c233 = c233  - 0.5d0*gup22 *t23 

      c311 = c311  - 0.5d0*gup23 *t11 
      c312 = c312  - 0.5d0*gup23 *t12 
      c313 = c313  - 0.5d0*gup23 *t13 
      c322 = c322  - 0.5d0*gup23 *t21 
      c323 = c323  - 0.5d0*gup23 *t22 
      c333 = c333  - 0.5d0*gup23 *t23 


return
end
