fileio_gnuplot.f90

Go to the documentation of this file.

!#############################################################################
!#                                                                           #
!# fosite - 3D hydrodynamical simulation program                             #
!# module: fileio_gnuplot.f90                                                #
!#                                                                           #
!# Copyright (C) 2008-2017                                                   #
!# Tobias Illenseer <tillense@astrophysik.uni-kiel.de>                       #
!# Björn Sperling   <sperling@astrophysik.uni-kiel.de>                       #
!# Jannes Klee      <jklee@astrophysik.uni-kiel.de>                          #
!#                                                                           #
!# This program is free software; you can redistribute it and/or modify      #
!# it under the terms of the GNU General Public License as published by      #
!# the Free Software Foundation; either version 2 of the License, or (at     #
!# your option) any later version.                                           #
!#                                                                           #
!# This program is distributed in the hope that it will be useful, but       #
!# WITHOUT ANY WARRANTY; without even the implied warranty of                #
!# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, GOOD TITLE or        #
!# NON INFRINGEMENT.  See the GNU General Public License for more            #
!# details.                                                                  #
!#                                                                           #
!# You should have received a copy of the GNU General Public License         #
!# along with this program; if not, write to the Free Software               #
!# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.                 #
!#                                                                           #
!#############################################################################
!----------------------------------------------------------------------------!
!----------------------------------------------------------------------------!
MODULE fileio_gnuplot_mod
  USE fileio_base_mod
  USE geometry_base_mod
  USE mesh_base_mod
  USE physics_base_mod
  USE timedisc_base_mod
  USE fluxes_base_mod
  USE sources_base_mod
  USE common_dict
  IMPLICIT NONE
  !--------------------------------------------------------------------------!
  PRIVATE
  ! exclude interface block from doxygen processing
!  INTERFACE Error
!     MODULE PROCEDURE Error_gnuplot, Error_common
!  END INTERFACE
  ! Private Attributes section starts here:
  INTEGER, PARAMETER      :: HLEN = 10000
  INTEGER, PARAMETER      :: DEFAULT_DECS = 5
  CHARACTER, PARAMETER    :: SP = achar(32)       
  CHARACTER, PARAMETER    :: LF = achar(10)       
  CHARACTER*2, PARAMETER  :: RECSEP = sp // sp    
  CHARACTER*2, PARAMETER  :: LINSEP = sp // lf    
  CHARACTER*2, PARAMETER  :: BLKSEP = lf // lf    
  CHARACTER(LEN=30), PARAMETER :: &
         header_string = "# Data output of fosite" // linsep
  CHARACTER(LEN=HLEN)  :: header_buf
  !--------------------------------------------------------------------------!
  TYPE, EXTENDS(fileio_base) :: fileio_gnuplot

  CONTAINS
    ! methods
    PROCEDURE :: initfileio_gnuplot
    PROCEDURE :: finalize
    PROCEDURE :: writeheader
!    PROCEDURE :: ReadHeader
!    PROCEDURE :: WriteTimestamp
!    PROCEDURE :: ReadTimestamp
    PROCEDURE :: writedataset
    PROCEDURE :: readdataset
    PROCEDURE :: openfile
    PROCEDURE :: closefile
    PROCEDURE :: adjusttimestep
    PROCEDURE :: inctime
    PROCEDURE :: getfilename
    PROCEDURE :: getfilestatus
  END TYPE

  !--------------------------------------------------------------------------!
  PUBLIC :: &
       ! types
       fileio_typ, &
       ! constants
       readonly, readend, replace, append, &
       ascii, bin, &
       file_exists
 !--------------------------------------------------------------------------!

