program write_unstrd
!
! This program sets up a 3-D unstructured grid,
! initializes cell-centered and vertex-centered arrays on
! the grid, and then calls routines in the Rocketeer API to write 
! out 2 HDF files, one with geometry data (unstr0.hdf) and one 
! (unstr1.hdf) that refers to the geometry data in unstr0.hdf.
!
! The grid came from the Silo User's Guide (LLNL).
!
! Written by Robert Fiedler, 12/13/2002.
!
      USE mod_hdf_writer_API

      implicit NONE
!
      character*(10) filename          ! HDF field data file
      character*(10) filegeom          ! HDF geometry data file
      character*(32) long_name
      integer ne,nn                    ! Number of elements, nodes
      parameter(ne=3,nn=16)
      REAL(MYREAL) xyz(nn,3)                   ! Cartesian coordinates
      integer conn(ne,8)               ! Connectivity
      REAL(MYREAL) dens(nn)                    ! Node centered scalar
      REAL(MYREAL) dsp(nn,3)                   ! nodal displacement
      REAL(MYREAL) s(nn,9)                     ! nodal tensor
      CHARACTER*(8) time                       ! Physical problem time
!
      integer i,l,indx
      REAL(MYREAL) x,y,z,r
      REAL(MYREAL) ax,ay,az,dens0,densm
      data ax,ay,az,dens0,densm /1.0,2.0,3.0,1.0,100.0/
!
! Set up unstructured grid from Silo User's Guide, Fig. 2.2.
!
      data xyz /0.,1.,2.,0.,1.,2.,0.,1.,2.,0.,1.,2.,3.,3.,3.,3.,  &
                1.,1.,1.,0.,0.,0.,1.,1.,1.,0.,0.,0.,.5,.5,.5,.5,  &
                1.,1.,1.,1.,1.,1.,0.,0.,0.,0.,0.,0.,.5,.5,.5,.5/
      data conn/ 1, 2, 3,  &
                 4, 5, 6,  &
                10,11,12,  &
                 7, 8, 9,  &
                 2, 3,13,  &
                 5, 6,13,  &
                11,12,13,  &
                 8, 9,13/

! Loop over the time index

      do indx=0,1

        write(filename,'("unstr",i1,".hdf")') indx
        filegeom = 'unstr0.hdf'
        write(time,'(F4.2)') real(indx)

! Make up some vector, scalar, and tensor data at nodes

        do i=1,nn
          dsp(i,1) = indx*xyz(i,1)*(-0.3)
          dsp(i,2) = indx*xyz(i,2)*( 0.1)
          dsp(i,3) = indx*xyz(i,3)*( 0.2)
          r = (xyz(i,1)/ax)**2 + (xyz(i,2)/ay)**2 + (xyz(i,3)/az)**2
          dens(i) = densm * exp(-r) * exp(-real(indx)) + dens0
          s(i,1) = 1. + indx*xyz(i,1)
          s(i,2) = 0.
          s(i,3) = 0.
          s(i,4) = 0.
          s(i,5) = 2. + indx*xyz(i,2)
          s(i,6) = 0.
          s(i,7) = 0.
          s(i,8) = 0.
          s(i,9) = 3. + indx*xyz(i,3)
        end do

! Call routines in the Rocketeer API to create a valid HDF
! file for Rocketeer.

! Write block header

        call w_block_header  &
         (filename, filegeom, 'block_001', time, 'solid', 3, 0,  &
          0, .false.)
!
! Write geometry data
!
! Note that if the geometry file name does not match the data
! file name, routine w_geometry_unstr does nothing.
!
        call w_geometry_unstr  &
         (filename, filegeom, 'block_001', 'solid', 'm',  &
          'x', 'y', 'z', nn, 0, ne, 8, conn, xyz(:,1), xyz(:,2), xyz(:,3))
!
! Write scalar variable (node centered)
!
        long_name = 'Density at time '//time

        call w_scalar_unstr  &
         (filename, long_name, 'd', 'kg/m^3',  &
          nn, 0, dens)
!
! Write vector variable (node centered)
!
        long_name = 'Displacement at time '//time

        call w_vector_unstr  &
         (filename, long_name, 'disp', 'm',  &
          nn, 0, dsp(:,1), dsp(:,2), dsp(:,3))
!
! Write tensor variable (cell centered)
!
        long_name = 'Stress at time '//time

        call w_tensor_unstr  &
         (filename, long_name, 'stress', 'mmm',  &
          nn, 0, s(:,1), s(:,2), S(:,3), s(:,4), s(:,5), s(:,6),  &
                 s(:,7), s(:,8), s(:,9))

       end do

       end