!
!  generate_grid.f90
!
      subroutine generate_grid(origin,x,y,z,lx,ly,lz,ni,nj,nk,rz,iz,dz,cz,ierr)        
!                                                                                
      implicit none                                                  
      real :: rz(*)
      integer :: iz(*)
      character (len=*) :: cz(*)
      real (selected_int_kind (5)) :: dz(*)
!
      integer :: ni,nj,nk,ierr                                                                                   
      real :: lx,ly,lz,x(ni,nj,nk),y(ni,nj,nk),z(ni,nj,nk),origin(3) 
      real :: xs,ys,zs                                              
!
      real :: dltx,dlty,dltz 
      integer :: nzp1,nyp1,nxp1,nzp2,nyp2,nxp2,i,j,k   
!                                                                                                                                                                                                                   
      dltx = lx/float(ni-1) 
      dlty = ly/float(nj-1)
      dltz = lz/float(nk-1)
!
      xs = origin(1) ; ys = origin(2) ; zs = origin(3)
      do k=1,nk                                                           
        do j=1,nj                                                           
          do i=1,ni                                                                           
            x(i,j,k) = xs + dltx*float(i-1)                                         
            y(i,j,k) = ys + dlty*float(j-1)                                        
            z(i,j,k) = zs + dltz*float(k-1)                                               
          enddo
        enddo
      enddo                                                                                                                                
!                                                                                                                                                                                              
      end subroutine generate_grid 