CONTAINS

  SUBROUTINE initfileio(this,Mesh,Physics,fmt,fmtname,fpath,filename,extension, &
       stoptime,dtwall,count,fcycles,sepfiles,unit)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base)  :: this
    CLASS(mesh_base)    :: Mesh
    CLASS(physics_base) :: Physics
    INTEGER             :: fmt
    CHARACTER(LEN=*)    :: fmtname
    CHARACTER(LEN=*)    :: fpath
    CHARACTER(LEN=*)    :: filename
    CHARACTER(LEN=*)    :: extension
    REAL                :: stoptime
    INTEGER             :: dtwall
    INTEGER             :: count
    INTEGER             :: fcycles
    LOGICAL             :: sepfiles
    INTEGER, OPTIONAL   :: unit
    !------------------------------------------------------------------------!
    INTENT(IN)    :: mesh,physics,fmt,fmtname,fpath,filename,extension,stoptime, &
         dtwall,count,fcycles,sepfiles,unit
    INTENT(INOUT) :: this
    !------------------------------------------------------------------------!
    ! basic FileIO initialization
    CALL this%InitFileIO_common(fmt,fmtname,fpath,filename,extension,fcycles,sepfiles,unit)

    this%stoptime = stoptime
    this%dtwall   = dtwall
    this%time     = 0.
    this%count    = count
    this%step     = 0
    ! count the number of output columns, i.e. fields per data point

#ifdef PARALLEL
    ! check data type extents in files
    ! first try to create a new dummy file
    CALL mpi_file_open(mpi_comm_world,getfilename(this),ior(ior(mpi_mode_rdwr,&
         mpi_mode_create),ior(mpi_mode_excl,mpi_mode_delete_on_close)),&
         mpi_info_null,this%handle,this%error)
    ! maybe file exists
    IF (this%error.NE.0) CALL mpi_file_open(mpi_comm_world,getfilename(this),&
         mpi_mode_rdonly,mpi_info_null,this%handle,this%error)
    ! then check the data type sizes
    ! extent of integer in file
    IF (this%error.EQ.0) CALL mpi_file_get_type_extent(this%handle,mpi_integer,&
         this%intext,this%error)
    ! extent of real in file
    IF (this%error.EQ.0) CALL mpi_file_get_type_extent(this%handle,default_mpi_real,&
         this%realext,this%error)
    IF (this%error.NE.0) CALL error(this,"InitFileIO","unable to check file properties")
    CALL mpi_file_close(this%handle,this%error)
#endif
  END SUBROUTINE initfileio


  SUBROUTINE initfileio_gnuplot(this,Mesh,Physics,IO,fmt,fpath,filename,stoptime,dtwall,&
       count,fcycles,unit,config)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base)              :: this
    CLASS(mesh_base)                :: Mesh
    CLASS(physics_base)             :: Physics
    TYPE(Dict_TYP),POINTER          :: IO
    TYPE(Dict_TYP),POINTER,OPTIONAL :: config
    INTEGER                         :: fmt
    CHARACTER(LEN=*)                :: fpath
    CHARACTER(LEN=*)                :: filename
    REAL                            :: stoptime
    INTEGER                         :: dtwall
    INTEGER                         :: count
    INTEGER                         :: fcycles
    INTEGER                         :: unit
    !------------------------------------------------------------------------!
    TYPE(Dict_TYP), POINTER                   :: node
    REAL, DIMENSION(:,:), POINTER             :: dummy2
    REAL, DIMENSION(:,:,:), POINTER           :: dummy3
    INTEGER                                   :: cartcoords
    INTEGER                                   :: depth
    INTEGER                                   :: err
#ifdef PARALLEL
    INTEGER                                   :: i
    INTEGER, DIMENSION(Mesh%IMAX-Mesh%IMIN+1) :: blocklen,indices
