output_vtk.F90 coverage: 80.00 %func 55.17 %block
1) module Output_VTK_module
2)
3) use Logging_module
4) use Output_Aux_module
5) use Output_Common_module
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) private
12)
13) #include "petsc/finclude/petscsys.h"
14)
15) PetscInt, parameter :: VTK_INTEGER = 0
16) PetscInt, parameter :: VTK_REAL = 1
17)
18) public :: OutputVTK
19)
20) contains
21)
22) ! ************************************************************************** !
23)
24) subroutine OutputVTK(realization_base)
25)
26) use Realization_Base_class, only : realization_base_type
27) use Discretization_module
28) use Grid_module
29) use Grid_Structured_module
30) use Option_module
31) use Field_module
32) use Patch_module
33) use String_module
34)
35) use Reaction_Aux_module
36) use Variables_module
37)
38) implicit none
39)
40) #include "petsc/finclude/petscvec.h"
41) #include "petsc/finclude/petscvec.h90"
42)
43) class(realization_base_type) :: realization_base
44)
45) PetscInt :: i, comma_count, quote_count
46) character(len=MAXSTRINGLENGTH) :: filename
47) character(len=MAXWORDLENGTH) :: word
48) character(len=MAXSTRINGLENGTH) :: string, string2
49) character(len=2) :: free_mol_char, tot_mol_char, sec_mol_char
50) type(grid_type), pointer :: grid
51) type(option_type), pointer :: option
52) type(discretization_type), pointer :: discretization
53) type(field_type), pointer :: field
54) type(patch_type), pointer :: patch
55) type(reaction_type), pointer :: reaction
56) type(output_option_type), pointer :: output_option
57) type(output_variable_type), pointer :: cur_variable
58) PetscReal, pointer :: vec_ptr(:)
59) Vec :: global_vec
60) Vec :: natural_vec
61) PetscErrorCode :: ierr
62)
63) discretization => realization_base%discretization
64) patch => realization_base%patch
65) grid => patch%grid
66) option => realization_base%option
67) field => realization_base%field
68) reaction => realization_base%reaction
69) output_option => realization_base%output_option
70)
71) ! open file
72) filename = OutputFilename(output_option,option,'vtk','')
73)
74) if (option%myrank == option%io_rank) then
75) option%io_buffer = '--> write vtk output file: ' // trim(filename)
76) call printMsg(option)
77) open(unit=OUTPUT_UNIT,file=filename,action="write")
78)
79) ! write header
80) write(OUTPUT_UNIT,'(''# vtk DataFile Version 2.0'')')
81) ! write title
82) write(OUTPUT_UNIT,'(''PFLOTRAN output'')')
83) write(OUTPUT_UNIT,'(''ASCII'')')
84) write(OUTPUT_UNIT,'(''DATASET UNSTRUCTURED_GRID'')')
85) endif
86)
87) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
88) option)
89) call DiscretizationCreateVector(discretization,ONEDOF,natural_vec,NATURAL, &
90) option)
91)
92) ! write out coordinates
93) call WriteVTKGrid(OUTPUT_UNIT,realization_base)
94)
95) if (option%myrank == option%io_rank) then
96) write(OUTPUT_UNIT,'(''CELL_DATA'',i8)') grid%nmax
97) endif
98)
99) cur_variable => output_option%output_snap_variable_list%first
100) do
101) if (.not.associated(cur_variable)) exit
102) call OutputGetVarFromArray(realization_base,global_vec,cur_variable%ivar, &
103) cur_variable%isubvar)
104) call DiscretizationGlobalToNatural(discretization,global_vec, &
105) natural_vec,ONEDOF)
106) word=trim(cur_variable%name)
107) call StringSwapChar(word," ","_")
108) if (cur_variable%iformat == 0) then
109) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base, &
110) word,natural_vec,VTK_REAL)
111) else
112) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base, &
113) word,natural_vec,VTK_INTEGER)
114) endif
115) cur_variable => cur_variable%next
116) enddo
117)
118) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
119) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
120)
121) if (option%myrank == option%io_rank) close(OUTPUT_UNIT)
122)
123) #if 1
124) if (output_option%print_vtk_vel_cent) then
125) call OutputVelocitiesVTK(realization_base)
126) endif
127) #endif
128)
129) #if 0
130) if (output_option%print_vtk_vel_cent) then
131) if (grid%structured_grid%nx > 1) then
132) call OutputFluxVelocitiesVTK(realization_base,LIQUID_PHASE, &
133) X_DIRECTION)
134) select case(option%iflowmode)
135) case(MPH_MODE,IMS_MODE,FLASH2_MODE,G_MODE)
136) call OutputFluxVelocitiesVTK(realization_base,GAS_PHASE, &
137) X_DIRECTION)
138) end select
139) endif
140) if (grid%structured_grid%ny > 1) then
141) call OutputFluxVelocitiesVTK(realization_base,LIQUID_PHASE, &
142) Y_DIRECTION)
143) select case(option%iflowmode)
144) case(MPH_MODE,IMS_MODE,FLASH2_MODE,G_MODE)
145) call OutputFluxVelocitiesVTK(realization_base,GAS_PHASE, &
146) Y_DIRECTION)
147) end select
148) endif
149) if (grid%structured_grid%nz > 1) then
150) call OutputFluxVelocitiesVTK(realization_base,LIQUID_PHASE, &
151) Z_DIRECTION)
152) select case(option%iflowmode)
153) case(MPH_MODE,IMS_MODE,FLASH2_MODE,G_MODE)
154) call OutputFluxVelocitiesVTK(realization_base,GAS_PHASE, &
155) Z_DIRECTION)
156) end select
157) endif
158) endif
159) #endif
160)
161) end subroutine OutputVTK
162)
163) #if 1
164)
165) ! ************************************************************************** !
166)
167) subroutine OutputVelocitiesVTK(realization_base)
168) !
169) ! Print velocities to Tecplot file in BLOCK format
170) !
171)
172) use Realization_Base_class, only : realization_base_type
173) use Discretization_module
174) use Grid_module
175) use Option_module
176) use Field_module
177) use Patch_module
178) use Variables_module
179)
180) implicit none
181)
182) #include "petsc/finclude/petscvec.h"
183) #include "petsc/finclude/petscvec.h90"
184)
185) class(realization_base_type) :: realization_base
186)
187) type(grid_type), pointer :: grid
188) type(option_type), pointer :: option
189) type(field_type), pointer :: field
190) type(discretization_type), pointer :: discretization
191) type(patch_type), pointer :: patch
192) type(output_option_type), pointer :: output_option
193) character(len=MAXSTRINGLENGTH) :: filename
194) character(len=MAXSTRINGLENGTH) :: string
195) character(len=MAXWORDLENGTH) :: word
196) Vec :: global_vec
197) Vec :: natural_vec
198) Vec :: global_vec_vx,global_vec_vy,global_vec_vz
199) PetscErrorCode :: ierr
200)
201) PetscReal, pointer :: vec_ptr(:)
202)
203) patch => realization_base%patch
204) grid => patch%grid
205) field => realization_base%field
206) option => realization_base%option
207) output_option => realization_base%output_option
208) discretization => realization_base%discretization
209)
210) ! open file
211) filename = OutputFilename(output_option,option,'vtk','vel')
212)
213) if (option%myrank == option%io_rank) then
214) option%io_buffer = '--> write vtk velocity output file: ' // &
215) trim(filename)
216) call printMsg(option)
217) open(unit=OUTPUT_UNIT,file=filename,action="write")
218)
219) ! write header
220) write(OUTPUT_UNIT,'(''# vtk DataFile Version 2.0'')')
221) ! write title
222) write(OUTPUT_UNIT,'(''PFLOTRAN output'')')
223) write(OUTPUT_UNIT,'(''ASCII'')')
224) write(OUTPUT_UNIT,'(''DATASET UNSTRUCTURED_GRID'')')
225)
226) endif
227)
228) ! write blocks
229) ! write out data sets
230) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
231) option)
232) call DiscretizationCreateVector(discretization,ONEDOF,natural_vec,NATURAL, &
233) option)
234) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vx)
235) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vy)
236) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vz)
237)
238) ! write out coordinates
239) call WriteVTKGrid(OUTPUT_UNIT,realization_base)
240)
241) if (option%myrank == option%io_rank) then
242) write(OUTPUT_UNIT,'(''CELL_DATA'',i8)') grid%nmax
243) endif
244)
245) word = 'Vlx'
246) call OutputGetCellCenteredVelocities(realization_base,global_vec_vx, &
247) global_vec_vy,global_vec_vz,LIQUID_PHASE)
248)
249) call DiscretizationGlobalToNatural(discretization,global_vec_vx,natural_vec,ONEDOF)
250) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
251)
252) word = 'Vly'
253) call DiscretizationGlobalToNatural(discretization,global_vec_vy,natural_vec,ONEDOF)
254) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
255)
256) word = 'Vlz'
257) call DiscretizationGlobalToNatural(discretization,global_vec_vz,natural_vec,ONEDOF)
258) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
259)
260) if (option%nphase > 1) then
261) word = 'Vgx'
262) call OutputGetCellCenteredVelocities(realization_base,global_vec_vx, &
263) global_vec_vy,global_vec_vz,GAS_PHASE)
264) call DiscretizationGlobalToNatural(discretization,global_vec_vx,natural_vec,ONEDOF)
265) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
266)
267) word = 'Vgy'
268) call DiscretizationGlobalToNatural(discretization,global_vec_vy,natural_vec,ONEDOF)
269) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
270)
271) word = 'Vgz'
272) call DiscretizationGlobalToNatural(discretization,global_vec_vz,natural_vec,ONEDOF)
273) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_REAL)
274) endif
275)
276) ! material id
277) word = 'Material_ID'
278) call OutputGetVarFromArray(realization_base,global_vec,MATERIAL_ID,ZERO_INTEGER)
279) call DiscretizationGlobalToNatural(discretization,global_vec,natural_vec,ONEDOF)
280) call WriteVTKDataSetFromVec(OUTPUT_UNIT,realization_base,word,natural_vec,VTK_INTEGER)
281)
282) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
283) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
284) call VecDestroy(global_vec_vx,ierr);CHKERRQ(ierr)
285) call VecDestroy(global_vec_vy,ierr);CHKERRQ(ierr)
286) call VecDestroy(global_vec_vz,ierr);CHKERRQ(ierr)
287)
288) if (option%myrank == option%io_rank) close(OUTPUT_UNIT)
289)
290) end subroutine OutputVelocitiesVTK
291) #endif
292)
293) ! ************************************************************************** !
294)
295) subroutine WriteVTKGrid(fid,realization_base)
296) !
297) ! Writes a grid in VTK format
298) !
299)
300) use Realization_Base_class, only : realization_base_type
301) use Discretization_module
302) use Grid_module
303) use Option_module
304) use Patch_module
305)
306) implicit none
307)
308) #include "petsc/finclude/petscvec.h"
309) #include "petsc/finclude/petscvec.h90"
310)
311) PetscInt :: fid
312) class(realization_base_type) :: realization_base
313)
314) type(discretization_type), pointer :: discretization
315) type(grid_type), pointer :: grid
316) type(option_type), pointer :: option
317) type(patch_type), pointer :: patch
318) PetscInt :: i, j, k, nx, ny, nz
319) PetscReal :: x, y, z
320) PetscInt :: nxp1Xnyp1, nxp1, nyp1, nzp1
321) PetscInt :: vertex_id
322) PetscErrorCode :: ierr
323)
324) 1000 format(es13.6,1x,es13.6,1x,es13.6)
325) 1001 format(i1,8(1x,i8))
326)
327) call PetscLogEventBegin(logging%event_output_grid_vtk,ierr);CHKERRQ(ierr)
328)
329) discretization => realization_base%discretization
330) patch => realization_base%patch
331) grid => patch%grid
332) option => realization_base%option
333)
334) if (realization_base%discretization%itype == STRUCTURED_GRID) then
335)
336) nx = grid%structured_grid%nx
337) ny = grid%structured_grid%ny
338) nz = grid%structured_grid%nz
339)
340) nxp1 = nx+1
341) nyp1 = ny+1
342) nzp1 = nz+1
343)
344) if (option%myrank == option%io_rank) then
345)
346) 1010 format("POINTS",1x,i12,1x,"float")
347) write(fid,1010) (nx+1)*(ny+1)*(nz+1)
348) do k=0,nz
349) if (k > 0) then
350) z = z + grid%structured_grid%dz_global(k)
351) else
352) z = discretization%origin_global(Z_DIRECTION)
353) endif
354) do j=0,ny
355) if (j > 0) then
356) y = y + grid%structured_grid%dy_global(j)
357) else
358) y = discretization%origin_global(Y_DIRECTION)
359) endif
360) x = discretization%origin_global(X_DIRECTION)
361) write(fid,1000) x,y,z
362) do i=1,nx
363) x = x + grid%structured_grid%dx_global(i)
364) write(fid,1000) x,y,z
365) enddo
366) enddo
367) enddo
368)
369) 1020 format('CELLS',1x,i12,1x,i12)
370) write(fid,1020) grid%nmax, grid%nmax*9
371) nxp1Xnyp1 = nxp1*nyp1
372) do k=0,nz-1
373) do j=0,ny-1
374) do i=0,nx-1
375) vertex_id = i+j*nxp1+k*nxp1Xnyp1
376) write(fid,1001) 8,vertex_id,vertex_id+1, &
377) vertex_id+nxp1+1,vertex_id+nxp1, &
378) vertex_id+nxp1Xnyp1,vertex_id+nxp1Xnyp1+1, &
379) vertex_id+nxp1Xnyp1+nxp1+1, &
380) vertex_id+nxp1Xnyp1+nxp1
381) enddo
382) enddo
383) enddo
384)
385) write(fid,'(a)') ""
386)
387) 1030 format('CELL_TYPES',1x,i12)
388) write(fid,1030) grid%nmax
389) do i=1,grid%nmax
390) write(fid,'(i2)') 12
391) enddo
392)
393) write(fid,'(a)') ""
394)
395) endif
396) endif
397)
398) call PetscLogEventEnd(logging%event_output_grid_vtk,ierr);CHKERRQ(ierr)
399)
400) end subroutine WriteVTKGrid
401)
402) ! ************************************************************************** !
403)
404) subroutine WriteVTKDataSetFromVec(fid,realization_base,dataset_name,vec,datatype)
405) !
406) ! Writes data from a Petsc Vec within a block
407) ! of a VTK file
408) !
409)
410) use Realization_Base_class, only : realization_base_type
411)
412) implicit none
413) #include "petsc/finclude/petscvec.h"
414) #include "petsc/finclude/petscvec.h90"
415)
416) PetscInt :: fid
417) class(realization_base_type) :: realization_base
418) Vec :: vec
419) character(len=MAXWORDLENGTH) :: dataset_name
420) PetscInt :: datatype
421) PetscErrorCode :: ierr
422)
423) PetscReal, pointer :: vec_ptr(:)
424)
425) call VecGetArrayF90(vec,vec_ptr,ierr);CHKERRQ(ierr)
426) call WriteVTKDataSet(fid,realization_base,dataset_name,vec_ptr,datatype, &
427) ZERO_INTEGER) ! 0 implies grid%nlmax
428) call VecRestoreArrayF90(vec,vec_ptr,ierr);CHKERRQ(ierr)
429)
430) end subroutine WriteVTKDataSetFromVec
431)
432) ! ************************************************************************** !
433)
434) subroutine WriteVTKDataSet(fid,realization_base,dataset_name,array,datatype, &
435) size_flag)
436) !
437) ! Writes data from an array within a block
438) ! of a VTK file
439) !
440)
441) use Realization_Base_class, only : realization_base_type
442) use Grid_module
443) use Option_module
444) use Patch_module
445)
446) implicit none
447)
448) PetscInt :: fid
449) class(realization_base_type) :: realization_base
450) PetscReal :: array(:)
451) character(len=MAXWORDLENGTH) :: dataset_name
452) PetscInt :: datatype
453) PetscInt :: size_flag ! if size_flag /= 0, use size_flag as the local size
454)
455) type(grid_type), pointer :: grid
456) type(option_type), pointer :: option
457) type(patch_type), pointer :: patch
458) PetscInt :: i
459) PetscInt :: max_proc, max_proc_prefetch
460) PetscMPIInt :: iproc_mpi, recv_size_mpi
461) PetscInt :: max_local_size
462) PetscMPIInt :: local_size_mpi
463) PetscInt :: istart, iend, num_in_array
464) PetscMPIInt :: status_mpi(MPI_STATUS_SIZE)
465) PetscInt, allocatable :: integer_data(:), integer_data_recv(:)
466) PetscReal, allocatable :: real_data(:), real_data_recv(:)
467) PetscErrorCode :: ierr
468)
469) 1001 format(10(es13.6,1x))
470) 1002 format(i3)
471)
472) patch => realization_base%patch
473) grid => patch%grid
474) option => realization_base%option
475)
476) call PetscLogEventBegin(logging%event_output_write_vtk,ierr);CHKERRQ(ierr)
477)
478) ! maximum number of initial messages
479) #define HANDSHAKE
480) max_proc = option%io_handshake_buffer_size
481) max_proc_prefetch = option%io_handshake_buffer_size / 10
482)
483) if (size_flag /= 0) then
484) call MPI_Allreduce(size_flag,max_local_size,ONE_INTEGER_MPI,MPIU_INTEGER, &
485) MPI_MAX,option%mycomm,ierr)
486) local_size_mpi = size_flag
487) else
488) ! if first time, determine the maximum size of any local array across
489) ! all procs
490) if (max_local_size_saved < 0) then
491) call MPI_Allreduce(grid%nlmax,max_local_size,ONE_INTEGER_MPI, &
492) MPIU_INTEGER,MPI_MAX,option%mycomm,ierr)
493) max_local_size_saved = max_local_size
494) if (OptionPrintToScreen(option)) print *, 'max_local_size_saved: ', &
495) max_local_size
496) endif
497) max_local_size = max_local_size_saved
498) local_size_mpi = grid%nlmax
499) endif
500)
501) ! transfer the data to an integer or real array
502) if (datatype == VTK_INTEGER) then
503) allocate(integer_data(max_local_size+10))
504) allocate(integer_data_recv(max_local_size))
505) do i=1,local_size_mpi
506) integer_data(i) = int(array(i))
507) enddo
508) else
509) allocate(real_data(max_local_size+10))
510) allocate(real_data_recv(max_local_size))
511) do i=1,local_size_mpi
512) real_data(i) = array(i)
513) enddo
514) endif
515)
516) ! communicate data to processor 0, round robin style
517) if (option%myrank == option%io_rank) then
518)
519) if (datatype == VTK_INTEGER) then
520) write(fid,'(''SCALARS '',a20,'' int 1'')') dataset_name
521) else
522) write(fid,'(''SCALARS '',a20,'' float 1'')') dataset_name
523) endif
524)
525) write(fid,'(''LOOKUP_TABLE default'')')
526)
527) if (datatype == VTK_INTEGER) then
528) ! This approach makes output files identical, regardless of processor
529) ! distribution. It is necessary when diffing files.
530) iend = 0
531) do
532) istart = iend+1
533) if (iend+10 > local_size_mpi) exit
534) iend = istart+9
535) write(fid,1002) integer_data(istart:iend)
536) enddo
537) ! shift remaining data to front of array
538) integer_data(1:local_size_mpi-iend) = integer_data(iend+1:local_size_mpi)
539) num_in_array = local_size_mpi-iend
540) else
541) iend = 0
542) do
543) istart = iend+1
544) if (iend+10 > local_size_mpi) exit
545) iend = istart+9
546) write(fid,1001) real_data(istart:iend)
547) enddo
548) ! shift remaining data to front of array
549) real_data(1:local_size_mpi-iend) = real_data(iend+1:local_size_mpi)
550) num_in_array = local_size_mpi-iend
551) endif
552) do iproc_mpi=1,option%mycommsize-1
553) #ifdef HANDSHAKE
554) if (option%io_handshake_buffer_size > 0 .and. &
555) iproc_mpi+max_proc_prefetch >= max_proc) then
556) max_proc = max_proc + option%io_handshake_buffer_size
557) call MPI_Bcast(max_proc,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
558) option%mycomm,ierr)
559) endif
560) #endif
561) call MPI_Probe(iproc_mpi,MPI_ANY_TAG,option%mycomm,status_mpi,ierr)
562) recv_size_mpi = status_mpi(MPI_TAG)
563) if (datatype == 0) then
564) call MPI_Recv(integer_data_recv,recv_size_mpi,MPIU_INTEGER,iproc_mpi, &
565) MPI_ANY_TAG,option%mycomm,status_mpi,ierr)
566) if (recv_size_mpi > 0) then
567) integer_data(num_in_array+1:num_in_array+recv_size_mpi) = &
568) integer_data_recv(1:recv_size_mpi)
569) num_in_array = num_in_array+recv_size_mpi
570) endif
571) iend = 0
572) do
573) istart = iend+1
574) if (iend+10 > num_in_array) exit
575) iend = istart+9
576) write(fid,1002) integer_data(istart:iend)
577) enddo
578) if (iend > 0) then
579) integer_data(1:num_in_array-iend) = integer_data(iend+1:num_in_array)
580) num_in_array = num_in_array-iend
581) endif
582) else
583) call MPI_Recv(real_data_recv,recv_size_mpi,MPI_DOUBLE_PRECISION,iproc_mpi, &
584) MPI_ANY_TAG,option%mycomm,status_mpi,ierr)
585) if (recv_size_mpi > 0) then
586) real_data(num_in_array+1:num_in_array+recv_size_mpi) = &
587) real_data_recv(1:recv_size_mpi)
588) num_in_array = num_in_array+recv_size_mpi
589) endif
590) iend = 0
591) do
592) istart = iend+1
593) if (iend+10 > num_in_array) exit
594) iend = istart+9
595) write(fid,1001) real_data(istart:iend)
596) enddo
597) if (iend > 0) then
598) real_data(1:num_in_array-iend) = real_data(iend+1:num_in_array)
599) num_in_array = num_in_array-iend
600) endif
601) endif
602) enddo
603) #ifdef HANDSHAKE
604) if (option%io_handshake_buffer_size > 0) then
605) max_proc = -1
606) call MPI_Bcast(max_proc,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
607) option%mycomm,ierr)
608) endif
609) #endif
610) ! Print the remaining values, if they exist
611) if (datatype == 0) then
612) if (num_in_array > 0) &
613) write(fid,1002) integer_data(1:num_in_array)
614) else
615) if (num_in_array > 0) &
616) write(fid,1001) real_data(1:num_in_array)
617) endif
618) write(fid,'(/)')
619) else
620) #ifdef HANDSHAKE
621) if (option%io_handshake_buffer_size > 0) then
622) do
623) if (option%myrank < max_proc) exit
624) call MPI_Bcast(max_proc,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
625) option%mycomm,ierr)
626) enddo
627) endif
628) #endif
629) if (datatype == VTK_INTEGER) then
630) call MPI_Send(integer_data,local_size_mpi,MPIU_INTEGER,option%io_rank, &
631) local_size_mpi,option%mycomm,ierr)
632) else
633) call MPI_Send(real_data,local_size_mpi,MPI_DOUBLE_PRECISION,option%io_rank, &
634) local_size_mpi,option%mycomm,ierr)
635) endif
636) #ifdef HANDSHAKE
637) if (option%io_handshake_buffer_size > 0) then
638) do
639) call MPI_Bcast(max_proc,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
640) option%mycomm,ierr)
641) if (max_proc < 0) exit
642) enddo
643) endif
644) #endif
645) #undef HANDSHAKE
646) endif
647)
648) if (datatype == VTK_INTEGER) then
649) deallocate(integer_data)
650) else
651) deallocate(real_data)
652) endif
653)
654) call PetscLogEventEnd(logging%event_output_write_vtk,ierr);CHKERRQ(ierr)
655)
656) end subroutine WriteVTKDataSet
657)
658) end module Output_VTK_module