fileio_xdmf.f90
Go to the documentation of this file.
1!#############################################################################
2!# #
3!# fosite - 3D hydrodynamical simulation program #
4!# module: fileio_xdmf.f90 #
5!# #
6!# Copyright (C) 2013-2024 #
7!# Manuel Jung <mjung@astrophysik.uni-kiel.de> #
8!# Jannes Klee <jklee@astrophysik.uni-kiel.de> #
9!# Tobias Illenseer <tillense@astrophysik.uni-kiel.de> #
10!# #
11!# This program is free software; you can redistribute it and/or modify #
12!# it under the terms of the GNU General Public License as published by #
13!# the Free Software Foundation; either version 2 of the License, or (at #
14!# your option) any later version. #
15!# #
16!# This program is distributed in the hope that it will be useful, but #
17!# WITHOUT ANY WARRANTY; without even the implied warranty of #
18!# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, GOOD TITLE or #
19!# NON INFRINGEMENT. See the GNU General Public License for more #
20!# details. #
21!# #
22!# You should have received a copy of the GNU General Public License #
23!# along with this program; if not, write to the Free Software #
24!# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. #
25!# #
26!#############################################################################
27!----------------------------------------------------------------------------!
51!----------------------------------------------------------------------------!
61 USE common_dict
62 !--------------------------------------------------------------------------!
63 PRIVATE
64 INTEGER, PARAMETER :: blk_indent = 2
65 CHARACTER, PARAMETER :: sp = achar(32)
66 CHARACTER, PARAMETER :: lf = achar(10)
67 CHARACTER(LEN=2), PARAMETER :: tb = repeat(sp,blk_indent)
68 !--------------------------------------------------------------------------!
70 TYPE, EXTENDS(fileio_binary) :: fileio_xdmf
71 TYPE(filehandle_fortran) :: xmffile
72 CHARACTER(LEN=14) :: endian_xdmf
73 CONTAINS
74 PROCEDURE :: initfileio_deferred => initfileio_xdmf
75 PROCEDURE :: iteratedict
76 PROCEDURE :: writexmf
77 PROCEDURE :: writenode_xdmf
78 PROCEDURE :: writekey_xdmf
79 PROCEDURE :: writeattribute
80 PROCEDURE :: writemeshxml
81 PROCEDURE :: writevector
82 PROCEDURE :: writedataitem
83 final :: finalize
84 END TYPE
85 !--------------------------------------------------------------------------!
86 PUBLIC :: &
87 ! types
89 ! constants
90 ! methods
91 !--------------------------------------------------------------------------!
92
93CONTAINS
98 SUBROUTINE initfileio_xdmf(this,Mesh,Physics,Timedisc,Sources,config,IO)
99 IMPLICIT NONE
100 !------------------------------------------------------------------------!
101 CLASS(fileio_xdmf), INTENT(INOUT) :: this
102 CLASS(mesh_base), INTENT(IN) :: Mesh
103 CLASS(physics_base), INTENT(IN) :: Physics
104 CLASS(timedisc_base), INTENT(IN) :: Timedisc
105 CLASS(sources_list), ALLOCATABLE, INTENT(IN) :: Sources
106 TYPE(dict_typ), INTENT(IN), POINTER &
107 :: config
108 TYPE(dict_typ), INTENT(IN), POINTER &
109 :: IO
110 !------------------------------------------------------------------------!
111 INTEGER :: idx(1)
112 !------------------------------------------------------------------------!
113 !\attention interface changed, however "xdmf" is not passed anymore
114 IF (.NOT.this%Initialized()) &
115 CALL this%InitFileio(mesh,physics,timedisc,sources,config,io,"xdmf","bin",textfile=.false.)
116
117 CALL this%fileio_binary%InitFileIO(mesh,physics,timedisc,sources,config,io)
118
119 CALL this%GetEndianness(this%endian_xdmf, '"Little"', '"Big"')
120
121 SELECT CASE (mesh%NDIMS)
122 CASE(2,3)
123 ! 2D & 3D -> do nothing
124 CASE DEFAULT
125 CALL this%Error("fileio_xdmf::InitFileIO_xdmf","only 2D and 3D mesh is supported")
126 END SELECT
127
128 ! create file handle for the xmf file; only on rank 0 in parallel mode
129#ifdef PARALLEL
130 IF (this%GetRank().EQ.0) &
131#endif
132 CALL this%xmffile%InitFilehandle(this%datafile%filename,this%datafile%path,"xmf",textfile=.false.,onefile=.true.,cycles=1)
133
134 CALL this%WriteXMF(mesh,io)
135 END SUBROUTINE initfileio_xdmf
136
137
139 RECURSIVE SUBROUTINE iteratedict(this,Mesh,config,offset,filename,path,indent)
140 IMPLICIT NONE
141 !------------------------------------------------------------------------!
142 CLASS(fileio_xdmf), INTENT(INOUT) :: this
143 CLASS(mesh_base), INTENT(IN) :: mesh
144 TYPE(dict_typ), INTENT(IN), POINTER &
145 :: config
146 CHARACTER(LEN=*), INTENT(IN) :: filename
147 CHARACTER(LEN=*), INTENT(IN), OPTIONAL &
148 :: path
149 INTEGER, INTENT(INOUT) :: offset
150 INTEGER :: indent
151 !------------------------------------------------------------------------!
152 CHARACTER(LEN=MAX_CHAR_LEN) :: str, key
153 TYPE(dict_typ),POINTER :: node
154 !------------------------------------------------------------------------!
155 IF(PRESENT(path)) THEN
156 str = path
157 ELSE
158 str = ""
159 ENDIF
160 node => config
161 DO WHILE(ASSOCIATED(node))
162 key = trim(str)//"/"//trim(getkey(node))
163
164 CALL this%WriteNode_xdmf(mesh,key,node,offset,filename,indent)
165
166 IF(haschild(node)) THEN
167 CALL this%IterateDict(mesh, getchild(node), offset, filename, key, indent)
168 END IF
169 node=>getnext(node)
170 END DO
171 END SUBROUTINE iteratedict
172
177 SUBROUTINE writekey_xdmf(this,offset,key,type,bytes,dims)
178 IMPLICIT NONE
179 !------------------------------------------------------------------------!
180 CLASS(fileio_xdmf), INTENT(INOUT) :: this
181 CHARACTER(LEN=*), INTENT(IN) :: key
182 INTEGER, INTENT(INOUT) :: offset
183 INTEGER, INTENT(IN) :: type
184 INTEGER, INTENT(IN) :: bytes
185 INTEGER,DIMENSION(5), OPTIONAL, INTENT(IN) &
186 :: dims
187 !------------------------------------------------------------------------!
188 INTEGER :: keylen,l
189 !------------------------------------------------------------------------!
190 keylen = len_trim(key)
191 IF(PRESENT(dims)) THEN
192 l = 3
193 IF(dims(4).GT.1) l = 4
194 IF(dims(5).GT.1) l = 5
195 ELSE
196 l = 0
197 END IF
198 offset = offset + (3+l) * this%intsize + keylen
199 END SUBROUTINE writekey_xdmf
200
201
206 SUBROUTINE writenode_xdmf(this,Mesh,key,node,offset,filename,indent)
207 IMPLICIT NONE
208 !------------------------------------------------------------------------!
209 CLASS(fileio_xdmf), INTENT(INOUT) :: this
210 CLASS(mesh_base), INTENT(IN) :: Mesh
211 CHARACTER(LEN=MAX_CHAR_LEN), INTENT(IN) &
212 :: key
213 TYPE(dict_typ), INTENT(IN), POINTER &
214 :: node
215 INTEGER , INTENT(INOUT) :: offset
216 CHARACTER(LEN=*), INTENT(IN) :: filename
217 INTEGER, INTENT(IN) :: indent
218 !------------------------------------------------------------------------!
219 REAL,DIMENSION(:,:,:), POINTER :: ptr3 => null()
220 REAL,DIMENSION(:,:,:,:), POINTER :: ptr4 => null()
221 REAL,DIMENSION(:,:,:,:,:),POINTER :: ptr5 => null()
222 INTEGER,DIMENSION(5) :: dims
223 INTEGER :: bytes
224 CHARACTER(LEN=1),DIMENSION(:),POINTER &
225 :: val
226 LOGICAL :: ref
227 INTEGER :: type
228 !------------------------------------------------------------------------!
229 ref = .false.
230 dims(:) = 1
231
232 type = getdatatype(node)
233 SELECT CASE(type)
234 CASE(dict_real_p)
235 bytes = this%realsize
236 type = dict_real
237 CASE(dict_real_threed)
238 CALL getattr(node,key,ptr3)
239 dims(1:3) = shape(ptr3)
240 CASE(dict_real_fourd)
241 CALL getattr(node,key,ptr4)
242 dims(1:4) = shape(ptr4)
243 CASE(dict_real_fived)
244 CALL getattr(node,key,ptr5)
245 dims(1:5) = shape(ptr5)
246 CASE DEFAULT
247 val => getdata(node)
248 IF(ASSOCIATED(val)) THEN
249 bytes = SIZE(val)
250 ELSE
251 bytes = 0
252 END IF
253 END SELECT
254
255 IF(product(dims).GT.1) THEN
256 CALL this%SetOutputDims(mesh,dims)
257 bytes = product(dims) * this%realsize
258 CALL this%WriteKey_xdmf(offset,key,type,bytes,dims)
259 ELSE
260 CALL this%WriteKey_xdmf(offset,key,type,bytes)
261 END IF
262
263 SELECT CASE(getdatatype(node))
264 CASE(dict_real_threed)
265 CALL this%WriteAttribute(mesh,trim(key),dims(1:3),filename,offset,ref,indent)
266! It is not clear, how other datatypes than DICT_REAL_TWOD should be mapped in
267! XDMF.
268! CASE(DICT_REAL_THREED)
269! CALL this%WriteAttribute(Mesh,TRIM(key),dims(1:3),filename,offset,ref)
270! CASE(DICT_REAL_FOURD)
271! CALL this%WriteAttribute(Mesh,TRIM(key),dims(1:4),filename,offset,ref)
272 END SELECT
273
274 offset = offset + bytes
275
276 END SUBROUTINE writenode_xdmf
277
278
281 SUBROUTINE writedataitem(this,dims,filename,offset,indent)
282 IMPLICIT NONE
283 !------------------------------------------------------------------------!
284 CLASS(fileio_xdmf), INTENT(INOUT) :: this
285 CHARACTER(LEN=*), INTENT(IN) :: dims
286 CHARACTER(LEN=*), INTENT(IN) :: filename
287 INTEGER, INTENT(IN) :: offset
288 INTEGER, INTENT(IN) :: indent
289 !------------------------------------------------------------------------!
290 CHARACTER(LEN=16) :: seek
291 !------------------------------------------------------------------------!
292 WRITE(seek,'(I16)',iostat=this%err) offset
293
294 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
295 repeat(tb,indent) // '<DataItem Dimensions=' // trim(dims) // sp &
296 // 'NumberType="Float" Precision=' // trim(this%realfmt)
297 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
298 sp // 'Format="Binary" Endian=' // trim(this%endian_xdmf) // &
299 sp // 'Seek="' // trim(adjustl(seek)) // '">' // lf
300 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
301 repeat(tb,indent+1) // trim(filename) // lf // &
302 repeat(tb,indent) // '</DataItem>' // lf
303 END SUBROUTINE writedataitem
304
307 SUBROUTINE writeattribute(this,Mesh,name,dims,filename,offset,ref,indent)
308 IMPLICIT NONE
309 !------------------------------------------------------------------------!
310 CLASS(fileio_xdmf), INTENT(INOUT) :: this
311 CLASS(mesh_base), INTENT(IN) :: Mesh
312 CHARACTER(LEN=*), INTENT(IN) :: name
313 CHARACTER(LEN=*), INTENT(IN) :: filename
314 INTEGER, INTENT(IN), DIMENSION(:) &
315 :: dims
316 INTEGER, INTENT(IN) :: offset
317 LOGICAL, INTENT(IN) :: ref
318 INTEGER, INTENT(IN) :: indent
319 !------------------------------------------------------------------------!
320 CHARACTER(LEN=32) :: dtype,center,dstr
321 INTEGER :: dsize, err
322 !------------------------------------------------------------------------!
323 !data size
324 dsize = SIZE(dims)
325 ! default (no mesh data)
326 center = "Grid"
327
328 ! check for mesh data
329 IF (dsize.GE.3) THEN
330 IF (all(dims(1:3).EQ.(/mesh%INUM,mesh%JNUM,mesh%KNUM/))) &
331 center = "Cell"
332 IF (all(dims(1:3).EQ.(/mesh%INUM+mesh%IP1,mesh%JNUM+mesh%JP1,mesh%KNUM+mesh%KP1/))) &
333 center = "Node"
334 IF (trim(center).NE."Grid") &
335 dsize = dsize - 3
336 END IF
337
338 SELECT CASE(dsize)
339 CASE(0)
340 dtype = "Scalar"
341 CASE(1)
342 dtype = "Vector"
343 CASE(2)
344 dtype = "Matrix"
345 CASE DEFAULT
346 CALL this%Error("fileio_xdmf::WriteAttribute","data type currently not supported")
347 END SELECT
348
349 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
350 repeat(tb,indent) // '<Attribute Name="' // trim(name) // '" ' &
351 // 'AttributeType="' // trim(dtype) // '" ' &
352 // 'Center="' // trim(center) // '">' // lf
353
354 CALL getdimsstr(dims,dstr)
355 CALL this%WriteDataItem(trim(dstr), trim(filename), offset,indent+1)
356
357 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
358 repeat(tb,indent) // '</Attribute>' // lf
359 END SUBROUTINE writeattribute
360
363 SUBROUTINE writevector(this,Mesh,name,dims,ref1,ref2,ref3,step)
364 IMPLICIT NONE
365 !------------------------------------------------------------------------!
366 CLASS(fileio_xdmf), INTENT(INOUT) :: this
367 CLASS(mesh_base), INTENT(IN) :: Mesh
368 CHARACTER(LEN=*), INTENT(IN) :: name
369 INTEGER, INTENT(IN), DIMENSION(3) &
370 :: dims
371 CHARACTER(LEN=*), INTENT(IN) :: ref1,ref2,ref3
372 CHARACTER(LEN=*), INTENT(IN) :: step
373 !------------------------------------------------------------------------!
374 CHARACTER(LEN=32) :: dstr
375 INTEGER :: err
376 !------------------------------------------------------------------------!
377 CALL getdimsstr(dims,dstr)
378 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
379 '<Attribute Name="' // trim(name) // '" ' &
380 // 'AttributeType="Vector" Center="Cell">' // lf&
381 // '<DataItem ItemType="Function" Function="JOIN($0, $1, $2)" '&
382 // 'Dimensions=' // trim(dstr) // '>' // lf &
383 ! first component
384 // '<DataItem Reference="/Xdmf/Domain/Grid'&
385 // "/Grid[@Name='step" // step // "']" &
386 // '/Attribute[@Name=' &
387 // "'" // trim(ref1) // "'" // ']/DataItem[1]"/>' // lf &
388 ! second component
389 // '<DataItem Reference="/Xdmf/Domain/Grid'&
390 // "/Grid[@Name='step" // step // "']" &
391 // '/Attribute[@Name=' &
392 // "'" // trim(ref2) // "'" // ']/DataItem[1]"/>' // lf &
393 ! third component
394 // '<DataItem Reference="/Xdmf/Domain/Grid'&
395 // "/Grid[@Name='step" // step // "']" &
396 // '/Attribute[@Name=' &
397 // "'" // trim(ref3) // "'" // ']/DataItem[1]"/>' // lf &
398 ! end
399 // '</DataItem>' // lf &
400 // '</Attribute>' // lf
401 END SUBROUTINE writevector
402
405 SUBROUTINE writemeshxml(this,Mesh,filename,indent)
406 IMPLICIT NONE
407 !------------------------------------------------------------------------!
408 CLASS(fileio_xdmf), INTENT(INOUT) :: this
409 CLASS(mesh_base), INTENT(IN) :: Mesh
410 CHARACTER(LEN=*), INTENT(IN) :: filename
411 INTEGER :: indent
412 !------------------------------------------------------------------------!
413 INTEGER :: dims(3)
414 CHARACTER(LEN=32) :: dstr
415 CHARACTER(LEN=14) :: gstr(3) = (/"'/mesh/grid_x'","'/mesh/grid_y'","'/mesh/grid_z'"/)
416 INTEGER :: i
417 !------------------------------------------------------------------------!
418 dims = (/ mesh%INUM+mesh%IP1, mesh%JNUM+mesh%JP1, mesh%KNUM+mesh%KP1 /)
419 CALL getdimsstr(dims(:),dstr)
420
421 SELECT CASE(mesh%NDIMS)
422 CASE(2)
423 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
424 repeat(tb,indent) // '<Topology TopologyType="2DSMesh" ' // 'NumberOfElements=' // trim(dstr) // '/>' // lf // &
425 repeat(tb,indent) // '<Geometry GeometryType="X_Y_Z">' // lf
426 CASE(3)
427 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
428 repeat(tb,indent) // '<Topology TopologyType="3DSMesh" ' // 'NumberOfElements=' // trim(dstr) // '/>' // lf // &
429 repeat(tb,indent) // '<Geometry GeometryType="X_Y_Z">' // lf
430 CASE DEFAULT
431 CALL this%Error("fileio_xdmf::WriteMeshXML","only 2D and 3D mesh is currently supported")
432 END SELECT
433
434 ! register all 3 coordinates even for 2D grids, because 2D curved surfaces may have 3D cartesian coordinates
435 ! e.g., surface of a sphere
436 DO i=1,3
437 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
438 repeat(tb,indent+1) // '<DataItem Reference="/Xdmf/Domain/Grid/Grid/Attribute[@Name='//gstr(i)//']/DataItem[1]"/>' // lf
439 END DO
440
441 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
442 repeat(tb,indent) // '</Geometry>' // lf
443 END SUBROUTINE writemeshxml
444
447 SUBROUTINE writexmf(this,Mesh,IO)
448 IMPLICIT NONE
449 !------------------------------------------------------------------------!
450 CLASS(fileio_xdmf), INTENT(INOUT) :: this
451 CLASS(mesh_base), INTENT(IN) :: Mesh
452 TYPE(dict_typ), INTENT(IN), POINTER &
453 :: IO
454 !------------------------------------------------------------------------!
455 INTEGER :: i, offset, err
456 REAL :: ftime
457 CHARACTER(LEN=4) :: step
458 CHARACTER(LEN=32) :: time
459 CHARACTER(LEN=256) :: filename
460 TYPE(dict_typ), POINTER :: meshIO => null()
461 !------------------------------------------------------------------------!
462 this%err = 0
463 IF (this%GetRank().EQ.0) THEN
464 ! write a xmf-file
465 CALL this%xmffile%OpenFile(replace,this%step)
466
467 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
468 '<?xml version="1.0" ?>' // lf &
469 // '<!DOCTYPE Xdmf SYSTEM "Xdmf.dtd" []>' // lf &
470 // '<Xdmf Version="2.0">' // lf &
471 // tb // '<Domain>' // lf &
472 // tb // tb // '<Grid Name="mesh" GridType="Collection" CollectionType="Temporal">' // lf
473
474! this works in paraview, but not in visit :(. So we have to use the simpler method for now
475! // '<Time TimeType="List">' // LF &
476! // '<DataItem Format="XML" NumberType="Float" Dimensions="11">' // LF &
477! // '0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0' // LF &
478! // '</DataItem>' // LF &
479! // '</Time>' // LF
480! CALL WriteDataItem_xdmf(this, '"1"', TRIM(filename)//'_'//step//'.bin')
481
482 CALL getattr(io,"mesh",meshio)
483
484 DO i=0,this%count
485 offset = 13
486 ftime = this%stoptime/(this%count)*i
487 WRITE(step,'(I4.4)') i
488 WRITE(time,'(ES16.9)') ftime
489 WRITE(filename, '(A)') trim(this%datafile%GetFilename(i))
490
491 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
492 repeat(tb,3) // '<Grid Name="step' // step // '" GridType="Uniform">' // lf // &
493 repeat(tb,4) // '<Time Value="' // trim(time) // '" />' // lf
494
495 CALL this%WriteMeshXML(mesh,filename,indent=5)
496! quick HACK to bind velocity components into vector structure
497! CALL this%WriteVector(Mesh,"/timedisc/velocity", &
498! (/Mesh%INUM,Mesh%JNUM,Mesh%KNUM/), &
499! "/timedisc/xvelocity","/timedisc/yvelocity", &
500! "/timedisc/zvelocity", step)
501 CALL this%IterateDict(mesh, io, offset, filename,indent=5)
502
503 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
504 repeat(tb,3) // '</Grid>' // lf
505 END DO
506
507 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
508 tb // tb // '</Grid>' // lf
509
510
511 WRITE(unit=this%xmffile%GetUnitNumber(),iostat=this%err) &
512 tb // '</Domain>' // lf &
513 // '</Xdmf>' // lf
514
515 IF(this%err.NE.0) CALL this%Error("fileio_xdmf::WriteXMF", "Can't write xmf-file")
516 CALL this%xmffile%CloseFile(this%step)
517 END IF
518 END SUBROUTINE writexmf
519
521 PURE SUBROUTINE getdimsstr(dims,res)
522 IMPLICIT NONE
523 !------------------------------------------------------------------------!
524 INTEGER,DIMENSION(3), INTENT(IN) :: dims
525 CHARACTER(LEN=*), INTENT(OUT) :: res
526 !------------------------------------------------------------------------!
527 INTEGER :: i,idx(3)
528 CHARACTER(LEN=8) :: dim_str
529 !------------------------------------------------------------------------!
530 ! reverse dimensional order (probably because of Fortran vs. C array order)
531 idx(1:3) = dims(3:1:-1)
532
533 res = ''
534 DO i=1,3
535 IF (idx(i).GT.1) THEN ! skip suppressed dimension in 2D case
536 WRITE(dim_str,'(I8)') idx(i)
537 ! append to resulting string
538 res = trim(res) // ' ' // trim(adjustl(dim_str))
539 END IF
540 END DO
541 res = '"' // trim(adjustl(res)) // '"'
542 END SUBROUTINE getdimsstr
543
544 SUBROUTINE finalize(this)
545 IMPLICIT NONE
546 !-----------------------------------------------------------------------!
547 TYPE(fileio_xdmf), INTENT(INOUT) :: this
548 !-----------------------------------------------------------------------!
549 ! do nothing specific for xdmf
550 ! finalizer of fileio_binary is called automatically
551 END SUBROUTINE
552
553END MODULE fileio_xdmf_mod
Dictionary for generic data types.
Definition: common_dict.f90:61
integer, parameter, public dict_real
Definition: common_dict.f90:95
pointer, public getdata(root)
Return the datatype of node 'root'.
logical function, public haschild(root)
Check if the node 'root' has one or more children.
integer, parameter, public dict_real_fourd
integer, parameter, public dict_real_threed
integer function, public getdatatype(root)
Return the datatype of node 'root'.
function, public getkey(root)
Get the key of pointer 'root'.
integer, parameter, public dict_real_fived
integer, parameter, public dict_real_p
type(dict_typ) function, pointer, public getnext(root)
Get the pointer to the next child.
type(logging_base), save this
type(dict_typ) function, pointer, public getchild(root)
Get the pointer to a direct child of the pointer 'root'.
Generic file I/O module.
Definition: fileio_base.f90:54
integer, parameter replace
read/write access replacing file
module for binary file I/O
character, parameter lf
line feed
subroutine finalize(this)
Closes the file I/O.
module for XDMF file I/O
Definition: fileio_xdmf.f90:52
subroutine writenode_xdmf(this, Mesh, key, node, offset, filename, indent)
Write the xdmf node.
recursive subroutine iteratedict(this, Mesh, config, offset, filename, path, indent)
Iterate the dictionary and run a Subroutine on every node.
subroutine writevector(this, Mesh, name, dims, ref1, ref2, ref3, step)
Writes the mesh to file.
subroutine writeattribute(this, Mesh, name, dims, filename, offset, ref, indent)
Writes description of data item in xml syntax.
pure subroutine getdimsstr(dims, res)
determines the string for the dimension attribute, i.e. array dimensions of mesh data
subroutine writexmf(this, Mesh, IO)
Main routine to write all data to xmf file.
integer, parameter blk_indent
block indentation
Definition: fileio_xdmf.f90:64
subroutine writedataitem(this, dims, filename, offset, indent)
Writes description of data item in xml syntax.
subroutine writekey_xdmf(this, offset, key, type, bytes, dims)
Write the xdmf key.
subroutine initfileio_xdmf(this, Mesh, Physics, Timedisc, Sources, config, IO)
Constructor for the xdmf file I/O.
Definition: fileio_xdmf.f90:99
character, parameter sp
space
Definition: fileio_xdmf.f90:65
subroutine writemeshxml(this, Mesh, filename, indent)
Writes the mesh to file.
character(len=2), parameter tb
Definition: fileio_xdmf.f90:67
base module for numerical flux functions
Definition: fluxes_base.f90:39
base class for geometrical properties
basic mesh module
Definition: mesh_base.f90:72
Basic physics module.
module to manage list of source terms
class for Fortran file handle
Definition: fileio_base.f90:99
FileIO gnuplot class.
Definition: fileio_xdmf.f90:70
mesh data structure
Definition: mesh_base.f90:122
container class to manage the list of source terms