#endif
    !------------------------------------------------------------------------!
    INTENT(IN)    :: mesh,physics,fmt,fpath,filename,stoptime,dtwall,count,fcycles,&
                     unit
    INTENT(INOUT) :: this
    !------------------------------------------------------------------------!
    CALL this%InitFileIO(mesh,physics,fmt,"GNUPLOT",fpath,filename,"dat",stoptime,&
         dtwall,count,fcycles,.false.,unit)

    CALL getattr(config, "/datafile/decimals", this%DECS, default_decs)
    ! compute length of character field for real number output
    ! and check if linebuffer is large enough
    ! flen = 1 (sign) + 1 (one digit) + 1 (decimal point) + decs (decimal places)
    !      + 1 (E character) + 1 (sign of exponent) + 2 (exponent digits) + 2 (spaces)
    this%FLEN = this%DECS + 9
    this%maxcols = len(this%linebuf)/this%FLEN-1

    ALLOCATE(this%output(this%maxcols),stat=err)
    IF (this%error.NE.0) &
       CALL this%Error("InitFileIO_gnuplot","memory allocation failed for this%output")

    ! check if cartesian coordinates are selected for output;
    ! default: curvilinear coordinates (0)
    CALL getattr(config, "/datafile/cartcoords", cartcoords, 0)
    depth = 1
    node => config
    CALL writeheaderstring(header_buf,node,depth)
    IF (cartcoords.EQ.0) THEN
       CALL getattr(io,"/mesh/bary_curv",dummy3)
    ELSE
       CALL getattr(io,"/mesh/bary_centers",dummy3)
    END IF
    ! pointer to sub-array: set lower bounds
    ! use special function for bound remapping
    IF (mesh%INUM.EQ.1) THEN
      ! use format string as temp
      WRITE (this%fmtstr,'(A5,I2,A1)')'(A1,A',this%FLEN-3,')'
      WRITE(this%linebuf,trim(this%fmtstr))'#','y'
      this%output(1)%val => dummy3(:,:,2)
      this%COLS = 1
    ELSE IF (mesh%JNUM.EQ.1) THEN
      WRITE (this%fmtstr,'(A5,I2,A1)')'(A1,A',this%FLEN-3,')'
      WRITE(this%linebuf,trim(this%fmtstr))'#','x'
      this%output(1)%val => dummy3(:,:,1)
      this%COLS = 1
    ELSE
      WRITE (this%fmtstr,'(A5,I2,A2,I2,A1)')'(A1,A',this%FLEN-3,',A',this%FLEN-1,')'
      WRITE(this%linebuf,trim(this%fmtstr))'#','x','y'
      this%output(1)%val => dummy3(:,:,1)
      this%output(2)%val => dummy3(:,:,2)
      this%COLS = 2
    END IF

    ! set output-pointer and count the number of output columns
    WRITE (this%fmtstr,'(A,I2,A1)')'(A',this%FLEN-1,')'
    node => io
    CALL this%GetOutputPointer(mesh,node,this%COLS)

    ! length of one output line
    this%linelen = this%COLS * this%FLEN
    IF (this%linelen.GT.len(this%linebuf)) &
       CALL this%Error("InitFileIO_gnuplot", &
          "linebuffer to small; reducing decimals or number of output fields may help")

    header_buf = trim(header_string) // trim(header_buf) // this%linebuf(1:this%linelen)

    ! local domain size
    this%inum = mesh%IMAX - mesh%IMIN + 1
    this%jnum = mesh%JMAX - mesh%JMIN + 1

