
      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 tensors(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, &
                   r11,r12,r13,r22,r23,r33, &
                   nx,ny,nz,h)
use hpfprocessors
use types

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

!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
!HPF$ ALIGN WITH *g11 :: r11,r12,r13,r22,r23,r33

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

!.......................................................................
INTERFACE

subroutine inverse(g11,g12,g13,g22,g23,g33, &
                   gup11,gup12,gup13,gup22,gup23,gup33)
use hpfprocessors
use types
real(SP), dimension(:,:,:), intent(IN) :: g11,g12,g13,g22,g23,g33
real(SP), dimension(:,:,:),intent(OUT) :: gup11,gup12,gup13,gup22,gup23,gup33
!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
end subroutine inverse

subroutine connect0(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(OUT) :: c111,c112,c113,c122,c123,c133
real(SP), dimension(:,:,:),intent(OUT) :: c211,c212,c213,c222,c223,c233
real(SP), dimension(:,:,:),intent(OUT) :: 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
end subroutine connect0

subroutine riccidd(r11,r12,r13,r22,r23,r33, &
                   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(OUT) :: r11,r12,r13,r22,r23,r33
real(SP), dimension(:,:,:),intent(IN)  :: c111,c112,c113,c122,c123,c133
real(SP), dimension(:,:,:),intent(IN)  :: c211,c212,c213,c222,c223,c233
real(SP), dimension(:,:,:),intent(IN)  :: c311,c312,c313,c322,c323,c333

!HPF$ DISTRIBUTE r11 *(*,*,BLOCK) ONTO *PROC
!HPF$ ALIGN WITH *r11 :: r12,r13,r22,r23,r33
!HPF$ ALIGN WITH *r11 :: c111,c112,c113,c122,c123,c133
!HPF$ ALIGN WITH *r11 :: c211,c212,c213,c222,c223,c233
!HPF$ ALIGN WITH *r11 :: c311,c312,c313,c322,c323,c333

integer(I4B)  nx, ny, nz
real(SP) h
end subroutine riccidd
END INTERFACE
!.......................................................................


call inverse(g11,g12,g13,g22,g23,g33, &
             gup11,gup12,gup13,gup22,gup23,gup33)
call connect0(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)
call riccidd(r11,r12,r13,r22,r23,r33, &
             c111,c112,c113,c122,c123,c133, &
             c211,c212,c213,c222,c223,c233, &
             c311,c312,c313,c322,c323,c333, &
             nx,ny,nz,h)

return
end subroutine tensors









