81 CHARACTER,
PARAMETER ::
lf = achar(10)
114 CHARACTER(LEN=32) :: buf
115 CHARACTER(LEN=32) :: endianness
116 CHARACTER(LEN=12) :: realfmt
118 INTEGER :: ip1,jp1,kp1
119 REAL,
DIMENSION(:,:,:,:),
POINTER :: &
121 REAL,
DIMENSION(:,:) ,
POINTER :: &
123 REAL,
DIMENSION(:,:,:,:),
POINTER :: &
126 CHARACTER(LEN=64),
DIMENSION(:),
POINTER :: &
157 TYPE(
dict_typ),
INTENT(IN),
POINTER :: config
158 TYPE(
dict_typ),
INTENT(IN),
POINTER :: IO
160 TYPE(
output_typ),
DIMENSION(:),
POINTER :: poutput
163 CHARACTER(LEN=MAX_CHAR_LEN),
DIMENSION(2) :: skip
166 INTEGER,
DIMENSION(:),
POINTER :: sendbuf,recvbuf
169 REAL,
DIMENSION(:,:,:,:,:),
POINTER :: corners
175 CALL this%InitFileIO(mesh,physics,timedisc,sources,config,io,
"VTK",
"vts",textfile=.false.)
178 IF (this%datafile%cycles.NE.this%count+1) &
179 CALL this%Error(
'fileio_vtk::InitFileIO',
'filecycles = count+1 required for VTK output')
182 IF (abs(mesh%xmax-mesh%xmin).LT.2*epsilon(mesh%xmin) &
183 .OR.(mesh%INUM.EQ.1.AND.mesh%ROTSYM.EQ.1))
THEN
189 IF (abs(mesh%ymax-mesh%ymin).LT.2*epsilon(mesh%ymin) &
190 .OR.(mesh%JNUM.EQ.1.AND.mesh%ROTSYM.EQ.2))
THEN
196 IF (abs(mesh%zmax-mesh%zmin).LT.2*epsilon(mesh%zmin) &
197 .OR.(mesh%KNUM.EQ.1.AND.mesh%ROTSYM.EQ.3))
THEN
207 ALLOCATE(this%binout(
maxcomp, &
208 mesh%IMIN:mesh%IMAX, &
209 mesh%JMIN:mesh%JMAX, &
210 mesh%KMIN:mesh%KMAX), &
212 mesh%IMIN:mesh%IMAX+this%IP1, &
213 mesh%JMIN:mesh%JMAX+this%JP1, &
214 mesh%KMIN:mesh%KMAX+this%KP1), &
215 this%output(this%MAXCOLS), &
216 this%tsoutput(this%MAXCOLS), &
217 this%bflux(physics%VNUM,6), &
220 CALL this%Error(
"fileio_vtk::InitFileio_vtk",
"Unable to allocate memory.")
224 this%realsize = storage_size(this%binout)
225 SELECT CASE(this%realsize)
228 WRITE(this%linebuf,fmt=
'(I4)',iostat=this%err) this%realsize
229 WRITE(this%realfmt,fmt=
'(A)',iostat=this%err)
'"Float' // trim(adjustl(this%linebuf)) //
'"'
231 this%realsize = this%realsize/8
233 CALL this%Error(
"fileio_vtk::InitFileio_vtk", &
234 "only single, double or quadruple precision real numbers (4,8,16 bytes) supported")
238 CALL this%GetEndianness(this%endianness,
'"LittleEndian"',
'"BigEndian"')
248 IF (this%GetRank() .EQ. 0 )
THEN
249 n = this%GetNumProcs()
250 ALLOCATE(this%extent(0:(n-1)),sendbuf(6),recvbuf(6*n),stat = this%err)
253 ALLOCATE(sendbuf(6),recvbuf(1),stat = this%err)
256 CALL this%Error(
"fileio_vtk::InitFileio_vtk",
"Unable to allocate memory.")
260 sendbuf(1) = mesh%IMIN
261 sendbuf(2) = mesh%IMAX+this%IP1
262 sendbuf(3) = mesh%JMIN
263 sendbuf(4) = mesh%JMAX+this%JP1
264 sendbuf(5) = mesh%KMIN
265 sendbuf(6) = mesh%KMAX+this%KP1
268 CALL mpi_gather(sendbuf, 6, mpi_integer, recvbuf, 6, mpi_integer, &
269 0, mpi_comm_world, this%err)
272 IF (this%GetRank() .EQ. 0 )
THEN
274 WRITE(this%extent(i),fmt=
'(6(I7))',iostat=this%err)&
275 recvbuf(i*6+1),recvbuf(i*6+2),recvbuf(i*6+3),recvbuf(i*6+4),recvbuf(i*6+5),recvbuf(i*6+6)
279 DEALLOCATE(sendbuf,recvbuf,stat=this%err)
281 CALL this%Error(
"fileio_vtk::InitFileio_vtk",
"Deallocation of sendbuf/recvbuf failed.")
285 IF (this%cartcoords)
THEN
286 corners => mesh%RemapBounds(mesh%cart%corners)
288 corners => mesh%RemapBounds(mesh%curv%corners)
292 DO k=mesh%KMIN,mesh%KMAX+mesh%KP1
293 DO j=mesh%JMIN,mesh%JMAX+mesh%JP1
294 DO i=mesh%IMIN,mesh%IMAX+mesh%IP1
295 this%vtkcoords(1:3,i,j,k) = corners(i,j,k,1,1:3)
301 IF (mesh%INUM.EQ.1.AND.this%IP1.EQ.1)
THEN
304 DO k=mesh%KMIN,mesh%KMAX+mesh%KP1
305 DO j=mesh%JMIN,mesh%JMAX+mesh%JP1
306 this%vtkcoords(1:3,i+1,j,k) = corners(i,j,k,2,1:3)
311 IF (mesh%JNUM.EQ.1.AND.this%JP1.EQ.1)
THEN
314 DO k=mesh%KMIN,mesh%KMAX+mesh%KP1
315 DO i=mesh%IMIN,mesh%IMAX+mesh%IP1
316 this%vtkcoords(1:3,i,j+1,k) = corners(i,j,k,3,1:3)
321 IF (mesh%KNUM.EQ.1.AND.this%KP1.EQ.1)
THEN
324 DO j=mesh%JMIN,mesh%JMAX+mesh%JP1
325 DO i=mesh%IMIN,mesh%IMAX+mesh%IP1
326 this%vtkcoords(1:3,i,j,k+1) = corners(i,j,k,5,1:3)
333 CALL getattr(io,
"/timedisc/time",dummy1)
334 IF (
ASSOCIATED(dummy1%p))
THEN
335 this%tsoutput(1)%val => dummy1%p
336 this%tsoutput(1)%key =
"TIME"
345 skip(1:2) = [
CHARACTER(LEN=MAX_CHAR_LEN) ::
"corners",
"time"]
346 CALL this%GetOutputlist(mesh,node,k,this%TSCOLS,skip)
349 poutput => this%output
350 ALLOCATE(this%output,source=poutput(1:k),stat=this%err)
352 CALL this%Error(
"fileio_gnuplot::InitFileIO",
"memory allocation failed for this%output")
356 DO k = 1,
SIZE(this%output)
357 this%output(k)%numbytes = 0
358 DO i = 1,
SIZE(this%output(k)%p)
359 this%output(k)%numbytes = this%output(k)%numbytes +
SIZE(this%output(k)%p(i)%val)*this%realsize
365 IF (this%GetRank().EQ.0)
THEN
366 CALL this%pvtsfile%InitFilehandle(this%datafile%filename,this%datafile%path,
"pvts",&
367 textfile=.true.,onefile=.false.,cycles=this%count+1)
369 CALL this%pvdfile%InitFilehandle(this%datafile%filename,this%datafile%path,
"pvd",&
370 textfile=.true.,onefile=.true.,cycles=1)
390 IF (this%GetRank().EQ.0)
THEN
392 CALL this%pvdfile%OpenFile(
replace,this%step)
395 WRITE(unit=this%pvdfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err) &
396 '<?xml version="1.0"?>'//
lf &
397 //
'<VTKFile type="Collection" version="0.1" byte_order=' &
398 //trim(this%endianness)//
'>'//
lf &
399 //repeat(
' ',2)//
'<Collection>'
403 ftime = this%stoptime*k/this%count
404 IF (this%err.EQ. 0) &
405 WRITE(unit=this%pvdfile%GetUnitNumber(),fmt=
'(A,ES12.5,A)',iostat=this%err) &
406 repeat(
' ',4)//
'<DataSet timestep="',ftime,
'" part="0" file="' // &
408 trim(this%pvtsfile%GetBasename(k))//&
409 "." // trim(this%pvtsfile%extension)//
'"/>'
411 trim(this%datafile%GetBasename(k))//&
412 "." // trim(this%datafile%extension)//
'"/>'
426 IF (this%err.EQ. 0) &
427 WRITE(unit=this%pvdfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err) &
428 repeat(
' ',2)//
'</Collection>'//
lf//
'</VTKFile>'
430 CALL this%pvdfile%CloseFile(this%step)
432 IF (this%err.NE.0)
CALL this%Error(
"fileio_vtk::WriteParaviewFile",
"cannot write pvd-file")
580 WRITE(this%linebuf,fmt=
'(A,6(I7),A)')
'<?xml version="1.0"?>' //
lf// &
581 '<VTKFile type="StructuredGrid" version="0.1" byte_order='&
582 //trim(this%endianness)//
'>'//
lf//&
583 ' <StructuredGrid WholeExtent="',&
584 mesh%IMIN,mesh%IMAX+this%IP1,mesh%JMIN,mesh%JMAX+this%JP1, &
585 mesh%KMIN,mesh%KMAX+this%KP1,
'">'//
lf// &
586 ' <FieldData>'//
lf// &
587 ' <InformationKey type="String" Name="fosite version" format="ascii">'//
lf// &
588 ' '//trim(version)//
lf// &
589 ' </InformationKey>'//
lf
590 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) trim(this%linebuf)
594 IF (this%GetRank() .EQ. 0 )
THEN
595 CALL this%pvtsfile%OpenFile(
replace,this%step)
596 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A,6(I7),A)',iostat=this%err)
'<?xml version="1.0"?>' //
lf// &
597 '<VTKFile type="PStructuredGrid" version="0.1" byte_order='&
598 //trim(this%endianness)//
'>'//
lf//&
599 ' <PStructuredGrid WholeExtent="',&
600 1,mesh%INUM+this%IP1,1,mesh%JNUM+this%JP1,1,mesh%KNUM+this%KP1,
'">'//
lf//&
601 ' <PFieldData>'//
lf// &
602 ' <InformationKey type="String" Name="fosite version" format="ascii">'//
lf// &
603 ' '//trim(version)//
lf//&
606 CALL this%Error(
"fileio_vtk::WriteHeader",
"Writing pvts file header failed")
612 WRITE(this%linebuf,fmt=
'(A,ES14.7,A)') &
613 ' <DataArray type='//trim(this%realfmt)//
' Name="'//trim(this%tsoutput(k)%key)// &
614 '" NumberOfTuples="1" format="ascii">'//
lf// &
615 repeat(
' ',8),this%tsoutput(k)%val,
lf//&
617 IF (this%err.EQ.0)
THEN
620 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) trim(this%linebuf) //
lf
623 IF (this%GetRank() .EQ. 0 )
THEN
624 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A)',advance=
'NO',iostat=this%err) trim(this%linebuf)
625 CALL this%pvtsfile%CloseFile(this%step)
631 IF (this%err.NE.0)
CALL this%Error(
"fileio_vtk::WriteHeader",
"cannot write (P)VTS file header")
665 CHARACTER(LEN=256) :: filename
681 CALL this%WriteParaviewFile()
693 IF (this%GetRank() .EQ. 0 )
THEN
694 CALL this%pvtsfile%OpenFile(
append,this%step)
695 IF(this%err .NE. 0)
CALL this%Error(
"fileio_vtk::WriteDataset_vtk",
"Can't open pvts file")
747 WRITE(this%linebuf,fmt=
'(A,(6(I7)),A,I8,A)')&
748 ' </FieldData>'//
lf// &
749 ' <Piece Extent="', &
750 mesh%IMIN,mesh%IMAX+this%IP1,mesh%JMIN,mesh%JMAX+this%JP1,&
751 mesh%KMIN,mesh%KMAX+this%KP1,
'">'//
lf// &
753 ' <DataArray type='//trim(this%realfmt)//
' NumberOfComponents="3"'//&
754 ' Name="Point" format="appended" offset="',offset,
'">'//
lf// &
755 ' </DataArray>'//
lf// &
759 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) trim(this%linebuf)
762 offset = offset +
SIZE(this%vtkcoords)*this%realsize + 4
765 IF (this%GetRank() .EQ. 0 )
THEN
766 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err)&
767 repeat(
' ',4)//
'</PFieldData>' &
768 //
lf//repeat(
' ',4)//
'<PPoints>' &
769 //
lf//repeat(
' ',6)//
'<DataArray type='//trim(this%realfmt) &
770 //
' NumberOfComponents="3" Name="Point"/>' &
771 //
lf//repeat(
' ',4)//
'</PPoints>' &
772 //
lf//repeat(
' ',4)//
'<PCellData>'
778 DO k = 1,
SIZE(this%output)
779 WRITE(this%linebuf,fmt=
'(A,I3,A,I8,A)')&
780 repeat(
' ',8)//
'<DataArray type='//this%realfmt//&
781 ' NumberOfComponents="',
SIZE(this%output(k)%p), &
782 '" Name="'//trim(this%output(k)%key)//&
783 '" format="appended" offset="',offset,
'"/>'//
lf
785 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) trim(this%linebuf)
788 IF (this%GetRank() .EQ. 0 )
THEN
789 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A,I3,A)',iostat=this%err)&
790 repeat(
' ',6)//
'<DataArray type='//this%realfmt//&
791 ' NumberOfComponents="',
SIZE(this%output(k)%p),&
792 '" Name="'//trim(this%output(k)%key)//
'"/>'
797 offset = offset + this%output(k)%numbytes + 4
801 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) &
802 repeat(
' ',6)//
'</CellData>' &
803 //
lf//repeat(
' ',4)//
'</Piece>' &
804 //
lf//repeat(
' ',2)//
'</StructuredGrid>' &
805 //
lf//repeat(
' ',2)//
'<GlobalData>'
808 IF (this%GetRank() .EQ. 0 )
THEN
809 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err) repeat(
' ',4)//
'</PCellData>'
811 SELECT TYPE(df=>this%datafile)
813 DO i=0,this%GetNumProcs()-1
814 filename = trim(df%path) // &
815 trim(df%filename) // &
816 trim(df%GetRankString(i)) // &
817 trim(df%GetStepString(this%step)) // &
818 "." // trim(this%datafile%extension)
819 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err) &
820 repeat(
' ',4)//
'<Piece Extent="'//trim(this%extent(i))&
821 //
'" Source="'//trim(filename)//
'"/>'
824 CALL this%Error(
"fileio_vtk::WriteDataset",
"datafile must be of type filehandle_vts in parallel mode")
827 WRITE(unit=this%pvtsfile%GetUnitNumber(),fmt=
'(A)',iostat=this%err)&
828 repeat(
' ',2)//
'</PStructuredGrid>' &
857 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err)&
858 lf//repeat(
' ',2)//
'</GlobalData>' &
859 //
lf//
' <AppendedData encoding="raw">' &
862 CALL this%Error(
"fileio_vtk::WriteDataset",
"writing data array specs failed")
865 IF (this%GetRank() .EQ. 0 )
THEN
866 CALL this%pvtsfile%CloseFile(this%step)
867 IF(this%err .NE. 0)
CALL this%Error(
"fileio_vtk::WriteDataset_vtk",
"Can't close pvts file")
878 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) &
879 int(
SIZE(this%vtkcoords)*this%realsize + 4,kind=4),this%vtkcoords(:,:,:,:)
881 CALL this%Error(
"fileio_vtk::WriteDataset",
"writing coordinates failed")
891 DO l = 1,
SIZE(this%output)
893 n =
SIZE(this%output(l)%p)
895 DO k=1,mesh%KMAX-mesh%KMIN+1
896 DO j=1,mesh%JMAX-mesh%JMIN+1
897 DO i=1,mesh%IMAX-mesh%IMIN+1
898 this%binout(m,mesh%IMIN+i-1,mesh%JMIN+j-1,mesh%KMIN+k-1) &
899 = this%output(l)%p(m)%val(i,j,k)
905 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err) &
906 int(this%output(l)%numbytes,kind=4), &
908 ((((this%binout(m,i,j,k),m=1,n),i=mesh%IMIN,mesh%IMAX), &
909 j=mesh%JMIN,mesh%JMAX),k=mesh%KMIN,mesh%KMAX)
912 CALL this%Error(
"fileio_vtk::WriteDataset",
"writing data")
914 WRITE(unit=this%datafile%GetUnitNumber(),iostat=this%err)
lf//repeat(
' ',2)//&
915 '</AppendedData>'//
lf//
'</VTKFile>'//
lf
927 IF (
ASSOCIATED(this%binout))
DEALLOCATE(this%binout)
928 IF (
ASSOCIATED(this%vtkcoords))
DEALLOCATE(this%vtkcoords)
929 IF (
ASSOCIATED(this%bflux))
DEALLOCATE(this%bflux)
930 NULLIFY(this%binout,this%vtkcoords,this%bflux)
935 SUBROUTINE initfilehandle(this,filename,path,extension,textfile,onefile,cycles,unit)
939 CHARACTER(LEN=*),
INTENT(IN) :: filename
940 CHARACTER(LEN=*),
INTENT(IN) :: path
941 CHARACTER(LEN=*),
OPTIONAL,
INTENT(IN) :: extension
942 LOGICAL,
OPTIONAL,
INTENT(IN) :: textfile
943 LOGICAL,
OPTIONAL,
INTENT(IN) :: onefile
944 INTEGER,
OPTIONAL,
INTENT(IN) :: cycles
945 INTEGER,
OPTIONAL,
INTENT(IN) :: unit
947 CALL this%filehandle_fortran%InitFilehandle(filename,path,extension,textfile,onefile,cycles,unit)
950 CALL this%Error(
"filehandle_vts::InitFilehandle", &
951 "we need at least 3 bytes for rank string: 2 for '-r' + 1 for process number (0..9)")
953 CALL this%Error(
"filehandle_vts::InitFilehandle", &
954 "number of processes for multiple file output exceeds limits")
962 INTEGER,
INTENT(IN) :: step
963 CHARACTER(LEN=256) :: fname
965 fname = trim(
this%path) // trim(
this%filename) // trim(
this%GetRankString()) &
966 // trim(
this%GetStepString(step))
974 INTEGER,
OPTIONAL,
INTENT(IN) :: rank
977 CHARACTER(LEN=16) :: fmtstr
981 IF (
PRESENT(rank))
THEN
Dictionary for generic data types.
type(logging_base), save this
integer, parameter append
read/write access at end
subroutine initfilehandle(this, filename, path, extension, textfile, onefile, cycles, unit)
basic initialization of Fortran file handle
integer, parameter replace
read/write access replacing file
character(len=256) function getbasename(this, step)
get the file name without path and extension but with step string appended; if thisonefile suppress t...
I/O for GNUPLOT readable tabular files.
subroutine finalize(this)
Closes the file I/O.
subroutine writeheader(this, Mesh, Physics, Header, IO)
Writes the configuration as a header to the file.
character, parameter lf
line feed
I/O for VTK files in XML format (vtkStructuredGrid)
subroutine writedataset_vtk(this, Mesh, Physics, Fluxes, Timedisc, Header, IO)
Reads the header (not yet implemented)
subroutine initfileio_vtk(this, Mesh, Physics, Timedisc, Sources, config, IO)
Constructor for the VTK file I/O.
integer, parameter maxcols
max of different output arrays
integer, parameter rank_str_len
length of multi process rank string
subroutine writeparaviewfile(this)
Write the Paraview global description file.
integer, parameter maxcomp
max. of allowed components 9 is a tensor (rank 2, dim 3)
character(len=rank_str_len) function getrankstring(this, rank)
convert process rank to string
base module for numerical flux functions
base class for geometrical properties
module to manage list of source terms
class for Fortran file handle
FileIO class for VTK output.
container class to manage the list of source terms