!#ifdef PARALLEL
!    ! create new data type handle for one line
!    CALL MPI_Type_contiguous(this%linelen,MPI_CHARACTER,this%basictype,this%error)
!    CALL MPI_Type_commit(this%basictype,this%error)
!
!    ! number of output blocks
!    this%blocknum = Mesh%IMAX - Mesh%IMIN + 1
!    ! size of the output buffer
!    this%bufsize  = Mesh%JMAX - Mesh%JMIN + 1
!
!    ! allocate memory for output buffer and displacement records
!    ALLOCATE(this%outbuf(this%linelen,Mesh%JMIN:Mesh%JMAX), &
!         STAT=err)
!    IF (this%error.NE.0) THEN
!       CALL Error(this,"InitFileIO_gnuplot","memory allocation failed for this%outbuf")
!    END IF
!
!    blocklen(:) = this%bufsize
!    DO i=Mesh%IMIN,Mesh%IMAX
!        indices(i-Mesh%IMIN+1) = (i-1)*Mesh%JNUM + Mesh%JMIN - 1
!    END DO
!
!    ! new file type for the staggered data
!    CALL MPI_Type_indexed(this%blocknum,blocklen,indices, &
!         this%basictype,this%filetype,this%error)
!    CALL MPI_Type_commit(this%filetype,this%error)
!#endif
    ! write the format string for one entry in the data file:
    ! FLEN-2 characters for the number and 2 for the separators
    WRITE (this%fmtstr,'(A3,I2,A,I2.2,A5)') '(ES', this%FLEN-2, '.', this%DECS,',A,A)'
    ! write format string for one output line
    WRITE (this%linefmt, '(A,I0,A)') "(A", this%linelen-1, ")"
  END SUBROUTINE initfileio_gnuplot


  RECURSIVE SUBROUTINE writeheaderstring(string,root,k,prefix)
  IMPLICIT NONE
  !------------------------------------------------------------------------!
    TYPE(dict_typ), POINTER    :: root,node,subnode
    CHARACTER(LEN=*)           :: string
    CHARACTER(LEN=*), OPTIONAL :: prefix
    CHARACTER(LEN=128)         :: buf
    !------------------------------------------------------------------------!
    INTEGER                    :: idummy, k
    LOGICAL                    :: ldummy
    CHARACTER(LEN=128)         :: cdummy
    REAL                       :: rdummy
    !------------------------------------------------------------------------!
    INTENT(INOUT)              :: string,k
    !------------------------------------------------------------------------!

    node => root
    DO WHILE(ASSOCIATED(node))
       SELECT CASE(getdatatype(node))
       CASE(dict_int)
          CALL getattr(node,getkey(node),idummy)
          WRITE(buf,'(A1,A25,I14,A)')'#',trim(getkey(node))//": ",idummy, linsep
          WRITE(string(k:),'(A)')buf
          k = k + len(trim(buf))
       CASE(dict_real)
          CALL getattr(node,getkey(node),rdummy)
          WRITE(buf,'(A1,A25,ES14.5,A)')'#',trim(getkey(node))//": ",rdummy, linsep
          WRITE(string(k:),'(A)')buf
          k = k + len(trim(buf))
       CASE(dict_char)
