fileio_binary.f90
Go to the documentation of this file.
1!#############################################################################
2!# #
3!# fosite - 3D hydrodynamical simulation program #
4!# module: fileio_binary.f90 #
5!# #
6!# Copyright (C) 2015-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!----------------------------------------------------------------------------!
59!----------------------------------------------------------------------------!
60#define HAVE_VTK
69 USE common_dict
70#ifdef PARALLEL
71#ifdef HAVE_MPI_MOD
72 USE mpi
73#endif
74#endif
75 IMPLICIT NONE
76#ifdef PARALLEL
77#ifdef HAVE_MPIF_H
78 include 'mpif.h'
79#endif
80#endif
81 !--------------------------------------------------------------------------!
82 PRIVATE
83 CHARACTER, PARAMETER :: lf = achar(10)
84 !--------------------------------------------------------------------------!
86 TYPE, EXTENDS(fileio_base) :: fileio_binary
88 CHARACTER(LEN=12) :: realfmt
89 CHARACTER(LEN=14) :: endianness
90 INTEGER :: realsize
91 INTEGER :: intsize
92#ifdef PARALLEL
93 LOGICAL :: first
94 INTEGER :: cbufsize, & !< size of corner output
95 clsizes(3)
96 INTEGER(KIND=MPI_OFFSET_KIND) :: offset
97#else
98 INTEGER :: offset
99#endif
100 CONTAINS
102 PROCEDURE :: initfileio_deferred => initfileio_binary
103 PROCEDURE :: writeheader
104 !PROCEDURE :: ReadHeader
105 !PROCEDURE :: WriteTimestamp
106 !PROCEDURE :: ReadTimestamp
107 PROCEDURE :: writedataset_deferred => writedataset_binary
108 !PROCEDURE :: ReadDataset
109 final :: finalize
110 !PRIVATE
111 PROCEDURE :: hasmeshdims
112 PROCEDURE :: hascornerdims
113 PROCEDURE :: setoutputdims
114 PROCEDURE :: writenode
115 PROCEDURE :: writekey
117 END TYPE
118 !--------------------------------------------------------------------------!
119 PUBLIC :: &
120 ! types
122 !--------------------------------------------------------------------------!
123
124CONTAINS
129 SUBROUTINE initfileio_binary(this,Mesh,Physics,Timedisc,Sources,config,IO)
130 IMPLICIT NONE
131 !------------------------------------------------------------------------!
132 CLASS(fileio_binary), INTENT(INOUT) :: this
133 CLASS(mesh_base), INTENT(IN) :: Mesh
134 CLASS(physics_base), INTENT(IN) :: Physics
135 CLASS(timedisc_base), INTENT(IN) :: Timedisc
136 CLASS(sources_list), ALLOCATABLE, INTENT(IN) :: Sources
137 TYPE(dict_typ), INTENT(IN), POINTER :: config,IO
138 !------------------------------------------------------------------------!
139#ifdef PARALLEL
140 INTEGER, DIMENSION(3) :: gsizes,lsizes,indices
141#endif
142 CHARACTER(LEN=1),DIMENSION(:),POINTER :: mold
143 REAL :: r
144 INTEGER :: i,err
145 !------------------------------------------------------------------------!
146 IF (.NOT.this%Initialized()) &
147 CALL this%InitFileio(mesh,physics,timedisc,sources,config,io,"binary","bin",textfile=.false.)
148
149 ! We mark the endianess similar to the TIFF format
150 ! See: http://en.wikipedia.org/wiki/Endianness#Endianness_in_files_and_byte_swap
151 ! II == Intel Format == Little Endian
152 ! MM == Motorola Format == Big Endian
153 CALL this%GetEndianness(this%endianness, 'II','MM')
154 this%realsize = SIZE(transfer(r, mold))
155 this%intsize = SIZE(transfer(i, mold))
156
157 SELECT CASE(this%realsize)
158 CASE(4,8)
159 WRITE(this%realfmt,'(A,I1,A)') '"',this%realsize,'"'
160 CASE DEFAULT
161 CALL this%Error("fileio_binary::InitFileIO_binary","Only single and double precision are allowed")
162 END SELECT
163
164 ! set local array dimensions (frequently used)
165 this%INUM = mesh%IMAX-mesh%IMIN+1
166 this%JNUM = mesh%JMAX-mesh%JMIN+1
167 this%KNUM = mesh%KMAX-mesh%KMIN+1
168
169#ifdef PARALLEL
170 ! create the data type for the distributed array of
171 ! coordinates and simulation data
172 gsizes(1:3) = (/ mesh%INUM, mesh%JNUM, mesh%KNUM /)
173 lsizes(1:3) = (/ this%INUM, this%JNUM, this%KNUM /)
174 indices(1:3) = (/ mesh%IMIN-1, mesh%JMIN-1, mesh%KMIN-1 /)
175 this%bufsize = product(lsizes)
176 this%err = 0
177 SELECT TYPE(df=>this%datafile)
178 CLASS IS(filehandle_mpi)
179 CALL mpi_type_create_subarray(3, gsizes, lsizes, indices, mpi_order_fortran,&
180 default_mpi_real,df%filetype,this%err)
181 IF (this%err.EQ.0) CALL mpi_type_commit(df%filetype,this%err)
182 END SELECT
183
184 ! create the data type for the distributed array of
185 ! mesh corner positions
186 this%clsizes(:) = lsizes(:)
187 WHERE ((/mesh%INUM,mesh%JNUM,mesh%KNUM/).GT.1)
188 ! add +1 to the global size if dimension has real extend,i.e. not one cell
189 gsizes(:) = gsizes(:)+1
190 WHERE ((/mesh%INUM,mesh%JNUM,mesh%KNUM/).EQ.(/mesh%IMAX,mesh%JMAX,mesh%KMAX/))
191 ! add +1 to the local size if process contains the outermost domain
192 this%clsizes(:) = lsizes(:)+1
193 END WHERE
194 END WHERE
195 this%cbufsize = product(this%clsizes)
196 SELECT TYPE(df=>this%datafile)
197 CLASS IS(filehandle_mpi)
198 IF (this%err.EQ.0) CALL mpi_type_create_subarray(3, gsizes, this%clsizes, indices, mpi_order_fortran,&
199 default_mpi_real,df%cfiletype,this%err)
200 IF (this%err.EQ.0) CALL mpi_type_commit(df%cfiletype,this%err)
201 END SELECT
202
203 IF (this%err.NE.0) &
204 CALL this%Error("fileio_binary::InitFileIO_binary","creating MPI data types failed")
205 this%first = .true.
206#endif
207 END SUBROUTINE initfileio_binary
208
215 SUBROUTINE writeheader(this,Mesh,Physics,Header,IO)
216 IMPLICIT NONE
217 !------------------------------------------------------------------------!
218 CLASS(fileio_binary), INTENT(INOUT) :: this
219 CLASS(mesh_base), INTENT(IN) :: Mesh
220 CLASS(physics_base), INTENT(IN) :: Physics
221 TYPE(dict_typ), POINTER :: Header,IO
222 !------------------------------------------------------------------------!
223 CHARACTER(LEN=6) :: magic = "FOSITE"
224 CHARACTER(LEN=1) :: version = achar(0)
225 CHARACTER(LEN=4) :: sizes
226 CHARACTER(LEN=13) :: sheader
227 !------------------------------------------------------------------------!
228 WRITE(sizes, '(I2,I2)') this%realsize, this%intsize
229 sheader = magic // this%endianness(1:2) // version // sizes
230 this%offset = 0
231 SELECT TYPE(df=>this%datafile)
232#ifndef PARALLEL
233 CLASS IS(filehandle_fortran)
234 WRITE (unit=df%GetUnitNumber(),iostat=this%err) sheader
235#else
236 CLASS IS(filehandle_mpi)
237 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
238 mpi_byte, 'native', mpi_info_null, this%err)
239 IF (this%GetRank().EQ.0) &
240 CALL mpi_file_write(df%GetUnitNumber(),sheader,len(sheader),mpi_byte, &
241 df%status,this%err)
242#endif
243 CLASS DEFAULT
244 CALL this%Error("fileio_binary::WriteHeader","unknown file handle")
245 END SELECT
246 this%offset = this%offset + len(sheader)
247 END SUBROUTINE writeheader
248
249! !> \public Reads the header of a file (not yet implemented)
250! !!
251! SUBROUTINE ReadHeader_binary(this,success)
252! IMPLICIT NONE
253! !------------------------------------------------------------------------!
254! TYPE(FileIO_TYP) :: this !< \param [in,out] this fileio type
255! LOGICAL :: success !< \param [out] success
256! !------------------------------------------------------------------------!
257! INTENT(OUT) :: success
258! INTENT(INOUT) :: this
259! !------------------------------------------------------------------------!
260! success = .FALSE.
261! END SUBROUTINE ReadHeader_binary
262!
263! !> \public Writes the timestep (not yet implemented)
264! !!
265! SUBROUTINE WriteTimestamp_binary(this,time)
266! IMPLICIT NONE
267! !------------------------------------------------------------------------!
268! TYPE(FileIO_TYP) :: this !< \param [in,out] this fileio type
269! REAL :: time !< \param [in] time
270! !------------------------------------------------------------------------!
271! INTENT(IN) :: time
272! INTENT(INOUT) :: this
273! !------------------------------------------------------------------------!
274! END SUBROUTINE WriteTimestamp_binary
275!
276! !> \public Writes the timestep (not yet implemented)
277! !!
278! SUBROUTINE ReadTimestamp_binary(this,time)
279! IMPLICIT NONE
280! !------------------------------------------------------------------------!
281! TYPE(FileIO_TYP) :: this !< \param [in,out] this fileio type
282! REAL :: time !< \param [out] time
283! !------------------------------------------------------------------------!
284! INTENT(OUT) :: time
285! INTENT(INOUT) :: this
286! !------------------------------------------------------------------------!
287! END SUBROUTINE ReadTimestamp_binary
288
289
298 SUBROUTINE writekey(this,key,type,bytes,dims)
299 IMPLICIT NONE
300 !------------------------------------------------------------------------!
301 CLASS(fileio_binary), INTENT(INOUT) :: this
302 CHARACTER(LEN=*) :: key
303 INTEGER :: type,bytes
304 INTEGER, DIMENSION(5), OPTIONAL :: dims
305 INTEGER :: bufsize
306 CHARACTER(LEN=1), DIMENSION(:), ALLOCATABLE :: buf
307 !------------------------------------------------------------------------!
308 INTEGER :: keylen,l,b,o
309 !------------------------------------------------------------------------!
310 INTENT(IN) :: key,type,bytes
311 !------------------------------------------------------------------------!
312 keylen = len_trim(key)
313 IF(PRESENT(dims)) THEN
314 l = 3
315 IF(dims(4).GT.1) l = 4
316 IF(dims(5).GT.1) l = 5
317 ELSE
318 l = 0
319 END IF
320 b = bytes + l * this%intsize
321 bufsize = (3+l) * this%intsize + keylen
322 ALLOCATE(buf(bufsize))
323 o = 1
324 CALL append(buf,o,transfer(keylen,buf))
325 CALL append(buf,o,transfer(trim(key),buf))
326 CALL append(buf,o,transfer(type,buf))
327 CALL append(buf,o,transfer(b,buf))
328 IF(l.GT.0) THEN
329 CALL append(buf,o,transfer(dims(1:l),buf))
330 END IF
331
332 SELECT TYPE(df=>this%datafile)
333#ifndef PARALLEL
334 CLASS IS(filehandle_fortran)
335 WRITE (unit=df%GetUnitNumber(),iostat=this%err) buf
336#else
337 CLASS IS(filehandle_mpi)
338 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
339 mpi_byte, 'native', mpi_info_null, this%err)
340 IF (this%GetRank().EQ.0) &
341 CALL mpi_file_write(df%GetUnitNumber(),buf,bufsize,mpi_byte, &
342 df%status,this%err)
343#endif
344 CLASS DEFAULT
345 CALL this%Error("fileio_binary::WriteKey","unknown file handle")
346 END SELECT
347
348 DEALLOCATE(buf)
349 this%offset = this%offset + bufsize
350
351 CONTAINS
352 SUBROUTINE append(buffer,i,d)
353 IMPLICIT NONE
354 !--------------------------------------------------------------------!
355 CHARACTER(LEN=1), DIMENSION(:) :: buffer,d
356 INTEGER :: i,s
357 !--------------------------------------------------------------------!
358 s = SIZE(d)
359 buffer(i:i+s-1) = d
360 o = o + s
361 END SUBROUTINE
362 END SUBROUTINE writekey
363
364
367 RECURSIVE SUBROUTINE writedataattributes(this,Mesh,config,path)
368 IMPLICIT NONE
369 !------------------------------------------------------------------------!
370 CLASS(fileio_binary), INTENT(INOUT) :: this
371 CLASS(mesh_base), INTENT(IN) :: mesh
372 TYPE(dict_typ), POINTER :: config
373 CHARACTER(LEN=*), OPTIONAL :: path
374 !------------------------------------------------------------------------!
375 CHARACTER(LEN=MAX_CHAR_LEN) :: str, key
376 TYPE(dict_typ), POINTER :: node
377 !------------------------------------------------------------------------!
378 INTENT(IN) :: path
379 !------------------------------------------------------------------------!
380 IF(PRESENT(path)) THEN
381 str = path
382 ELSE
383 str = ""
384 ENDIF
385 node => config
386 DO WHILE(ASSOCIATED(node))
387 key = trim(str)//"/"//trim(getkey(node))
388
389 CALL this%WriteNode(mesh,key,node)
390
391 IF(haschild(node)) THEN
392 CALL this%WriteDataAttributes(mesh, getchild(node), key)
393 END IF
394 node=>getnext(node)
395 END DO
396 END SUBROUTINE writedataattributes
397
398 SUBROUTINE writenode(this,Mesh,key,node)
399 IMPLICIT NONE
400 !------------------------------------------------------------------------!
401 CLASS(fileio_binary), INTENT(INOUT) :: this
402 CLASS(mesh_base), INTENT(IN) :: Mesh
403 CHARACTER(LEN=MAX_CHAR_LEN) :: key
404 TYPE(dict_typ), POINTER :: node
405 !------------------------------------------------------------------------!
406 TYPE(real_t) :: ptr0
407 TYPE(int_t) :: ptrint
408 REAL, DIMENSION(:,:), POINTER, CONTIGUOUS :: ptr2
409 REAL, DIMENSION(:,:,:), POINTER :: ptr3
410 REAL, DIMENSION(:,:,:,:), POINTER :: ptr4
411 REAL, DIMENSION(:,:,:,:,:), POINTER :: ptr5
412 INTEGER, DIMENSION(5) :: dims
413 INTEGER :: bytes,k,l
414 CHARACTER(LEN=1), DIMENSION(:), POINTER :: val
415#ifdef PARALLEL
416 INTEGER(KIND=MPI_OFFSET_KIND) :: omax,omin
417#endif
418 !------------------------------------------------------------------------!
419 INTENT(IN) :: key
420 !------------------------------------------------------------------------!
421 dims(:) = 1
422 NULLIFY(ptr2,ptr3,ptr4,ptr5)
423
424 SELECT CASE(getdatatype(node))
425 CASE(dict_real_twod)
426 CALL getattr(node,key,ptr2)
427 dims(1:2) = shape(ptr2)
428 CASE(dict_real_threed)
429 CALL getattr(node,key,ptr3)
430 dims(1:3) = shape(ptr3)
431 CASE(dict_real_fourd)
432 CALL getattr(node,key,ptr4)
433 dims(1:4) = shape(ptr4)
434 CASE(dict_real_fived)
435 CALL getattr(node,key,ptr5)
436 dims(1:5) = shape(ptr5)
437 END SELECT
438
439 IF(product(dims).GT.1) THEN
440 CALL this%SetOutputDims(mesh,dims)
441 bytes = product(dims(1:3)) * this%realsize
442 CALL this%WriteKey(key,getdatatype(node),&
443 product(dims)*this%realsize,dims)
444 DO l=1,dims(5)
445 DO k=1,dims(4)
446 IF(ASSOCIATED(ptr4)) THEN
447 ptr3 => ptr4(:,:,:,k)
448 ELSE IF(ASSOCIATED(ptr5)) THEN
449 ptr3 => ptr5(:,:,:,k,l)
450 ELSE IF(ASSOCIATED(ptr2)) THEN
451 ptr3(1:dims(1),1:dims(2),1:1) => ptr2
452 END IF
453 SELECT TYPE(df=>this%datafile)
454#ifndef PARALLEL
455 CLASS IS(filehandle_fortran)
456 WRITE (unit=df%GetUnitNumber(),iostat=this%err) ptr3
457#else
458 CLASS IS(filehandle_mpi)
459 IF(this%HasMeshDims(mesh,shape(ptr3))) THEN
460 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,default_mpi_real,&
461 df%filetype, 'native', mpi_info_null, this%err)
462 CALL mpi_file_write_all(df%GetUnitNumber(),ptr3,this%bufsize,&
463 default_mpi_real,df%status,this%err)
464
465 ELSE IF(this%HasCornerDims(mesh,shape(ptr3))) THEN
466 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,default_mpi_real,&
467 df%cfiletype, 'native', mpi_info_null, this%err)
468 CALL mpi_file_write_all(df%GetUnitNumber(),ptr3(1:this%clsizes(1),1:this%clsizes(2),1:this%clsizes(3)),&
469 this%cbufsize,default_mpi_real,df%status,this%err)
470 ELSE
471 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
472 mpi_byte, 'native', mpi_info_null, this%err)
473 IF(this%GetRank().EQ.0) THEN
474 CALL mpi_file_write(df%GetUnitNumber(),ptr3,bytes,mpi_byte, &
475 df%status,this%err)
476 END IF
477 END IF
478#endif
479 CLASS DEFAULT
480 CALL this%Error("fileio_binary::WriteNode","unknown file handle")
481 END SELECT
482 this%offset = this%offset + bytes
483 END DO
484 END DO
485 ELSE
486 SELECT CASE(getdatatype(node))
487 CASE(dict_real_p)
488 CALL getattr(node,key,ptr0)
489 bytes = this%realsize
490 CALL this%WriteKey(key,dict_real,bytes)
491 SELECT TYPE(df=>this%datafile)
492#ifndef PARALLEL
493 CLASS IS(filehandle_fortran)
494 WRITE (unit=df%GetUnitNumber(),iostat=this%err) ptr0%p
495#else
496 CLASS IS(filehandle_mpi)
497 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
498 mpi_byte, 'native', mpi_info_null, this%err)
499 IF(this%GetRank().EQ.0) &
500 CALL mpi_file_write(df%GetUnitNumber(),ptr0%p,bytes,mpi_byte, &
501 df%status,this%err)
502#endif
503 CLASS DEFAULT
504 CALL this%Error("fileio_binary::WriteNode","unknown file handle")
505 END SELECT
506 CASE(dict_int_p)
507 CALL getattr(node,key,ptrint)
508 bytes = this%intsize
509 CALL this%WriteKey(key,dict_int,bytes)
510 SELECT TYPE(df=>this%datafile)
511#ifndef PARALLEL
512 CLASS IS(filehandle_fortran)
513 WRITE (unit=df%GetUnitNumber(),iostat=this%err) ptrint%p
514#else
515 CLASS IS(filehandle_mpi)
516 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
517 mpi_byte, 'native', mpi_info_null, this%err)
518 IF(this%GetRank().EQ.0) &
519 CALL mpi_file_write(df%GetUnitNumber(),ptrint%p,bytes,mpi_byte, &
520 df%status,this%err)
521#endif
522 CLASS DEFAULT
523 CALL this%Error("fileio_binary::WriteNode","unknown file handle")
524 END SELECT
525 CASE DEFAULT
526 val => getdata(node)
527 IF(ASSOCIATED(val)) THEN
528 bytes = SIZE(val)
529 ELSE
530 bytes = 0
531 END IF
532 CALL this%WriteKey(key,getdatatype(node),bytes)
533 IF(bytes.GT.0) THEN
534 SELECT TYPE(df=>this%datafile)
535#ifndef PARALLEL
536 CLASS IS(filehandle_fortran)
537 WRITE (unit=df%GetUnitNumber(),iostat=this%err) val
538#else
539 CLASS IS(filehandle_mpi)
540 CALL mpi_file_set_view(df%GetUnitNumber(),this%offset,mpi_byte,&
541 mpi_byte, 'native', mpi_info_null, this%err)
542 IF(this%GetRank().EQ.0) THEN
543 CALL mpi_file_write(df%GetUnitNumber(),val,bytes,mpi_byte, &
544 df%status,this%err)
545 END IF
546#endif
547 CLASS DEFAULT
548 CALL this%Error("fileio_binary::WriteNode","unknown file handle")
549 END SELECT
550 END IF
551 END SELECT
552 this%offset = this%offset + bytes
553 END IF
554#ifdef PARALLEL
555 IF(this%first) THEN
556 ! If MPI is used and this is the first output of the run,
557 ! it is checked, if the offsets of the different nodes are still
558 ! in sync. They can get easily out of sync if an output array is not
559 ! a mesh or corner array, but still has a different size on some nodes.
560 ! This easily happens, e.g. if one forgets to specify the subarray limits.
561 ! output_array(Mesh%IMIN:Mesh%IMAX,Mesh%JMIN:Mesh%JMAX)
562 ! The array would without them be interpreted as generic 2D array in
563 ! this case and therefore the calculated size could be different on
564 ! different nodes.
565 CALL mpi_allreduce(this%offset,omax,1,mpi_offset,mpi_max,&
566 mesh%comm_cart,this%err)
567 CALL mpi_allreduce(this%offset,omin,1,mpi_offset,mpi_min,&
568 mesh%comm_cart,this%err)
569 IF((this%offset.NE.omax).OR.(this%offset.NE.omin)) &
570 CALL this%Error("WriteNode_binary",&
571 "The offsets on different nodes are not in sync anymore." // lf &
572 //"The last key was '" // trim(key) // "'.")
573 END IF
574#endif
575 END SUBROUTINE writenode
576
577 FUNCTION hasmeshdims(this,Mesh,dims) RESULT(res)
578 IMPLICIT NONE
579 !------------------------------------------------------------------------!
580 CLASS(fileio_binary), INTENT(INOUT) :: this
581 CLASS(mesh_base), INTENT(IN) :: mesh
582 INTEGER, DIMENSION(:), INTENT(IN) :: dims
583 LOGICAL :: res
584 !------------------------------------------------------------------------!
585 IF(SIZE(dims).GE.3) THEN
586 res = all(dims(1:3).EQ.(/mesh%IMAX-mesh%IMIN+1, &
587 mesh%JMAX-mesh%JMIN+1, &
588 mesh%KMAX-mesh%KMIN+1/))
589 ELSE
590 res = .false.
591 END IF
592 END FUNCTION hasmeshdims
593
594
595 FUNCTION hascornerdims(this,Mesh,dims) RESULT(res)
596 IMPLICIT NONE
597 !------------------------------------------------------------------------!
598 CLASS(fileio_binary), INTENT(INOUT) :: this
599 CLASS(mesh_base), INTENT(IN) :: mesh
600 INTEGER, DIMENSION(:) :: dims
601 LOGICAL :: res
602 !------------------------------------------------------------------------!
603 INTENT(IN) :: dims
604 !------------------------------------------------------------------------!
605 IF(SIZE(dims).GE.3) THEN
606 ! remark: xP1 = 1 if xNUM > 1 and 0 if xNUM = 1 for x=I,J,K
607 res = all(dims(1:3).EQ.(/mesh%IMAX-mesh%IMIN+mesh%IP1+1, &
608 mesh%JMAX-mesh%JMIN+mesh%JP1+1, &
609 mesh%KMAX-mesh%KMIN+mesh%KP1+1/))
610 ELSE
611 res = .false.
612 END IF
613 END FUNCTION hascornerdims
614
615
616 SUBROUTINE setoutputdims(this,Mesh,dims)
617 IMPLICIT NONE
618 !------------------------------------------------------------------------!
619 CLASS(fileio_binary), INTENT(INOUT) :: this
620 CLASS(mesh_base), INTENT(IN) :: Mesh
621 INTEGER, DIMENSION(:), INTENT(INOUT):: dims
622 !------------------------------------------------------------------------!
623 IF(this%HasMeshDims(mesh,dims)) THEN
624 dims(1:3) = (/mesh%INUM,mesh%JNUM,mesh%KNUM/)
625 ELSE IF(this%HasCornerDims(mesh,dims)) THEN
626 dims(1:3) = (/mesh%INUM,mesh%JNUM,mesh%KNUM/)
627 WHERE ((/mesh%INUM,mesh%JNUM,mesh%KNUM/).GT.1)
628 ! add +1 to the global size to write one additional data point
629 ! in each dimension which has extend larger than 1
630 dims(1:3) = dims(1:3) + 1
631 END WHERE
632 END IF
633 END SUBROUTINE setoutputdims
634
635
638 SUBROUTINE writedataset_binary(this,Mesh,Physics,Fluxes,Timedisc,Header,IO)
639 IMPLICIT NONE
640 !------------------------------------------------------------------------!
641 CLASS(fileio_binary), INTENT(INOUT) :: this
642 CLASS(mesh_base), INTENT(IN) :: Mesh
643 CLASS(physics_base), INTENT(INOUT) :: Physics
644 CLASS(fluxes_base), INTENT(IN) :: Fluxes
645 CLASS(timedisc_base), INTENT(IN) :: Timedisc
646 TYPE(dict_typ), POINTER :: Header,IO
647 !------------------------------------------------------------------------!
648 CALL this%WriteDataAttributes(mesh,io)
649
650#ifdef PARALLEL
651 this%first = .false.
652#endif
653 END SUBROUTINE writedataset_binary
654
655! !> \public Reads the data arrays from file (not yet implemented)
656! !!
657! SUBROUTINE ReadDataset_binary(this,Mesh,Physics,Timedisc)
658! IMPLICIT NONE
659! !------------------------------------------------------------------------!
660! TYPE(FileIO_TYP) :: this !< \param [in,out] this fileio type
661! TYPE(Mesh_TYP) :: Mesh !< \param [in] mesh mesh type
662! TYPE(Physics_TYP) :: Physics !< \param [in] physics physics type
663! TYPE(Timedisc_TYP):: Timedisc !< \param [in,out] timedisc timedisc type
664! !------------------------------------------------------------------------!
665! INTENT(IN) :: Mesh,Physics
666! INTENT(INOUT) :: this,Timedisc
667! !------------------------------------------------------------------------!
668! END SUBROUTINE ReadDataset_binary
669
672 SUBROUTINE finalize(this)
673 IMPLICIT NONE
674 !------------------------------------------------------------------------!
675 TYPE(fileio_binary), INTENT(INOUT) :: this
676 !------------------------------------------------------------------------!
677#ifdef PARALLEL
678 SELECT TYPE(df=>this%datafile)
679 CLASS IS(filehandle_mpi)
680 CALL mpi_type_free(df%cfiletype,this%err)
681 CALL mpi_type_free(df%filetype,this%err)
682 END SELECT
683#endif
684 CALL this%Finalize_base()
685 END SUBROUTINE finalize
686
687END MODULE fileio_binary_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, parameter, public dict_int_p
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_twod
Definition: common_dict.f90:99
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
integer, parameter, public dict_int
Definition: common_dict.f90:94
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 append
read/write access at end
module for binary file I/O
recursive subroutine writedataattributes(this, Mesh, config, path)
Writes data attributes to a file.
subroutine initfileio_binary(this, Mesh, Physics, Timedisc, Sources, config, IO)
Constructor for the binary file I/O.
subroutine writekey(this, key, type, bytes, dims)
Writes key structure This subroutine writes the the key, data type and data sizes....
subroutine writeheader(this, Mesh, Physics, Header, IO)
Write the file header The header is written in ASCII and is 13 Byte long. First a "magic" identifier ...
subroutine setoutputdims(this, Mesh, dims)
character, parameter lf
line feed
subroutine finalize(this)
Closes the file I/O.
subroutine writenode(this, Mesh, key, node)
logical function hasmeshdims(this, Mesh, dims)
logical function hascornerdims(this, Mesh, dims)
subroutine writedataset_binary(this, Mesh, Physics, Fluxes, Timedisc, Header, IO)
Writes all desired data arrays to a file.
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
class for MPI file handle
mesh data structure
Definition: mesh_base.f90:122
container class to manage the list of source terms