!          CALL GetAttr(node,GetKey(node),cdummy)
!          WRITE(buf,'(A1,A25,A,A)')'#',TRIM(GetKey(node))//": ",TRIM(cdummy), LINSEP
!          WRITE(string(k:),'(A)')buf
!          k = k + LEN(TRIM(buf))
       CASE(dict_bool)
          CALL getattr(node,getkey(node),ldummy)
          WRITE(buf,'(A1,A25,L14,A)')'#',trim(getkey(node))//": ",ldummy, linsep
          WRITE(string(k:),'(A)')buf
          k = k + len(trim(buf))
       END SELECT
       IF(haschild(node)) THEN
          IF (present(prefix)) THEN
             WRITE(buf,'(A)')'#  ['//trim(prefix)//'/'//trim(getkey(node))//']' // linsep
          ELSE
             WRITE(buf,'(A)')'#  ['//trim(getkey(node))//']' // linsep
          END IF
          WRITE(string(k:),'(A)')buf
          k = k + len(trim(buf))
          IF (present(prefix)) THEN
             buf = trim(prefix)//'/'//trim(getkey(node))
          ELSE
             buf = trim(getkey(node))
          END IF
          subnode => getchild(node)
          CALL writeheaderstring(string,subnode,k,trim(buf))
       END IF
       node => getnext(node)
    END DO
  END SUBROUTINE writeheaderstring

  RECURSIVE SUBROUTINE getoutputpointer(this,Mesh,node,k)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base)      :: this
    CLASS(mesh_base)        :: mesh
    TYPE(dict_typ), POINTER :: node
    INTEGER                 :: k
    !------------------------------------------------------------------------!
    TYPE(dict_typ), POINTER           :: dir
    REAL, DIMENSION(:,:), POINTER     :: dummy2
    REAL, DIMENSION(:,:,:), POINTER   :: dummy3
    REAL, DIMENSION(:,:,:,:), POINTER :: dummy4
    INTEGER, DIMENSION(4)             :: dims
    INTEGER                           :: dim3,dim4,i,j
    !------------------------------------------------------------------------!
    INTENT(IN)    :: mesh
    INTENT(INOUT) :: this,k
    !------------------------------------------------------------------------!
    ! reset error code
    this%error = 0
    DO WHILE(ASSOCIATED(node))
      ! check for directory and exclude any coordinates (these are handled elsewhere)
      IF(haschild(node)) THEN
      ! recursion
        CALL getattr(node,getkey(node),dir)
        CALL getoutputpointer(this,mesh,dir,k)
      ELSE IF (hasdata(node).AND. .NOT.(getkey(node).EQ."bary_curv".OR. &
        getkey(node).EQ."bary_centers".OR.getkey(node).EQ."corners")) THEN
        WRITE(this%linebuf(k*this%FLEN:),trim(this%fmtstr))trim(getkey(node))
      ! value found!
        SELECT CASE(getdatatype(node))
        CASE(dict_real_twod)
          CALL getattr(node,getkey(node),dummy2)
          dims(1:2) = shape(dummy2)
          IF((dims(1).EQ.mesh%IMAX-mesh%IMIN+1).AND.&
             (dims(2).EQ.mesh%JMAX-mesh%JMIN+1)) THEN
            k = k+1
            IF (k .GT. this%maxcols) THEN
               this%error = 1
               EXIT
            END IF
            this%output(k)%val=> dummy2
          END IF
        CASE(dict_real_threed)
          CALL getattr(node,getkey(node),dummy3)
          dims(1:3) = shape(dummy3)
          IF((dims(1).EQ.mesh%IMAX-mesh%IMIN+1).AND.&
             (dims(2).EQ.mesh%JMAX-mesh%JMIN+1)) THEN
            dim3 = SIZE(dummy3, dim = 3)
            IF (k+dim3 .GT. this%maxcols) THEN
               this%error = 1
               EXIT
            END IF
            DO i=k+1, k+dim3
              this%output(i)%val => dummy3(:,:,i-k)
            END DO
            k = k+dim3
          END IF
        CASE(dict_real_fourd)
          CALL getattr(node,getkey(node),dummy4)
          dims(1:4) = shape(dummy4)
          IF((dims(1).EQ.mesh%IMAX-mesh%IMIN+1).AND.&
             (dims(2).EQ.mesh%JMAX-mesh%JMIN+1)) THEN
            dim3 = SIZE(dummy4, dim = 3)
            dim4 = SIZE(dummy4, dim = 4)
            IF (k+dim3*dim4 .GT. this%maxcols) THEN
               this%error = 1
               EXIT
            END IF
            DO j=1,dim4
              DO i=1,dim3
              this%output(k+j+(i-1)*dim4)%val => dummy4(:,:,i,j)
              END DO
            END DO
            k = k+dim3*dim4
          END IF
        CASE DEFAULT
          !do nothing (wrong type)
        END SELECT
      END IF
      node=>getnext(node)
    END DO
    IF (this%error.NE.0) &
         CALL error(this,"GetOutputPointer_gnuplot","number of output fields exceeds upper limit")
  END SUBROUTINE getoutputpointer


 PURE SUBROUTINE adjusttimestep(this,time,dt,dtcause)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    REAL               :: time
    REAL               :: dt
    INTEGER            :: dtcause
    !------------------------------------------------------------------------!
    INTENT(IN)         :: this
    INTENT(INOUT)      :: time,dt,dtcause
    !------------------------------------------------------------------------!
    IF ((time+dt)/this%time.GT.1.0) THEN
       dt = this%time - time
       dtcause = dtcause_fileio
    ELSE IF((time+1.5*dt)/this%time.GT.1.0) THEN
       dt = 0.5*(this%time - time)
       dtcause = dtcause_fileio
    END IF
  END SUBROUTINE adjusttimestep

  PURE SUBROUTINE inctime(this)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base), INTENT(INOUT) :: this
    !------------------------------------------------------------------------!
    this%time = this%time + abs(this%stoptime) / this%count
    this%step = this%step + 1
  END SUBROUTINE inctime

  SUBROUTINE openfile(this,action,fformat)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    INTEGER            :: action
    CHARACTER(LEN=*)   :: fformat
    !------------------------------------------------------------------------!
#ifdef PARALLEL
    INTEGER(KIND=MPI_OFFSET_KIND) :: offset
#endif
    !------------------------------------------------------------------------!
    INTENT(IN)       :: action,fformat
    INTENT(INOUT)    :: this
    !------------------------------------------------------------------------!
    SELECT CASE(action)
    CASE(readonly)
#ifdef PARALLEL
       CALL mpi_file_open(mpi_comm_world,getfilename(this),mpi_mode_rdonly, &
            mpi_info_null,this%handle,this%error)
        this%offset = 0
        CALL mpi_file_seek(this%handle,this%offset,mpi_seek_set,this%error)
#else
       OPEN(this%unit,file=getfilename(this),form=fformat,status="OLD", &
            action="READ",position="REWIND",iostat=this%error)
       !REWIND (UNIT=this%unit,IOSTAT=this%error)
#endif
    CASE(readend)
#ifdef PARALLEL
       CALL mpi_file_open(mpi_comm_world,getfilename(this),ior(mpi_mode_rdonly,&
            mpi_mode_append),mpi_info_null,this%handle,this%error)
       ! opening in append mode doesn't seem to work for pvfs2, hence ...
       offset = 0
       CALL mpi_file_seek(this%handle,offset,mpi_seek_end,this%error)
       CALL mpi_file_sync(this%handle,this%error)
#else
       OPEN(this%unit,file=getfilename(this),form=fformat,status="OLD", &
            action="READ",position="APPEND",iostat=this%error)
#endif
    CASE(replace)
#ifdef PARALLEL
       CALL mpi_file_delete(getfilename(this),mpi_info_null,this%error)
       CALL mpi_file_open(mpi_comm_world,getfilename(this),ior(mpi_mode_wronly,&
            mpi_mode_create),mpi_info_null,this%handle,this%error)
#else
       OPEN(this%unit,file=getfilename(this),form=fformat,status="REPLACE",&
            action="WRITE",position="REWIND",iostat=this%error)
#endif
    CASE(append)
#ifdef PARALLEL
       CALL mpi_file_open(mpi_comm_world,getfilename(this),ior(mpi_mode_rdwr,&
            mpi_mode_append),mpi_info_null,this%handle,this%error)
       ! opening in append mode doesn't seem to work for pvfs2, hence ...
       offset = 0
       CALL mpi_file_seek(this%handle,offset,mpi_seek_end,this%error)
       CALL mpi_file_sync(this%handle,this%error)
#else
       OPEN(this%unit,file=getfilename(this),form=fformat,status="OLD",&
            action="READWRITE",position="APPEND",iostat=this%error)
#endif
    CASE DEFAULT
       CALL error(this,"OpenFile","Unknown access mode.")
    END SELECT
  END SUBROUTINE openfile


  SUBROUTINE openfile_gnuplot(this,action)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    INTEGER            :: action
    !------------------------------------------------------------------------!
    INTENT(IN)         :: action
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
    CALL this%OpenFile(action,ascii)
  END SUBROUTINE openfile_gnuplot


  SUBROUTINE closefile_gnuplot(this)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    !------------------------------------------------------------------------!
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
#ifdef PARALLEL
    CALL mpi_file_close(this%handle,this%error)
#else
    CLOSE(this%unit,iostat=this%error)
#endif
  END SUBROUTINE closefile_gnuplot

  SUBROUTINE writeheader(this)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    !------------------------------------------------------------------------!
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
    IF (getrank(this).EQ.0) THEN
#ifdef PARALLEL
       CALL mpi_file_write(this%handle,trim(header_buf),len(trim(header_buf)), &
            mpi_character,this%status,this%error)
#else
       WRITE (this%unit,fmt='(A)',iostat=this%error) trim(header_buf) !(1:HLEN-1)
#endif
    END IF
  END SUBROUTINE writeheader

  SUBROUTINE readheader_gnuplot(this,success)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    LOGICAL            :: success
    !------------------------------------------------------------------------!
    INTENT(OUT)        :: success
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
#ifdef PARALLEL
    IF (getrank(this).EQ.0) THEN
    END IF
#else
#endif
    success = .false.
  END SUBROUTINE readheader_gnuplot

  SUBROUTINE writetimestamp_gnuplot(this,time)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    REAL               :: time
    !------------------------------------------------------------------------!
    INTENT(IN)         :: time
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
    IF (getrank(this).EQ.0) THEN
#ifdef PARALLEL
#else
#endif
    END IF
!!$    CALL Warning(this,"WriteTimestamp_gnuplot",&
!!$         "function is not implemented")
  END SUBROUTINE writetimestamp_gnuplot

  SUBROUTINE readtimestamp_gnuplot(this,time)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    REAL               :: time
    !------------------------------------------------------------------------!
    INTENT(OUT)        :: time
    INTENT(INOUT)      :: this
    !------------------------------------------------------------------------!
    IF (getrank(this).EQ.0) THEN
#ifdef PARALLEL
#else
#endif
    END IF
    time = 0.0
!!$    CALL Warning(this,"ReadTimestamp_gnuplot",&
!!$         "function is not implemented")
  END SUBROUTINE readtimestamp_gnuplot

  SUBROUTINE writedataset(this,Mesh)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base) :: this
    CLASS(mesh_base)   :: Mesh
    !------------------------------------------------------------------------!
    INTEGER             :: i,j,k,l
#ifdef PARALLEL
    INTEGER(KIND=MPI_OFFSET_KIND) :: offset
    INTEGER                       :: request
#endif
    !------------------------------------------------------------------------!
    INTENT(IN)    :: mesh
    INTENT(INOUT) :: this
    !------------------------------------------------------------------------!

    IF (ASSOCIATED(timedisc%w)) THEN
      IF (mesh%FARGO.EQ.3.AND.mesh%SN_shear) THEN
        CALL physics%AddBackgroundVelocityX(mesh,timedisc%w,timedisc%pvar,timedisc%cvar)
      ELSE
        CALL physics%AddBackgroundVelocityY(mesh,timedisc%w,timedisc%pvar,timedisc%cvar)
      END IF
    END IF


#ifdef PARALLEL
    ! be sure to write at the end by getting the offset from the file's size
    CALL mpi_file_get_size(this%handle,offset,this%error)
    ! very importan
    CALL mpi_barrier(mpi_comm_world,this%error)
    ! write _one_ line feed at the beginning of each time step
    IF (getrank(this).EQ.0) THEN
       CALL mpi_file_write_at(this%handle, offset, lf, 1, mpi_character, &
            this%status, this%error)
    END IF
    ! add the initial line feed and the general offset (depends on Mesh%IMIN)
    offset = offset + 1
    ! create the file view
    CALL mpi_file_set_view(this%handle,offset,this%basictype,this%filetype, &
         'native',mpi_info_null,this%error)
#else
    ! write _one_ line feed at the beginning of each time step
    WRITE (this%unit,fmt='(A)',advance='NO') lf
#endif

    DO k=1,this%cols
      ! trim the data for gnuplot output
      WHERE (abs(this%output(k)%val(:,:)).LT.max(tiny(this%output(k)%val(:,:)),1.0d-99))
        this%output(k)%val(:,:) = 0.0e+00
      END WHERE
    END DO

    ! start i,j from 1 to rank local size, because RemapBounds has not been
    ! used/cannot be used.
    DO i=1,this%inum
       DO j=1,this%jnum
          ! write positions to line buffer
          DO k=1,this%COLS-1
             WRITE (this%linebuf((k-1)*this%FLEN+1:k*this%FLEN),trim(this%fmtstr)) &
                   this%output(k)%val(i,j), recsep
          END DO

          IF ((j.EQ.mesh%JNUM).AND.((mesh%JNUM.GT.1).OR.(mesh%INUM.EQ.i))) THEN
             ! finish the block
             WRITE (this%linebuf((this%COLS-1)*this%FLEN+1:this%linelen),trim(this%fmtstr)) &
                  this%output(this%COLS)%val(i,j), blksep
          ELSE
             ! finish the line
             WRITE (this%linebuf((this%COLS-1)*this%FLEN+1:this%linelen),trim(this%fmtstr)) &
                  this%output(this%COLS)%val(i,j), linsep
          END IF

#ifdef PARALLEL
          ! write line buffer to output buffer
          DO k=1,this%linelen
             this%outbuf(k,j-1+mesh%JMIN) = this%linebuf(k:k)
          END DO
#else
          ! write line buffer to output file
          WRITE (this%unit,fmt=trim(this%linefmt),advance='YES') this%linebuf(1:this%linelen-1)
#endif
       END DO
#ifdef PARALLEL
       !*****************************************************************!
       ! This collective call doesn't work for pvfs2 -> bug in ROMIO ?
!!$       CALL MPI_File_write_all(this%handle,this%binout,this%bufsize,&
!!$            this%basictype, this%status, this%error)
       !*****************************************************************!
       ! so we use these two commands instead
       CALL mpi_file_iwrite(this%handle,this%outbuf,this%bufsize,this%basictype,&
            request,this%error)
       CALL mpi_wait(request,this%status,this%error)
#endif
    END DO
  END SUBROUTINE writedataset

  SUBROUTINE readdataset_gnuplot(this,Mesh,Physics,Timedisc)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base)   :: this
    CLASS(mesh_base)     :: Mesh
    CLASS(physics_base)  :: Physics
    CLASS(timedisc_base) :: Timedisc
    !------------------------------------------------------------------------!
    INTENT(IN)           :: mesh,physics,timedisc
    INTENT(INOUT)        :: this
    !------------------------------------------------------------------------!
  END SUBROUTINE readdataset_gnuplot


  SUBROUTINE error(this,modproc,msg)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_base), INTENT(INOUT) :: this
    CHARACTER(LEN=*),  INTENT(IN)     :: modproc
    CHARACTER(LEN=*),  INTENT(IN)     :: msg
    !------------------------------------------------------------------------!
    IF (initialized(this)) &
         CALL this%CloseFile_gnuplot()
    CALL this%Error_fileio(modproc,msg)
  END SUBROUTINE error

  SUBROUTINE finalize(this)
    IMPLICIT NONE
    !------------------------------------------------------------------------!
    CLASS(fileio_gnuplot) :: this
    !------------------------------------------------------------------------!
#ifdef PARALLEL
    DEALLOCATE(this%outbuf)
#endif
    DEALLOCATE(this%output)
    CALL this%Finalize_base()
  END SUBROUTINE finalize

END MODULE fileio_gnuplot_mod