MFIX  2016-1
vtk_out.f
Go to the documentation of this file.
1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2 ! C
3 ! Module name: WRITE_DBG_VTU_AND_VTP_FILES C
4 ! Purpose: Writes the cell and particle data in VTK format C
5 ! for debug regions only. C
6 ! C
7 ! Author: Jeff Dietiker Date: 22-Jul-15 C
8 ! Reviewer: Date: C
9 ! C
10 ! Revision Number # Date: ##-###-## C
11 ! Author: # C
12 ! Purpose: # C
13 ! C
14 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
16 
17  use discretelement, only: discrete_element
18  use vtp, only: write_vtp_file
19  use vtk, only: dimension_vtk
20 
21  IMPLICIT NONE
22  INTEGER :: LC
23 
24  DO lc = 1, dimension_vtk
25  CALL write_vtu_file(lc,1)
26  IF(discrete_element) CALL write_vtp_file(lc,1)
27  ENDDO
28 
29  END SUBROUTINE write_dbg_vtu_and_vtp_files
30 
31 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
32 ! C
33 ! Module name: WRITE_VTU_FILE C
34 ! Purpose: Writes the cell data grid in VTK format (Unstructured VTU) C
35 ! Binary format C
36 ! C
37 ! Author: Jeff Dietiker Date: 21-Feb-08 C
38 ! Reviewer: Date: C
39 ! C
40 ! Revision Number # Date: ##-###-## C
41 ! Author: # C
42 ! Purpose: # C
43 ! C
44 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
45  SUBROUTINE write_vtu_file(LCV,MODE)
46 
47  USE compar
48  USE constant
49  USE cutcell
50  USE discretelement, Only : discrete_element
51  USE fldvar
52  USE functions
53  USE geometry
54  USE indices
55  USE mfix_pic
56  USE mpi_utility
57  USE output
58  USE parallel
59  USE parallel_mpi
60  USE param
61  USE param1
62  USE pgcor
63  USE pgcor
64  USE physprop
65  USE pscor
66  USE quadric
67  USE run
68  USE rxns
69  USE scalars
70  USE sendrecv
71  USE stl
72  USE toleranc
73  USE visc_s
74  USE vtk
75  USE vtp
76 
77  IMPLICIT NONE
78  INTEGER :: I,J,K,L,M,NN,R,IJK,LCV
79 
80  DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: FACET_COUNT_DES, NEIGHBORING_FACET
81 
82  INTEGER :: SPECIES_COUNTER,LT
83 
84  CHARACTER (LEN=32) :: SUBM,SUBN,SUBR
85  CHARACTER (LEN=64) :: VAR_NAME
86 
87  DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: DP_BC_ID, IJK_ARRAY
88 
89  INTEGER :: PASS
90  INTEGER :: WRITE_HEADER = 1
91  INTEGER :: WRITE_DATA = 2
92  INTEGER :: MODE ! MODE = 0 : Write regular VTK region file
93  ! MODE = 1 : Write debug VTK region file (VTK_DBG_FILE = .TRUE.)
94 
95  ! There is nothing to write if we are not in adefined vtk region
96  vtk_region = lcv
97  IF(.NOT.vtk_defined(vtk_region)) RETURN
98  IF(vtk_data(lcv)/='C') RETURN
99  IF(mode==0.AND.(vtk_dbg_file(lcv))) RETURN
100  IF(mode==1.AND.(.NOT.vtk_dbg_file(lcv))) RETURN
101 
102 ! Location of U-momentum cells for original (uncut grid)
103  IF (do_i) THEN
104  xg_e(1) = zero
105  DO i = imin1, imax2
106  xg_e(i) = xg_e(i-1) + dx(i)
107  END DO
108  ENDIF
109 
110 ! Location of V-momentum cells for original (uncut grid)
111  IF (do_j) THEN
112  yg_n(1) = zero
113  DO j = jmin1, jmax2
114  yg_n(j) = yg_n(j-1) + dy(j)
115  END DO
116  ENDIF
117 
118 ! Location of W-momentum cells for original (uncut grid)
119  IF (do_k) THEN
120  zg_t(1) = zero
121  DO k = kmin1, kmax2
122  zg_t(k) = zg_t(k-1) + dz(k)
123  END DO
124  ELSE
125  zg_t = zero
126  ENDIF
127 
128 
129  CALL setup_vtk_region
130 
131  CALL open_vtu_file_bin(mode)
132 
133  IF(mode==0) CALL open_pvd_file
134 
135  CALL clean_geometry
136 
137  DO pass=write_header,write_data
138 
139 
140  CALL write_geometry_in_vtu_bin(pass)
141 
142  IF(vtk_ep_g(vtk_region)) &
143  CALL write_scalar_in_vtu_bin('EP_G',ep_g,pass)
144 
145  IF(vtk_p_g(vtk_region)) &
146  CALL write_scalar_in_vtu_bin('P_G',p_g,pass)
147 
148  IF(vtk_p_star(vtk_region)) &
149  CALL write_scalar_in_vtu_bin('P_S',p_s,pass)
150 
151  IF(vtk_vel_g(vtk_region)) &
152  CALL write_vector_in_vtu_bin('Gas_Velocity',u_g,v_g,w_g,pass)
153 
154  IF(vtk_u_g(vtk_region)) &
155  CALL write_scalar_in_vtu_bin('U_G',p_g,pass)
156 
157  IF(vtk_v_g(vtk_region)) &
158  CALL write_scalar_in_vtu_bin('V_G',p_g,pass)
159 
160  IF(vtk_w_g(vtk_region)) &
161  CALL write_scalar_in_vtu_bin('W_G',p_g,pass)
162 
163  DO m = 1,mmax
164  IF(vtk_vel_s(vtk_region,m)) THEN
165  WRITE(subm,*)m
166  CALL write_vector_in_vtu_bin('Solids_Velocity_'//adjustl(subm),u_s(:,m),v_s(:,m),w_s(:,m),pass)
167  ENDIF
168  IF(vtk_u_s(vtk_region,m)) THEN
169  WRITE(subm,*)m
170  CALL write_scalar_in_vtu_bin('U_s_'//adjustl(subm),u_s(:,m),pass)
171  ENDIF
172  IF(vtk_v_s(vtk_region,m)) THEN
173  WRITE(subm,*)m
174  CALL write_scalar_in_vtu_bin('V_s_'//adjustl(subm),v_s(:,m),pass)
175  ENDIF
176  IF(vtk_w_s(vtk_region,m)) THEN
177  WRITE(subm,*)m
178  CALL write_scalar_in_vtu_bin('W_s_'//adjustl(subm),w_s(:,m),pass)
179  ENDIF
180  END DO
181 
182  DO m = 1,mmax
183  IF(vtk_rop_s(vtk_region,m)) THEN
184  WRITE(subm,*)m
185  CALL write_scalar_in_vtu_bin('Solids_density_'//adjustl(subm),rop_s(:,m),pass)
186  ENDIF
187  END DO
188 
189  IF(vtk_t_g(vtk_region)) &
190  CALL write_scalar_in_vtu_bin('Gas_temperature',t_g,pass)
191 
192  DO m = 1,mmax
193  IF(vtk_t_s(vtk_region,m)) THEN
194  WRITE(subm,*)m
195  CALL write_scalar_in_vtu_bin('Solids_temperature_'//adjustl(subm),t_s(:,m),pass)
196  ENDIF
197  END DO
198 
199 
200  species_counter = 0
201  DO nn = 1,nmax(0)
202  IF(vtk_x_g(vtk_region,nn)) THEN
203  WRITE(subn,*)nn
204  IF(use_rrates) THEN
205  species_counter = species_counter + 1
206  var_name = adjustl(species_name(species_counter))
207  lt = len_trim(adjustl(species_name(species_counter)))
208  ELSE
209  var_name = adjustl(species_alias_g(nn))
210  lt = len_trim(adjustl(species_alias_g(nn)))
211  ENDIF
212  var_name = var_name(1:lt)//'_Gas_mass_fractions_'//adjustl(subn)
213  CALL write_scalar_in_vtu_bin(var_name,x_g(:,nn),pass)
214  ENDIF
215  END DO
216 
217  DO m = 1, mmax
218  WRITE(subm,*)m
219  DO nn = 1,nmax(m)
220  IF(vtk_x_s(vtk_region,m,nn)) THEN
221  WRITE(subn,*)nn
222  IF(use_rrates) THEN
223  species_counter = species_counter + 1
224  var_name = adjustl(species_name(species_counter))
225  lt = len_trim(adjustl(species_name(species_counter)))
226  ELSE
227  var_name = adjustl(species_alias_s(m,nn))
228  lt = len_trim(adjustl(species_alias_s(m,nn)))
229  ENDIF
230  var_name = var_name(1:lt)//'_Solids_mass_fractions_'//trim(adjustl(subm))//'_'//adjustl(subn)
231  CALL write_scalar_in_vtu_bin(var_name,x_s(:,m,nn),pass)
232  ENDIF
233  END DO
234  END DO
235 
236  DO m = 1,mmax
237  IF(vtk_theta_m(vtk_region,m)) THEN
238  WRITE(subm,*)m
239  CALL write_scalar_in_vtu_bin('Granular_temperature_'//adjustl(subm),theta_m(:,m),pass)
240  ENDIF
241  END DO
242 
243  DO nn = 1,nscalar
244  IF(vtk_scalar(vtk_region,nn)) THEN
245  WRITE(subn,*)nn
246  var_name = 'Scalar_'//adjustl(subn)
247  CALL write_scalar_in_vtu_bin(var_name,scalar(:,nn),pass)
248  ENDIF
249  END DO
250 
251  DO r = 1,nrr
252  IF(vtk_rrate(vtk_region,r)) THEN
253  WRITE(subr,*)r
254  var_name = 'RRates_'//adjustl(subr)
255  CALL write_scalar_in_vtu_bin(var_name,reactionrates(:, r),pass)
256  ENDIF
257  END DO
258 
259  IF(k_epsilon) THEN
260  IF(vtk_k_turb_g(vtk_region)) &
261  CALL write_scalar_in_vtu_bin('K_Turb_G',k_turb_g,pass)
262  IF(vtk_e_turb_g(vtk_region)) &
263  CALL write_scalar_in_vtu_bin('E_Turb_G',e_turb_g,pass)
264  ENDIF
265 
266 
268  CALL calc_vorticity
269  ENDIF
270 
272  CALL write_scalar_in_vtu_bin('VORTICITY_MAG',vorticity,pass)
273  IF(vtk_lambda_2(vtk_region)) &
274  CALL write_scalar_in_vtu_bin('LAMBDA_2',lambda2,pass)
275 
276 
278  CALL write_scalar_in_vtu_bin('PARTITION',partition,pass)
279 
280 
281  IF(vtk_bc_id(vtk_region)) THEN
282  Allocate(dp_bc_id(dimension_3))
283  dp_bc_id = dble(bc_id)
284  CALL write_scalar_in_vtu_bin('BC_ID',dp_bc_id,pass)
285  DeAllocate(dp_bc_id)
286  ENDIF
287 
288 
289  IF(vtk_dwall(vtk_region)) &
290  CALL write_scalar_in_vtu_bin('DISTANCE_TO_WALL',dwall,pass)
291 
292  IF(vtk_ijk(vtk_region)) THEN
293  Allocate(ijk_array(dimension_3))
294  DO ijk = ijkstart3, ijkend3
295  ijk_array(ijk) = dble(ijk)
296  ENDDO
297  CALL write_scalar_in_vtu_bin('IJK',ijk_array,pass)
298  DeAllocate(ijk_array)
299  ENDIF
300 
301  IF(vtk_ijk(vtk_region)) &
302  CALL write_vector_in_vtu_bin('Scalar normal',normal_s(:,1),normal_s(:,2),normal_s(:,3),pass)
303 
304  DO nn=1,15
305  IF(vtk_debug(vtk_region,nn)) THEN
306  WRITE(subn,*)nn
307  var_name = 'DEBUG_'//adjustl(subn)
308  CALL write_scalar_in_vtu_bin(var_name,debug_cg(:,nn),pass)
309  ENDIF
310  ENDDO
311 
312 
313 
314  ENDDO ! PASS LOOP, EITHER HEADER OR DATA
315 
316 
317  CALL close_vtu_file_bin(mode)
318  IF(mode==0) CALL update_and_close_pvd_file
319 
320 #ifdef MPI
321  call mpi_barrier(mpi_comm_world,mpierr)
322 #endif
323 
324 ! Update Frames
325  IF (mype == pe_io.AND.time_dependent_filename) THEN
326  OPEN(convert='BIG_ENDIAN',unit = vtu_frame_unit, file = trim(vtu_frame_filename))
327  DO l = 1, dimension_vtk
328  IF(vtk_defined(l)) WRITE(vtu_frame_unit,*) l,frame(l)
329  ENDDO
330  CLOSE(vtu_frame_unit)
331  ENDIF
332 
333  IF (full_log.AND.mype == pe_io) WRITE(*,20)' DONE.'
334 
335 20 FORMAT(a,1x/)
336  RETURN
337 
338  END SUBROUTINE write_vtu_file
339 
340 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
341 ! C
342 ! Module name: OPEN_VTU_FILE C
343 ! Purpose: Open a vtu file and writes the header C
344 ! Binary format C
345 ! C
346 ! Author: Jeff Dietiker Date: 21-Feb-08 C
347 ! Reviewer: Date: C
348 ! C
349 ! Revision Number # Date: ##-###-## C
350 ! Author: # C
351 ! Purpose: # C
352 ! C
353 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
354  SUBROUTINE open_vtu_file_bin(MODE)
356  USE cdist
357  USE compar
358  USE constant
359  USE cutcell
360  USE exit, only: mfix_exit
361  USE fldvar
362  USE functions
363  USE geometry
364  USE indices
365  USE mpi_utility
366  USE output
367  USE parallel
368  USE param
369  USE param1
370  USE quadric
371  USE run
372  USE sendrecv
373  USE toleranc
374  USE vtk
375 
376  IMPLICIT NONE
377  LOGICAL :: VTU_FRAME_FILE_EXISTS
378  INTEGER :: ISTAT,BUFF1,BUFF2,L
379  INTEGER :: MODE ! MODE = 0 : Write regular VTK region file
380  ! MODE = 1 : Write debug VTK region file (VTK_DBG_FILE = .TRUE.)
381 
382 #ifdef MPI
383  call mpi_barrier(mpi_comm_world,mpierr)
384 #endif
385 
386 ! Only open the file from head node when not using distributed I/O
387  IF (mype /= pe_io.AND.(.NOT.bdist_io)) RETURN
388 
389  IF(trim(vtu_dir)/='.') CALL create_dir(trim(vtu_dir))
390 
391  IF(time_dependent_filename) THEN
392  INQUIRE(file=vtu_frame_filename,exist=vtu_frame_file_exists)
393  IF(vtu_frame_file_exists) THEN
394  OPEN(convert='BIG_ENDIAN',unit = vtu_frame_unit, file = trim(vtu_frame_filename))
395  DO l = 1, dimension_vtk
396  IF(vtk_defined(l)) THEN
397  READ(vtu_frame_unit,*)buff1,buff2
398  frame(l)=buff2
399  ENDIF
400  ENDDO
401  CLOSE(vtu_frame_unit)
402  ENDIF
403  IF(reset_frame_at_time_zero.AND.time==zero) THEN
404  DO l = 1, dimension_vtk
405  IF(l==vtk_region) frame(l)=-1
406  ENDDO
407  ENDIF
408  DO l = 1, dimension_vtk
409  IF(l==vtk_region) frame(l) = frame(l) + 1
410  ENDDO
411  ENDIF
412 
413  IF (bdist_io.AND.number_of_vtk_cells>0) THEN
414 
415 
416 ! For distributed I/O, define the file name for each processor
417  IF(time_dependent_filename.AND.mode==0) THEN
419  ELSE
420  WRITE(vtu_filename,25) trim(vtk_filebase(vtk_region)),mype
421  ENDIF
422  ELSE
423  IF(mype.EQ.pe_io) THEN
424  IF(time_dependent_filename.AND.mode==0) THEN
426  ELSE
427  WRITE(vtu_filename,35) trim(vtk_filebase(vtk_region))
428  ENDIF
429  END IF
430  END IF
431 
432 ! Add the VTU directory path if necessary
433 
434  IF(trim(vtu_dir)/='.') vtu_filename='./'//trim(vtu_dir)//'/'//vtu_filename
435 
436 ! Echo
437  IF (full_log) THEN
438  IF (.NOT.bdist_io) THEN
439  WRITE(*,10,advance='NO')' WRITING VTU FILE : ', trim(vtu_filename),' .'
440  ELSE
441  IF(mype==pe_io) WRITE(*,15,advance='NO')' EACH PROCESOR IS WRITING ITS OWN VTU FILE.'
442  ENDIF
443  ENDIF
444 
445 ! Open File
446 ! OPEN(CONVERT='BIG_ENDIAN',UNIT = VTU_UNIT, FILE = TRIM(VTU_FILENAME),FORM='BINARY',IOSTAT=ISTAT)
447 
448 
449  IF(number_of_vtk_cells>0) THEN
450 
451  OPEN(convert='BIG_ENDIAN',unit = vtu_unit, &
452  file = trim(vtu_filename), &
453  form = 'UNFORMATTED', & ! works with gfortran 4.3.4 and ifort 10.1 but may not be supported by all compilers
454  ! use 'BINARY' if 'UNFORMATTED' is not supported
455  access = 'STREAM', & ! works with gfortran 4.3.4 and ifort 10.1 but may not be supported by all compilers
456  ! use 'SEQUENTIAL' if 'STREAM' is not supported
457  action = 'WRITE', &
458  iostat=istat)
459 
460 
461  IF (istat /= 0) THEN
462  IF(dmp_log) WRITE(unit_log, 1001) vtu_filename, vtu_unit,vtu_dir
463  IF(full_log.AND.mype == pe_io) WRITE(*, 1001) vtu_filename, vtu_unit,vtu_dir
464  CALL mfix_exit(mype)
465  ENDIF
466 
467 
468  1001 FORMAT(/1x,70('*')//, ' From: OPEN_VTU_FILE',/,' Message: ', &
469  'Error opening vtu file. Terminating run.',/10x, &
470  'File name: ',a,/10x, &
471  'DES_UNIT : ',i4, /10x, &
472  'PLEASE VERIFY THAT VTU_DIR EXISTS: ', a, &
473  /1x,70('*')/)
474 
475 
476  ! Write file Header
477  buffer='<?xml version="1.0"?>'
478  WRITE(vtu_unit)trim(buffer)//end_rec
479 
480 
481  WRITE(buffer,110)'<!-- Time =',time,' sec. -->'
482  WRITE(vtu_unit)trim(buffer)//end_rec
483 
484  buffer='<VTKFile type="UnstructuredGrid" version="0.1" byte_order="BigEndian">'
485  WRITE(vtu_unit)trim(buffer)//end_rec
486 
487  buffer=' <UnstructuredGrid>'
488  WRITE(vtu_unit)trim(buffer)//end_rec
489 
490  ENDIF
491 ! For distributed I/O, open .pvtu file that combines all *.vtu files for a given FRAME
492 ! this is a simple ASCII file
493 
494  IF (mype == pe_io.AND.bdist_io) THEN
495 
496  IF(time_dependent_filename.AND.mode==0) THEN
497  WRITE(pvtu_filename,40) trim(vtk_filebase(vtk_region)),frame(vtk_region)
498  ELSE
499  WRITE(pvtu_filename,45) trim(vtk_filebase(vtk_region))
500  ENDIF
501 
502  IF(trim(vtu_dir)/='.') pvtu_filename='./'//trim(vtu_dir)//'/'//pvtu_filename
503 
504  OPEN(convert='BIG_ENDIAN',unit = pvtu_unit, file = trim(pvtu_filename))
505 
506  WRITE(pvtu_unit,100) '<?xml version="1.0"?>'
507  WRITE(pvtu_unit,110) '<!-- Time =',time,' sec. -->'
508  WRITE(pvtu_unit,120) '<VTKFile type="PUnstructuredGrid"',&
509  ' version="0.1" byte_order="BigEndian">'
510 
511  WRITE(pvtu_unit,100) ' <PUnstructuredGrid GhostLevel="0">'
512  WRITE(pvtu_unit,100) ' <PPoints>'
513  WRITE(pvtu_unit,100) ' <PDataArray type="Float32" Name="coordinates" NumberOfComponents="3" &
514  &format="appended" offset=" 0" />'
515  WRITE(pvtu_unit,100) ' </PPoints>'
516  WRITE(pvtu_unit,100) ''
517  WRITE(pvtu_unit,100) ' <PCellData Scalars="scalars">'
518 
519  ENDIF
520 
521 100 FORMAT(a)
522 110 FORMAT(a,e14.7,a)
523 120 FORMAT(a,a)
524 10 FORMAT(/1x,3a)
525 15 FORMAT(/1x,a)
526 20 FORMAT(a,"_",i4.4,"_",i5.5,".vtu")
527 25 FORMAT(a,"_",i5.5,".vtu")
528 30 FORMAT(a,"_",i4.4,".vtu")
529 35 FORMAT(a,".vtu")
530 40 FORMAT(a,"_",i4.4,".pvtu")
531 45 FORMAT(a,".pvtu")
532 
533  RETURN
534 
535  END SUBROUTINE open_vtu_file_bin
536 
537 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
538 ! C
539 ! Module name: WRITE_GEOMETRY_IN_VTU_BIN C
540 ! Purpose: Write Geometry and connectivity in a vtu file C
541 ! Binary format C
542 ! C
543 ! Author: Jeff Dietiker Date: 21-Feb-08 C
544 ! Reviewer: Date: C
545 ! C
546 ! Revision Number # Date: ##-###-## C
547 ! Author: # C
548 ! Purpose: # C
549 ! C
550 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
551  SUBROUTINE write_geometry_in_vtu_bin(PASS)
553  use, INTRINSIC :: iso_c_binding
554  USE param
555  USE param1
556  USE parallel
557  USE constant
558  USE run
559  USE toleranc
560  USE geometry
561  USE indices
562  USE compar
563  USE mpi_utility
564  USE sendrecv
565  USE quadric
566  USE cutcell
567  USE fldvar
568  USE vtk
569  USE cdist
570  USE functions
571 
572  IMPLICIT NONE
573 
574  INTEGER :: IJK,L
575  INTEGER :: OFFSET
576 
577  INTEGER :: CELL_TYPE
578 
579  REAL(c_float) :: float
580  INTEGER(c_int) :: int
581 
582  INTEGER :: nbytes_xyz,nbytes_connectivity,nbytes_offset,nbytes_type
583  INTEGER :: offset_xyz,offset_connectivity,offset_offset,offset_type
584 
585  INTEGER :: PASS
586  INTEGER :: WRITE_HEADER = 1
587  INTEGER :: WRITE_DATA = 2
588 
589 
590 ! First a series of tags is written for the geometry (PASS=WRITE_HEADER)
591 ! - Coordinates
592 ! - Connectivity
593 ! - Connectivity offset
594 ! - cell types
595 !
596 
597 ! Since the data is appended (i.e., written after all tags), the
598 ! offset, in number of bytes must be specified. The offset includes
599 ! the size of the data for each field, plus the size of the integer
600 ! that stores the number of bytes. this is why the offset of a field
601 ! equals the offset of the previous field plus sizeof(int) plus the
602 ! number of bytes of the field.
603 
604 ! Next, the actual data is written for the geometry (PASS=WRITE_DATA)
605 ! The DATA is converted to single precision to save memory.
606 
607  IF (mype == pe_io.AND.(.NOT.bdist_io)) THEN
608 ! The number of points and number of VTK cells is now computed in
609 ! SETUP_VTK_REGION
610 
611 ! Number of bytes of each field
612  nbytes_xyz = number_of_points * 3 * sizeof(float)
613 
614  nbytes_connectivity = 0
615  DO ijk = 1,ijkmax3
616  IF (belongs_to_vtk_subdomain(ijk)) THEN
617  nbytes_connectivity = nbytes_connectivity + global_number_of_nodes(ijk)
618  ENDIF
619  END DO
620  nbytes_connectivity = nbytes_connectivity * sizeof(int)
621 
622 
623  nbytes_offset = number_of_vtk_cells * sizeof(int)
624 
625  nbytes_type = number_of_vtk_cells * sizeof(int)
626 
627 
628 ! Offset of each field
629  offset_xyz = 0
630  offset_connectivity = offset_xyz + sizeof(int) + nbytes_xyz
631  offset_offset = offset_connectivity + sizeof(int) + nbytes_connectivity
632  offset_type = offset_offset + sizeof(int) + nbytes_offset
633 
634 
635  IF(pass==write_header) THEN
636 
637  WRITE(buffer,100)' <Piece NumberOfPoints="',number_of_points,'" NumberOfCells="',number_of_vtk_cells,'" >'
638  WRITE(vtu_unit)trim(buffer)//end_rec
639 
640  WRITE(buffer,110)' <Points>'
641  WRITE(vtu_unit)trim(buffer)//end_rec
642 
643  WRITE(buffer,100)' <DataArray type="Float32" Name="coordinates" NumberOfComponents="3" &
644  &format="appended" offset="',offset_xyz,'" />'
645  WRITE(vtu_unit)trim(buffer)//end_rec
646 
647  WRITE(buffer,110)' </Points>'
648  WRITE(vtu_unit)trim(buffer)//end_rec
649 
650  WRITE(buffer,110)' <Cells>'
651  WRITE(vtu_unit)trim(buffer)//end_rec
652 
653  WRITE(buffer,100)' <DataArray type="Int32" Name="connectivity" format="appended" offset="', &
654  offset_connectivity,'" />'
655  WRITE(vtu_unit)trim(buffer)//end_rec
656 
657  WRITE(buffer,100)' <DataArray type="Int32" Name="offsets" format="appended" offset="',offset_offset,'" />'
658  WRITE(vtu_unit)trim(buffer)//end_rec
659 
660  WRITE(buffer,100)' <DataArray type="Int32" Name="types" format="appended" offset="',offset_type,'" />'
661  WRITE(vtu_unit)trim(buffer)//end_rec
662 
663  WRITE(buffer,110)' </Cells>'
664  WRITE(vtu_unit)trim(buffer)//end_rec
665 
666  vtu_offset = offset_type + sizeof(int) + nbytes_type ! Store offset for first variable to be written
667 
668  WRITE(buffer,110)' <CellData>' ! Preparing CellData
669  WRITE(vtu_unit)trim(buffer)//end_rec
670 
671 
672 
673 
674  ELSEIF(pass==write_data) THEN
675 
676  WRITE(buffer,110)' </CellData>'
677  WRITE(vtu_unit)trim(buffer)//end_rec
678 
679  WRITE(buffer,110)' </Piece>'
680  WRITE(vtu_unit)trim(buffer)//end_rec
681 
682 
683  WRITE(buffer,110)' </UnstructuredGrid>'
684  WRITE(vtu_unit)trim(buffer)//end_rec
685 
686 
687  WRITE(buffer,110)' <AppendedData encoding="raw">'
688  WRITE(vtu_unit)trim(buffer)//end_rec
689 
690 
691 ! Starting raw binary data with an underscore
692 
693  WRITE(buffer,110)'_'
694  WRITE(vtu_unit)trim(buffer)
695 
696 ! X,Y,Z coordinates
697  WRITE(vtu_unit) nbytes_xyz, (global_coords_of_points(1:3,l), l = 1,number_of_points)
698 
699 ! Connectivity
700  WRITE(vtu_unit) nbytes_connectivity
701 
702  DO ijk = 1,ijkmax3
703  IF (belongs_to_vtk_subdomain(ijk)) THEN
705  ENDIF
706  END DO
707 
708 ! Offsets
709  WRITE(vtu_unit) nbytes_offset
710 
711  offset = 0
712  DO ijk = 1,ijkmax3
713  IF (belongs_to_vtk_subdomain(ijk)) THEN
714  offset = offset + global_number_of_nodes(ijk)
715  WRITE(vtu_unit) offset
716  ENDIF
717  END DO
718 
719 ! Types
720  WRITE(vtu_unit)nbytes_type
721 
722  IF(no_k) THEN
723  cell_type = 7
724  ELSE
725  cell_type = 41
726  ENDIF
727 
728  DO ijk = 1,ijkmax3
729  IF (belongs_to_vtk_subdomain(ijk)) WRITE(vtu_unit) cell_type
730  END DO
731 
732 
733  ENDIF
734 
735 
736  ELSEIF(bdist_io) THEN
737 
738 ! For distributed I/O, it works the same as above, except, the data is local to each processor
739 ! First compute local number of cells and points
740 
741 ! The number of points and number of VTK cells is now computed in
742 ! SETUP_VTK_REGION
743 
744 ! Number of bytes of each field
745  nbytes_xyz = number_of_points * 3 * sizeof(float)
746 
747  nbytes_connectivity = 0
748  DO ijk = 1,ijkend3
749  IF (belongs_to_vtk_subdomain(ijk)) THEN
750  nbytes_connectivity = nbytes_connectivity + number_of_nodes(ijk)
751  ENDIF
752  END DO
753  nbytes_connectivity = nbytes_connectivity * sizeof(int)
754 
755 
756  nbytes_offset = number_of_vtk_cells * sizeof(int)
757 
758  nbytes_type = number_of_vtk_cells * sizeof(int)
759 
760 
761 ! Offset of each field
762  offset_xyz = 0
763  offset_connectivity = offset_xyz + sizeof(int) + nbytes_xyz
764  offset_offset = offset_connectivity + sizeof(int) + nbytes_connectivity
765  offset_type = offset_offset + sizeof(int) + nbytes_offset
766 
767 
768  IF(pass==write_header) THEN
769 
770  WRITE(buffer,100)' <Piece NumberOfPoints="',number_of_points,'" NumberOfCells="',number_of_vtk_cells,'" >'
771  WRITE(vtu_unit)trim(buffer)//end_rec
772 
773  WRITE(buffer,110)' <Points>'
774  WRITE(vtu_unit)trim(buffer)//end_rec
775 
776  WRITE(buffer,100)' <DataArray type="Float32" Name="coordinates" NumberOfComponents="3" &
777  &format="appended" offset="',offset_xyz,'" />'
778  WRITE(vtu_unit)trim(buffer)//end_rec
779 
780  WRITE(buffer,110)' </Points>'
781  WRITE(vtu_unit)trim(buffer)//end_rec
782 
783  WRITE(buffer,110)' <Cells>'
784  WRITE(vtu_unit)trim(buffer)//end_rec
785 
786  WRITE(buffer,100)' <DataArray type="Int32" Name="connectivity" format="appended" offset="', &
787  offset_connectivity,'" />'
788  WRITE(vtu_unit)trim(buffer)//end_rec
789 
790  WRITE(buffer,100)' <DataArray type="Int32" Name="offsets" format="appended" offset="',offset_offset,'" />'
791  WRITE(vtu_unit)trim(buffer)//end_rec
792 
793  WRITE(buffer,100)' <DataArray type="Int32" Name="types" format="appended" offset="',offset_type,'" />'
794  WRITE(vtu_unit)trim(buffer)//end_rec
795 
796  WRITE(buffer,110)' </Cells>'
797  WRITE(vtu_unit)trim(buffer)//end_rec
798 
799  vtu_offset = offset_type + sizeof(int) + nbytes_type ! Store offset for first variable to be written
800 
801  WRITE(buffer,110)' <CellData>' ! Preparing CellData
802  WRITE(vtu_unit)trim(buffer)//end_rec
803 
804 
805 
806 
807  ELSEIF(pass==write_data) THEN
808 
809  WRITE(buffer,110)' </CellData>'
810  WRITE(vtu_unit)trim(buffer)//end_rec
811 
812  WRITE(buffer,110)' </Piece>'
813  WRITE(vtu_unit)trim(buffer)//end_rec
814 
815 
816  WRITE(buffer,110)' </UnstructuredGrid>'
817  WRITE(vtu_unit)trim(buffer)//end_rec
818 
819 
820  WRITE(buffer,110)' <AppendedData encoding="raw">'
821  WRITE(vtu_unit)trim(buffer)//end_rec
822 
823 
824 ! Starting raw binary data with an underscore
825 
826  WRITE(buffer,110)'_'
827  WRITE(vtu_unit)trim(buffer)
828 
829 ! X,Y,Z coordinates
830  WRITE(vtu_unit) nbytes_xyz, (coords_of_points(l,1:3), l = 1,number_of_points)
831 
832 ! Connectivity
833  WRITE(vtu_unit) nbytes_connectivity
834 
835  DO ijk = 1,ijkend3
836  IF (belongs_to_vtk_subdomain(ijk)) THEN
837  WRITE(vtu_unit) (cleaned_connectivity(ijk,l)-1,l=1,number_of_nodes(ijk))
838  ENDIF
839  END DO
840 
841 ! Offsets
842  WRITE(vtu_unit) nbytes_offset
843 
844  offset = 0
845  DO ijk = 1,ijkend3
846  IF (belongs_to_vtk_subdomain(ijk)) THEN
847  offset = offset + number_of_nodes(ijk)
848  WRITE(vtu_unit) offset
849  ENDIF
850  END DO
851 
852 ! Types
853  WRITE(vtu_unit)nbytes_type
854 
855  IF(no_k) THEN
856  cell_type = 7
857  ELSE
858  cell_type = 41
859  ENDIF
860 
861  DO ijk = 1,ijkend3
862  IF (belongs_to_vtk_subdomain(ijk)) WRITE(vtu_unit) cell_type
863  END DO
864 
865 
866  ENDIF
867 
868 
869  ENDIF
870 
871 
872 100 FORMAT(a,i12,a,i12,a)
873 110 FORMAT(a)
874 
875  RETURN
876 
877  END SUBROUTINE write_geometry_in_vtu_bin
878 
879 
880 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
881 ! C
882 ! Module name: WRITE_SCALAR_IN_VTU_BIN C
883 ! Purpose: Write Scalar variable in a vtu file C
884 ! Binary format C
885 ! C
886 ! Author: Jeff Dietiker Date: 21-Feb-08 C
887 ! Reviewer: Date: C
888 ! C
889 ! Revision Number # Date: ##-###-## C
890 ! Author: # C
891 ! Purpose: # C
892 ! C
893 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
894  SUBROUTINE write_scalar_in_vtu_bin(VAR_NAME,VAR,PASS)
896  use, INTRINSIC :: iso_c_binding
897  USE param
898  USE param1
899  USE parallel
900  USE constant
901  USE run
902  USE toleranc
903  USE geometry
904  USE indices
905  USE compar
906  USE mpi_utility
907  USE sendrecv
908  USE quadric
909  USE cutcell
910  USE fldvar
911  USE vtk
912  USE cdist
913  USE output
914  USE functions
915 
916  IMPLICIT NONE
917  INTEGER :: I,IJK
918 
919  CHARACTER (*) :: VAR_NAME
920  DOUBLE PRECISION, DIMENSION(DIMENSION_3) :: VAR
921  DOUBLE PRECISION, ALLOCATABLE :: GLOBAL_VAR(:)
922  DOUBLE PRECISION, DIMENSION(DIMENSION_3) :: TMP_VAR
923 
924  REAL(c_float) :: float
925 
926  INTEGER :: nbytes_scalar
927 
928  INTEGER :: PASS
929  INTEGER :: WRITE_HEADER = 1
930  INTEGER :: WRITE_DATA = 2
931 
932  IF (.NOT.bdist_io) THEN
933 
934 ! For each scalar, write a tag, with corresponding offset
935 
936  nbytes_scalar = number_of_vtk_cells * sizeof(float)
937 
938  IF(pass==write_header) THEN
939 ! For each scalar, write a tag, with corresponding offset
940 
941  DO i = 1,len_trim(var_name)
942  IF(var_name(i:i) == ' ') var_name(i:i) = '_'
943  ENDDO
944 
945  WRITE(buffer,90)' <DataArray type="Float32" Name="', &
946  trim(var_name),'" format="appended" offset="',vtu_offset,'" />'
947  WRITE(vtu_unit)trim(buffer)//end_rec
948 
949  vtu_offset = vtu_offset + sizeof(float) + nbytes_scalar
950 
951 
952  ELSEIF(pass==write_data) THEN
953 ! and write the data, always preceded by its size in number of bytes
954 
955  IF (mype == pe_io) THEN
956  allocate (global_var(ijkmax3))
957  ELSE
958  allocate (global_var(1))
959  ENDIF
960 
961  IF(re_indexing) THEN
962  CALL unshift_dp_array(var,tmp_var)
963  CALL gather (tmp_var,global_var,root)
964  ELSE
965  CALL gather (var,global_var,root)
966  ENDIF
967 
968  IF (mype /= pe_io) RETURN
969 
970 
971  WRITE(vtu_unit) nbytes_scalar
972 
973  DO ijk = 1,ijkmax3
974  IF (belongs_to_vtk_subdomain(ijk)) WRITE(vtu_unit) REAL(global_var(ijk))
975  ENDDO
976 
977 
978  Deallocate (global_var)
979 
980  ENDIF
981 
982 
983  ELSE ! BDIST_IO=.TRUE.
984 
985 
986  nbytes_scalar = number_of_vtk_cells * sizeof(float)
987 
988  IF(pass==write_header) THEN
989 ! For each scalar, write a tag, with corresponding offset
990 
991  DO i = 1,len_trim(var_name)
992  IF(var_name(i:i) == ' ') var_name(i:i) = '_'
993  ENDDO
994 
995  WRITE(buffer,90)' <DataArray type="Float32" Name="', &
996  trim(var_name),'" format="appended" offset="',vtu_offset,'" />'
997  WRITE(vtu_unit)trim(buffer)//end_rec
998 
999 
1000  vtu_offset = vtu_offset + sizeof(float) + nbytes_scalar
1001 
1002 
1003  ELSEIF(pass==write_data) THEN
1004 ! and write the data, always preceded by its size in number of bytes
1005 
1006  WRITE(vtu_unit) nbytes_scalar
1007 
1008  DO ijk = 1,ijkend3
1009  IF (belongs_to_vtk_subdomain(ijk)) WRITE(vtu_unit) REAL(var(ijk))
1010  ENDDO
1011 
1012  ENDIF
1013 
1014 
1015  IF (mype == pe_io) THEN ! Update pvtu file with variable name
1016  WRITE(pvtu_unit,90) ' <DataArray type="Float32" Name="', &
1017  trim(var_name),'" format="appended" offset="',vtu_offset,'" />'
1018  ENDIF
1019 
1020 
1021  ENDIF
1022 
1023 
1024  IF (pass==write_data.AND.full_log.AND.mype == pe_io) WRITE(*,10,advance='NO')'.'
1025 
1026 10 FORMAT(a)
1027 90 FORMAT(a,a,a,i12,a)
1028 
1029  RETURN
1030 
1031  END SUBROUTINE write_scalar_in_vtu_bin
1032 
1033 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1034 ! C
1035 ! Module name: WRITE_VECTOR_IN_VTU C
1036 ! Purpose: Write Vector variable in a vtu file C
1037 ! Binary format C
1038 ! C
1039 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1040 ! Reviewer: Date: C
1041 ! C
1042 ! Revision Number # Date: ##-###-## C
1043 ! Author: # C
1044 ! Purpose: # C
1045 ! C
1046 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1047  SUBROUTINE write_vector_in_vtu_bin(VAR_NAME,VARX,VARY,VARZ,PASS)
1049  use, INTRINSIC :: iso_c_binding
1050  USE param
1051  USE param1
1052  USE parallel
1053  USE constant
1054  USE run
1055  USE toleranc
1056  USE geometry
1057  USE indices
1058  USE compar
1059  USE mpi_utility
1060  USE sendrecv
1061  USE quadric
1062  USE cutcell
1063  USE fldvar
1064  USE vtk
1065  USE cdist
1066  USE output
1067  USE functions
1068 
1069  IMPLICIT NONE
1070  INTEGER :: IJK
1071 
1072  CHARACTER (*) :: VAR_NAME
1073  DOUBLE PRECISION, DIMENSION(DIMENSION_3) :: VARX,VARY,VARZ
1074  DOUBLE PRECISION, ALLOCATABLE :: GLOBAL_VARX(:),GLOBAL_VARY(:),GLOBAL_VARZ(:)
1075  DOUBLE PRECISION, DIMENSION(DIMENSION_3) :: TMP_VAR
1076 
1077  REAL(c_float) :: float
1078 
1079  INTEGER :: nbytes_vector
1080 
1081  INTEGER :: PASS
1082  INTEGER :: WRITE_HEADER = 1
1083  INTEGER :: WRITE_DATA = 2
1084 
1085  IF (.NOT.bdist_io) THEN
1086 
1087  nbytes_vector = number_of_vtk_cells * 3 * sizeof(float)
1088 
1089  IF(pass==write_header) THEN
1090 ! For each vector, write a tag, with corresponding offset
1091 
1092  WRITE(buffer,90)' <DataArray type="Float32" Name="', &
1093  trim(var_name),'" NumberOfComponents="3" format="appended" offset="',vtu_offset,'" />'
1094  WRITE(vtu_unit)trim(buffer)//end_rec
1095 
1096 
1097  vtu_offset = vtu_offset + sizeof(float) + nbytes_vector
1098 
1099 
1100  ELSEIF(pass==write_data) THEN
1101 ! and write the data, always preceded by its size in number of bytes
1102 
1103  IF (mype == pe_io) THEN
1104  allocate (global_varx(ijkmax3))
1105  allocate (global_vary(ijkmax3))
1106  allocate (global_varz(ijkmax3))
1107  ELSE
1108  allocate (global_varx(1))
1109  allocate (global_vary(1))
1110  allocate (global_varz(1))
1111  ENDIF
1112 
1113  IF(re_indexing) THEN
1114  CALL unshift_dp_array(varx,tmp_var)
1115  call gather (tmp_var,global_varx,root)
1116 
1117  CALL unshift_dp_array(vary,tmp_var)
1118  call gather (tmp_var,global_vary,root)
1119 
1120  CALL unshift_dp_array(varz,tmp_var)
1121  call gather (tmp_var,global_varz,root)
1122 
1123  ELSE
1124  call gather (varx,global_varx,root)
1125  call gather (vary,global_vary,root)
1126  call gather (varz,global_varz,root)
1127  ENDIF
1128 
1129 
1130  IF (mype /= pe_io) RETURN
1131 
1132 
1133  WRITE(vtu_unit) nbytes_vector
1134 
1135  DO ijk = 1,ijkmax3
1136  IF (belongs_to_vtk_subdomain(ijk)) THEN
1137  WRITE(vtu_unit) REAL(GLOBAL_VARX(IJK)),REAL(GLOBAL_VARY(IJK)),REAL(global_varz(ijk))
1138  ENDIF
1139  ENDDO
1140 
1141 
1142  Deallocate (global_varx)
1143  Deallocate (global_vary)
1144  Deallocate (global_varz)
1145 
1146  ENDIF
1147 
1148 
1149  ELSE ! BDIST_IO=.TRUE.
1150 
1151 
1152  nbytes_vector = number_of_vtk_cells * 3 * sizeof(float)
1153 
1154  IF(pass==write_header) THEN
1155 ! For each vector, write a tag, with corresponding offset
1156 
1157 
1158  WRITE(buffer,90)' <DataArray type="Float32" Name="', &
1159  trim(var_name),'" NumberOfComponents="3" format="appended" offset="',vtu_offset,'" />'
1160  WRITE(vtu_unit)trim(buffer)//end_rec
1161 
1162 
1163  vtu_offset = vtu_offset + sizeof(float) + nbytes_vector
1164 
1165 
1166  ELSEIF(pass==write_data) THEN
1167 ! and write the data, always preceded by its size in number of bytes
1168 
1169  WRITE(vtu_unit) nbytes_vector
1170 
1171  DO ijk = 1,ijkend3
1172  IF (belongs_to_vtk_subdomain(ijk)) THEN
1173  WRITE(vtu_unit) REAL(VARX(IJK)),REAL(VARY(IJK)),REAL(varz(ijk))
1174  ENDIF
1175  ENDDO
1176 
1177  ENDIF
1178 
1179 
1180  IF (mype == pe_io) THEN ! Update pvtu file with variable name
1181  WRITE(pvtu_unit,90)' <DataArray type="Float32" Name="', &
1182  trim(var_name),'" NumberOfComponents="3" format="appended" offset="',vtu_offset,'" />'
1183  ENDIF
1184 
1185  ENDIF
1186 
1187 
1188  IF (pass==write_data.AND.full_log.AND.mype == pe_io) WRITE(*,10,advance='NO')'.'
1189 
1190 10 FORMAT(a)
1191 90 FORMAT(a,a,a,i12,a)
1192 
1193  RETURN
1194 
1195  END SUBROUTINE write_vector_in_vtu_bin
1196 
1197 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1198 ! C
1199 ! Module name: CLOSE_VTU_FILE_BIN C
1200 ! Purpose: Close a vtu file C
1201 ! Binary format C
1202 ! C
1203 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1204 ! Reviewer: Date: C
1205 ! C
1206 ! Revision Number # Date: ##-###-## C
1207 ! Author: # C
1208 ! Purpose: # C
1209 ! C
1210 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1211  SUBROUTINE close_vtu_file_bin(MODE)
1213  USE compar
1214  Use run
1215  USE vtk
1216  use cdist
1217  USE mpi_utility
1218 
1219  IMPLICIT NONE
1220 
1221  INTEGER:: N
1222  CHARACTER (LEN=32) :: VTU_NAME
1223  INTEGER, DIMENSION(0:numPEs-1) :: ALL_VTK_CELL_COUNT
1224  INTEGER :: IERR
1225  INTEGER :: MODE ! MODE = 0 : Write regular VTK region file
1226  ! MODE = 1 : Write debug VTK region file (VTK_DBG_FILE = .TRUE.)
1227 
1228 
1229  IF (mype /= pe_io.AND.(.NOT.bdist_io)) RETURN
1230 
1231  IF(number_of_vtk_cells>0) THEN
1232 
1233 ! Write last tags and close the vtu file
1234  WRITE(buffer,110)' </AppendedData>'
1235  WRITE(vtu_unit)end_rec//trim(buffer)//end_rec
1236 
1237  WRITE(buffer,110)'</VTKFile>'
1238  WRITE(vtu_unit)trim(buffer)//end_rec
1239 
1240  CLOSE(vtu_unit)
1241 
1242  ENDIF
1243 
1244 ! Update pvtu file and close
1245 
1246  IF(bdist_io) CALL allgather_1i (number_of_vtk_cells,all_vtk_cell_count,ierr)
1247 
1248  IF (mype == pe_io.AND.bdist_io) THEN
1249  WRITE(pvtu_unit,100) ' </PCellData>'
1250 
1251  DO n = 0,numpes-1
1252  IF(all_vtk_cell_count(n)>0) THEN
1253  IF(time_dependent_filename.AND.mode==0) THEN
1254  WRITE(vtu_name,20) trim(vtk_filebase(vtk_region)),frame(vtk_region),n
1255  ELSE
1256  WRITE(vtu_name,25) trim(vtk_filebase(vtk_region)),n
1257  ENDIF
1258 
1259  WRITE(pvtu_unit,110) ' <Piece Source="',trim(vtu_name),'"/>'
1260  ENDIF
1261  ENDDO
1262 
1263 
1264  WRITE(pvtu_unit,100) ' </PUnstructuredGrid>'
1265  WRITE(pvtu_unit,100) '</VTKFile>'
1266  CLOSE(pvtu_unit)
1267  ENDIF
1268 
1269 
1270 20 FORMAT(a,"_",i4.4,"_",i5.5,".vtu")
1271 25 FORMAT(a,"_",i5.5,".vtu")
1272 
1273 100 FORMAT(a)
1274 110 FORMAT(a,a,a)
1275 
1276  RETURN
1277 
1278  END SUBROUTINE close_vtu_file_bin
1279 
1280 
1281 
1282 
1283 
1284 
1285 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1286 ! C
1287 ! Module name: OPEN_PVD_FILE C
1288 ! Purpose: Open a PVD file and writes the header C
1289 ! C
1290 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1291 ! Reviewer: Date: C
1292 ! C
1293 ! Revision Number # Date: ##-###-## C
1294 ! Author: # C
1295 ! Purpose: # C
1296 ! C
1297 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1298  SUBROUTINE open_pvd_file
1300  USE compar
1301  USE constant
1302  USE cutcell
1303  USE exit, only: mfix_exit
1304  USE fldvar
1305  USE functions
1306  USE geometry
1307  USE indices
1308  USE output
1309  USE parallel
1310  USE param
1311  USE param1
1312  USE quadric
1313  USE run
1314  USE sendrecv
1315  USE toleranc
1316  USE vtk
1317 
1318  IMPLICIT NONE
1319  LOGICAL :: PVD_EXISTS
1320 
1321  IF (mype /= pe_io) RETURN
1322 
1323  pvd_filename = trim(vtk_filebase(vtk_region)) // '.pvd'
1324 
1325 ! First, check if the file already exists.
1326 
1327  INQUIRE(file=pvd_filename,exist=pvd_exists)
1328 
1329 ! The first time this subroutine is executed, properly initialize the pvd file
1330 
1331  IF(.NOT.pvd_file_initialized(vtk_region)) THEN
1332 
1333  IF(run_type == 'NEW'.OR.run_type=='RESTART_2')THEN
1334  ! For a new or RESTART_2 run, the pvd file should not exist, and is created with appropriate header
1335  IF (.NOT.pvd_exists) THEN
1336  OPEN(convert='BIG_ENDIAN',unit = pvd_unit, file = trim(pvd_filename))
1337  WRITE(pvd_unit,100) '<?xml version="1.0"?>'
1338  WRITE(pvd_unit,100) '<VTKFile type="Collection" version="0.1" byte_order="LittleEndian">'
1339  WRITE(pvd_unit,100) '<Collection>'
1340 ! CALL UPDATE_AND_CLOSE_PVD_FILE
1342  ELSE ! If the pvd file exists, print error message and exits
1343  WRITE(*,1002) trim(pvd_filename)
1344  WRITE(unit_log, 1002) trim(pvd_filename)
1345  CALL mfix_exit(mype)
1346  ENDIF
1347  ELSE
1348  ! For a restart_1 run, the pvd file must exist
1349  IF (.NOT.pvd_exists) THEN
1350  ! If the pvd file does not exist, print error message and exits
1351  WRITE(*,1003) trim(pvd_filename)
1352  WRITE(unit_log, 1003) trim(pvd_filename)
1353  CALL mfix_exit(mype)
1354  ELSE
1355  ! If it already exists, go to the bottom of the file and prepare to append data (remove last two lines)
1356  OPEN(convert='BIG_ENDIAN',unit=pvd_unit,file = trim(pvd_filename),position="APPEND",status='OLD')
1357  backspace(pvd_unit)
1358  backspace(pvd_unit)
1360  ENDIF
1361  ENDIF
1362  ELSE
1363  ! When properly initialized, open the file and go to the
1364  ! bottom of the file and prepare to append data (remove last two lines)
1365  OPEN(convert='BIG_ENDIAN',unit=pvd_unit,file = trim(pvd_filename),position="APPEND",status='OLD')
1366  backspace(pvd_unit)
1367  backspace(pvd_unit)
1368  ENDIF
1369 
1370 
1371 100 FORMAT(a)
1372 
1373 1002 FORMAT(/1x,70('*')/,' From: OPEN_PVD_FILE',/,' Message: ', &
1374  a,' already exists in the run directory.',/10x, &
1375  'This is not allowed for a new run.',/10x, &
1376  'Terminating run.',/1x,70('*')/)
1377 
1378 1003 FORMAT(/1x,70('*')/,' From: OPEN_PVD_FILE',/,' Message: ', &
1379  a,' is missing from the the run directory,',/10x, &
1380  ' and must be present for a restart run.',/10x, &
1381  'Terminating run.',/1x,70('*')/)
1382 
1383  RETURN
1384 
1385  END SUBROUTINE open_pvd_file
1386 
1387 
1388 
1389 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1390 ! C
1391 ! Module name: UPDATE_AND_CLOSE_PVD_FILE C
1392 ! Purpose: Updates and close a pvd file C
1393 ! C
1394 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1395 ! Reviewer: Date: C
1396 ! C
1397 ! Revision Number # Date: ##-###-## C
1398 ! Author: # C
1399 ! Purpose: # C
1400 ! C
1401 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1402  SUBROUTINE update_and_close_pvd_file
1404  USE compar
1405  use cdist
1406  USE run
1407  USE vtk
1408 
1409  IMPLICIT NONE
1410 
1411  CHARACTER (LEN=255) :: FILENAME
1412  CHARACTER (LEN=5) :: EXT
1413 
1414  IF (mype /= pe_io) RETURN
1415 
1416  IF(vtk_data(vtk_region)=='C') THEN
1417  ext = '.pvtu'
1418  ELSEIF(vtk_data(vtk_region)=='P') THEN
1419  ext = '.pvtp'
1420  ENDIF
1421 
1422  IF(.NOT.bdist_io) THEN
1423  filename=vtu_filename
1424  ELSE
1425  IF(time_dependent_filename) THEN
1426  WRITE(filename,40) trim(vtk_filebase(vtk_region)),frame(vtk_region),ext
1427  ELSE
1428  WRITE(filename,45) trim(vtk_filebase(vtk_region)),ext
1429  ENDIF
1430  IF(trim(vtu_dir)/='.') filename='./'//trim(vtu_dir)//'/'//filename
1431  ENDIF
1432 
1433 ! Write the data to the file
1434  WRITE(pvd_unit,100)&
1435  '<DataSet timestep="',time,'" ',& ! simulation time
1436  'group="" part="0" ',& ! necessary file data
1437  'file="',trim(filename),'"/>' ! file name of vtp
1438 
1439 ! Write the closing tags
1440  WRITE(pvd_unit,110)'</Collection>'
1441  WRITE(pvd_unit,110)'</VTKFile>'
1442 
1443  CLOSE(pvd_unit)
1444 
1445 ! 40 FORMAT(A,"_",I4.4,".pvtu")
1446 ! 45 FORMAT(A,".pvtu")
1447 40 FORMAT(a,"_",i4.4,a5)
1448 45 FORMAT(a,a5)
1449 100 FORMAT(6x,a,e14.7,5a)
1450 110 FORMAT(a)
1451 
1452  RETURN
1453 
1454  END SUBROUTINE update_and_close_pvd_file
1455 
1456 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1457 ! C
1458 ! Module name: WRITE_CUT_SURFACE_VTK C
1459 ! Purpose: Writes the cut cell surface in VTK format C
1460 ! C
1461 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1462 ! Reviewer: Date: C
1463 ! C
1464 ! Revision Number # Date: ##-###-## C
1465 ! Author: # C
1466 ! Purpose: # C
1467 ! C
1468 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1469  SUBROUTINE write_cut_surface_vtk
1471  USE compar
1472  USE constant
1473  USE cut_cell_preproc, only: eval_f
1474  USE cutcell
1475  USE exit, only: mfix_exit
1476  USE fldvar
1477  USE functions
1478  USE geometry
1479  USE indices
1480  USE parallel
1481  USE param
1482  USE param1
1483  USE polygon
1484  USE quadric
1485  USE run
1486  USE sendrecv
1487  USE stl
1488  USE toleranc
1489  USE vtk
1490 
1491  IMPLICIT NONE
1492 
1493  INTEGER :: L,IJK,NODE
1494  INTEGER :: POINT_ID,POLY_COUNT,FACE_ID,Q_ID
1495  INTEGER :: N_CUT_FACE_NODES
1496 
1497  INTEGER NUMBER_OF_FACES
1498  INTEGER NUMBER_OF_SURFACE_POINTS
1499 
1500  DOUBLE PRECISION, DIMENSION(3,15) :: COORD_CUT_FACE_NODES
1501  DOUBLE PRECISION, DIMENSION(3) :: NORMAL
1502 
1503  INTEGER, DIMENSION(DIMENSION_MAX_CUT_CELL,6) :: FACE_CONNECTIVITY
1504  INTEGER, DIMENSION(DIMENSION_MAX_CUT_CELL) :: NUMBER_OF_CUT_FACE_POINTS
1505 
1506  DOUBLE PRECISION, DIMENSION(DIMENSION_MAX_CUT_CELL) :: X_FACE_POINT
1507  DOUBLE PRECISION, DIMENSION(DIMENSION_MAX_CUT_CELL) :: Y_FACE_POINT
1508  DOUBLE PRECISION, DIMENSION(DIMENSION_MAX_CUT_CELL) :: Z_FACE_POINT
1509 
1510  DOUBLE PRECISION :: X_COPY,Y_COPY,Z_COPY,F_COPY
1511 
1512  LOGICAL :: CLIP_FLAG
1513 
1514  CHARACTER (LEN=255) :: FILENAME
1515 
1516  LOGICAL :: CORNER_POINT
1517  INTEGER :: NODE_OF_CORNER, IERROR
1518 
1519  IF(mype/=0) RETURN
1520 
1521 !======================================================================
1522 ! Set-up connectivity for each cell, i.e., regular cells and cut cells
1523 !======================================================================
1524 
1525  poly_count = 0
1526 
1527  number_of_surface_points = 0
1528 
1529  number_of_faces = 0
1530 
1531  DO ijk = 1,ijkmax3
1532 
1533  IF(global_cut_cell_at(ijk)) THEN
1534 
1535 !======================================================================
1536 ! Filter the connectivity to identify nodes belonging to cut face
1537 !======================================================================
1538 
1539 
1540  number_of_faces = number_of_faces + 1
1541 
1542  n_cut_face_nodes = 0
1543 
1544  CALL get_global_cell_node_coordinates(ijk,'SCALAR')
1545 
1546  DO l = 1, global_number_of_nodes(ijk)
1547  IF(global_connectivity(ijk,l)>ijkmax3) THEN ! One of the new point
1551  corner_point = .false.
1552  ELSE ! An existing point
1553  DO node = 1,8
1554  corner_point = .true.
1555  IF(global_connectivity(ijk,l) == ijk_of_node(node)) THEN
1556  node_of_corner = node
1557  x_copy = x_node(node)
1558  y_copy = y_node(node)
1559  z_copy = z_node(node)
1560 
1561  IF (global_snap(ijk_of_node(node))) THEN ! One of the snapped corner point which now belongs to the cut face
1562  n_cut_face_nodes = n_cut_face_nodes + 1
1563  coord_cut_face_nodes(1,n_cut_face_nodes) = x_copy
1564  coord_cut_face_nodes(2,n_cut_face_nodes) = y_copy
1565  coord_cut_face_nodes(3,n_cut_face_nodes) = z_copy
1566  ENDIF
1567  ENDIF
1568  END DO
1569 
1570  ENDIF
1571 
1572  IF(corner_point) THEN
1573  q_id = 1
1574 
1575  CALL eval_f('QUADRIC',x_copy,y_copy,z_copy,q_id,f_copy,clip_flag)
1576 
1577  CALL eval_f('POLYGON',x_copy,y_copy,z_copy,n_polygon,f_copy,clip_flag)
1578 
1579  CALL eval_f('USR_DEF',x_copy,y_copy,z_copy,n_usr_def,f_copy,clip_flag)
1580 
1581  IF(use_stl.OR.use_msh) f_copy = global_f_at(ijk_of_node(node_of_corner))
1582 
1583 ! CALL EVAL_STL_FCT_AT('SCALAR',IJK,NODE_OF_CORNER,F_COPY,CLIP_FLAG,BCID2)
1584  ELSE
1585  f_copy = zero
1586  ENDIF
1587 
1588  IF (abs(f_copy) < tol_f ) THEN ! belongs to cut face
1589  n_cut_face_nodes = n_cut_face_nodes + 1
1590  coord_cut_face_nodes(1,n_cut_face_nodes) = x_copy
1591  coord_cut_face_nodes(2,n_cut_face_nodes) = y_copy
1592  coord_cut_face_nodes(3,n_cut_face_nodes) = z_copy
1593  ENDIF
1594 
1595  END DO
1596 
1597  CALL reorder_polygon(n_cut_face_nodes,coord_cut_face_nodes,normal,ierror)
1598 
1599  number_of_cut_face_points(number_of_faces) = n_cut_face_nodes
1600  poly_count = poly_count + n_cut_face_nodes + 1
1601  DO node = 1,n_cut_face_nodes
1602  number_of_surface_points = number_of_surface_points + 1
1603 
1604  IF(number_of_surface_points>=dimension_max_cut_cell) THEN
1605  WRITE(*,3000) 'ERROR IN SUBROUTINE WRITE_3DCUT_SURFACE_VTK:'
1606  WRITE(*,3000) 'NUMBER_OF_SURFACE_POINTS>=DIMENSION_MAX_CUT_CELL:'
1607  WRITE(*,3000) 'INCREASE VALUE OF FAC_DIM_MAX_CUT_CELL.'
1608  WRITE(*,3010) 'CURRENT VALUE OF FAC_DIM_MAX_CUT_CELL =',fac_dim_max_cut_cell
1609  WRITE(*,3020) 'CURRENT VALUE OF DIMENSION_MAX_CUT_CELL =',dimension_max_cut_cell
1610  WRITE(*,3000) 'MFiX will exit now.'
1611  CALL mfix_exit(mype)
1612  ENDIF
1613 
1614  x_face_point(number_of_surface_points) = coord_cut_face_nodes(1,node)
1615  y_face_point(number_of_surface_points) = coord_cut_face_nodes(2,node)
1616  z_face_point(number_of_surface_points) = coord_cut_face_nodes(3,node)
1617  face_connectivity(number_of_faces,node) = number_of_surface_points
1618  ENDDO
1619 
1620  ENDIF
1621 
1622  END DO
1623 
1624  filename= trim(run_name) // '_boundary.vtk'
1625  filename = trim(filename)
1626  OPEN(convert='BIG_ENDIAN',unit = 123, file = filename)
1627  WRITE(123,1001)'# vtk DataFile Version 2.0'
1628  WRITE(123,1001)'3D CUT-CELL SURFACE'
1629  WRITE(123,1001)'ASCII'
1630 
1631  IF(no_k) THEN ! 2D GEOMETRY
1632  WRITE(123,1001)'DATASET UNSTRUCTURED_GRID'
1633  ELSE ! 3D GEOMETRY
1634  WRITE(123,1001)'DATASET POLYDATA'
1635  ENDIF
1636 
1637  WRITE(123,1010)'POINTS ',number_of_surface_points,' float'
1638 
1639  DO point_id = 1,number_of_surface_points
1640  WRITE(123,1020) x_face_point(point_id),y_face_point(point_id),z_face_point(point_id)
1641  ENDDO
1642 
1643  IF(no_k) THEN ! 2D GEOMETRY
1644 
1645  WRITE(123,1030)'CELLS ',number_of_faces,poly_count
1646  DO face_id = 1 , number_of_faces
1647  WRITE(123,1040) number_of_cut_face_points(face_id),(face_connectivity(face_id,l)-1,&
1648  l=1,number_of_cut_face_points(face_id))
1649  ENDDO
1650  WRITE(123,1030)'CELL_TYPES ',number_of_faces
1651  DO face_id = 1 , number_of_faces
1652  WRITE(123,1040) 3
1653  ENDDO
1654 
1655  ELSE ! 3D GEOMETRY
1656 
1657  WRITE(123,1030)'POLYGONS ',number_of_faces,poly_count
1658  DO face_id = 1 , number_of_faces
1659  WRITE(123,1040) number_of_cut_face_points(face_id),(face_connectivity(face_id,l)-1,&
1660  l=1,number_of_cut_face_points(face_id))
1661  ENDDO
1662 
1663  ENDIF
1664 
1665 1001 FORMAT(a)
1666 1010 FORMAT(a,i8,a)
1667 1020 FORMAT(3(e16.8,2x))
1668 1030 FORMAT(a,2(i8,2x))
1669 1040 FORMAT(20(i8,2x))
1670 3000 FORMAT(1x,a)
1671 3010 FORMAT(1x,a,f8.4)
1672 3020 FORMAT(1x,a,i8)
1673 3030 FORMAT(1x,a,a)
1674  CLOSE (123)
1675 
1676 
1677  WRITE(*,3030)'WROTE BOUNDARY IN VTK FILE : ',filename
1678  RETURN
1679 
1680 
1681  END SUBROUTINE write_cut_surface_vtk
1682 
1683 
1684 
1685 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1686 ! C
1687 ! Module name: GATHER_DATA C
1688 ! Purpose: Gather data from all processes in preparation of C
1689 ! Writing VTK files C
1690 ! C
1691 ! Author: Jeff Dietiker Date: 21-Feb-08 C
1692 ! Reviewer: Date: C
1693 ! C
1694 ! Revision Number # Date: ##-###-## C
1695 ! Author: # C
1696 ! Purpose: # C
1697 ! C
1698 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1699  SUBROUTINE gather_data
1701  USE param
1702  USE param1
1703  USE parallel
1704  USE constant
1705  USE run
1706  USE toleranc
1707  USE geometry
1708  USE indices
1709  USE compar
1710  USE mpi_utility
1711  USE sendrecv
1712  USE quadric
1713  USE cutcell
1714  USE fldvar
1715  USE vtk
1716  USE functions
1717 
1718  IMPLICIT NONE
1719 
1720  INTEGER :: IJK,I,J,K,L
1721  INTEGER :: IJK_OFFSET
1722 
1723  INTEGER :: iproc,IERR
1724  INTEGER, DIMENSION(0:numPEs-1) :: disp,rcount
1725  INTEGER, DIMENSION(:,:), ALLOCATABLE :: SHIFTED_CONNECTIVITY
1726 
1727 !======================================================================
1728 ! parallel processing
1729 !======================================================================
1730 
1731  CALL allgather_1i (number_of_new_points,rcount,ierr)
1732 
1733  IF (mype == 0) THEN
1734  ijk_offset = 0
1735  ELSE
1736  ijk_offset = 0
1737  DO iproc=0,mype-1
1738  ijk_offset = ijk_offset + rcount(iproc)
1739  ENDDO
1740  ENDIF
1741 
1742  CALL allgather_1i (ijk_offset,disp,ierr)
1743 
1744  IF(.NOT.global_var_allocated) THEN
1745 
1746  IF (mype == pe_io) THEN
1747  allocate (global_i_of(ijkmax3))
1748  allocate (global_j_of(ijkmax3))
1749  allocate (global_k_of(ijkmax3))
1750  allocate (global_connectivity(ijkmax3,15))
1751  allocate (global_number_of_nodes(ijkmax3))
1752  allocate (global_interior_cell_at(ijkmax3))
1753  allocate (global_blocked_cell_at(ijkmax3))
1754  allocate (global_standard_cell_at(ijkmax3))
1755  allocate (global_cut_cell_at(ijkmax3))
1756  allocate (global_snap(ijkmax3))
1757  allocate (global_x_new_point(ijkmax3))
1758  allocate (global_y_new_point(ijkmax3))
1759  allocate (global_z_new_point(ijkmax3))
1760  allocate (global_f_at(ijkmax3))
1761 
1762  ELSE
1763  allocate (global_i_of(1))
1764  allocate (global_j_of(1))
1765  allocate (global_k_of(1))
1766  allocate (global_connectivity(1,15))
1767  allocate (global_number_of_nodes(1))
1768  allocate (global_interior_cell_at(1))
1769  allocate (global_blocked_cell_at(1))
1770  allocate (global_standard_cell_at(1))
1771  allocate (global_cut_cell_at(1))
1772  allocate (global_snap(1))
1773  allocate (global_x_new_point(1))
1774  allocate (global_y_new_point(1))
1775  allocate (global_z_new_point(1))
1776  allocate (global_f_at(1))
1777  ENDIF
1778 
1779  global_var_allocated = .true.
1780 
1781  ENDIF
1782 
1783  IF(numpes==1) THEN ! Serial run
1787  ELSE !Parallel run
1788  call gatherv_1d( x_new_point, number_of_new_points, global_x_new_point, rcount, disp, pe_io, ierr )
1789  call gatherv_1d( y_new_point, number_of_new_points, global_y_new_point, rcount, disp, pe_io, ierr )
1790  call gatherv_1d( z_new_point, number_of_new_points, global_z_new_point, rcount, disp, pe_io, ierr )
1791  ENDIF
1792 
1794 
1795  Allocate( shifted_connectivity(dimension_3,15) )
1796 
1797  shifted_connectivity = connectivity
1798 
1799  WHERE (shifted_connectivity > ijkend3)
1800  shifted_connectivity = shifted_connectivity - ijkend3 + ijkmax3 + disp(mype)
1801  END WHERE
1802 
1803  DO ijk = ijkstart3,ijkend3
1804  DO l=1,number_of_nodes(ijk)
1805  IF(connectivity(ijk,l) <= ijkend3) THEN
1806  i = i_of(connectivity(ijk,l))
1807  j = j_of(connectivity(ijk,l))
1808  k = k_of(connectivity(ijk,l))
1809  shifted_connectivity(ijk,l) = funijk_gl(i,j,k)
1810  ENDIF
1811  ENDDO
1812  ENDDO
1813 
1814 
1815  global_interior_cell_at = .false.
1816  global_blocked_cell_at = .false.
1817  global_cut_cell_at = .false.
1818  call gather (i_of,global_i_of,root)
1819  call gather (j_of,global_j_of,root)
1820  call gather (k_of,global_k_of,root)
1821  call gather (shifted_connectivity,global_connectivity,root)
1827  call gather (snap,global_snap,root)
1828  call gather (f_at,global_f_at,root)
1829 
1830  Deallocate( shifted_connectivity )
1831 
1832  IF (mype == pe_io) THEN
1833 
1834  poly_counter = 0
1835 
1836  number_of_cells = 0
1837 
1839 
1841 
1843 
1844  DO ijk = 1, ijkmax3
1845 
1846  IF (global_interior_cell_at(ijk)) THEN
1847 
1849 
1853 
1855 
1856  ENDIF
1857 
1858  END DO
1859 
1860 
1862 
1863  ENDIF
1864 
1865  RETURN
1866 
1867 
1868  END SUBROUTINE gather_data
1869 
1870 
1871 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
1872 ! C
1873 ! Module name: SETUP_VTK_NO_CUTCELL C
1874 ! Purpose: Setup VTK data for the regular grid (no cut cells) C
1875 ! This i scalled when CARTESIAN_GRID is .FALSE. C
1876 ! and WRITE>VTK_FILES is .TRUE. . C
1877 ! C
1878 ! Author: Jeff Dietiker Date: 14-Jan-15 C
1879 ! Reviewer: Date: C
1880 ! C
1881 ! Revision Number # Date: ##-###-## C
1882 ! Author: # C
1883 ! Purpose: # C
1884 ! C
1885 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1886  SUBROUTINE setup_vtk_no_cutcell
1888  USE param
1889  USE param1
1890  USE parallel
1891  USE constant
1892  USE run
1893  USE toleranc
1894  USE geometry
1895  USE indices
1896  USE compar
1897  USE mpi_utility
1898  USE sendrecv
1899  USE quadric
1900  USE cutcell
1901  USE fldvar
1902  USE vtk
1903  USE functions
1904 
1905  IMPLICIT NONE
1906 
1907  INTEGER :: IJK,I,J,K,L,NODE
1908 
1909  INTEGER, DIMENSION(:,:), ALLOCATABLE :: SHIFTED_CONNECTIVITY
1910 
1911 ! Only a few arrays need to be allocated here simce we do not need
1912 ! all Cartesian-grid arrays
1913  IF(.NOT.ALLOCATED(xg_e)) Allocate( xg_e(0:dimension_i) )
1914  IF(.NOT.ALLOCATED(yg_n)) Allocate( yg_n(0:dimension_j) )
1915  IF(.NOT.ALLOCATED(zg_t)) Allocate( zg_t(0:dimension_k) )
1916 
1917  IF(.NOT.ALLOCATED(interior_cell_at)) THEN
1918  Allocate( interior_cell_at(dimension_3) )
1919  interior_cell_at = .false.
1920  ENDIF
1921 
1922  IF(.NOT.ALLOCATED(blocked_cell_at)) THEN
1923  Allocate( blocked_cell_at(dimension_3) )
1924  blocked_cell_at = .false.
1925  ENDIF
1926 
1927  IF(.NOT.ALLOCATED(standard_cell_at)) THEN
1928  Allocate( standard_cell_at(dimension_3) )
1929  standard_cell_at = .true.
1930  ENDIF
1931 
1932  IF(.NOT.ALLOCATED(cut_cell_at)) THEN
1933  Allocate( cut_cell_at(dimension_3) )
1934  cut_cell_at = .false.
1935  ENDIF
1936 
1937  IF(.NOT.ALLOCATED(number_of_nodes)) Allocate( number_of_nodes(dimension_3) )
1938  number_of_nodes= 0
1939 
1940  IF(.NOT.ALLOCATED(connectivity)) Allocate( connectivity(dimension_3,15) )
1941 
1942 
1943 ! This is a shoter version of Get_cut_cell_Flags
1944  DO ijk = ijkstart3, ijkend3
1945 
1946 !======================================================================
1947 ! Get coordinates of eight nodes
1948 !======================================================================
1949 
1950  CALL get_cell_node_coordinates(ijk,'SCALAR')
1951 
1952  i = i_of(ijk)
1953  j = j_of(ijk)
1954  k = k_of(ijk)
1955 
1956  IF(no_k) THEN ! 2D case
1957 
1958  interior_cell_at(ijk) = ( (i >= istart1 ).AND.(i <= iend1 ) &
1959  .AND.(j >= jstart1 ).AND.(j <= jend1 ) )
1960 
1961  ELSE ! 3D case
1962 
1963  interior_cell_at(ijk) = ( (i >= istart1 ).AND.(i <= iend1 ) &
1964  .AND.(j >= jstart1 ).AND.(j <= jend1 ) &
1965  .AND.(k >= kstart1 ).AND.(k <= kend1 ) )
1966 
1967  ENDIF
1968 
1969 
1970  IF(interior_cell_at(ijk)) THEN
1971 
1972 ! Set up the connectivity: 4 nodes in 2D, 8 nodes in 3D
1973  IF(no_k) THEN
1974  number_of_nodes(ijk) = 4
1975  connectivity(ijk,1) = ijk_of_node(5)
1976  connectivity(ijk,2) = ijk_of_node(6)
1977  connectivity(ijk,3) = ijk_of_node(8)
1978  connectivity(ijk,4) = ijk_of_node(7)
1979  ELSE
1980  number_of_nodes(ijk) = 8
1981  DO node = 1,8
1982  connectivity(ijk,node) = ijk_of_node(node)
1983  END DO
1984  ENDIF
1985 
1986 ! If obstacles are defined, they will be flagged as blocked cells
1987 ! and will not be visible in the VTK files
1988  IF(wall_at(ijk)) THEN
1989  blocked_cell_at(ijk) = .true.
1990  cut_cell_at(ijk) = .false.
1991  standard_cell_at(ijk) = .false.
1992  ENDIF
1993 
1994  ENDIF
1995 
1996  ENDDO
1997 
1998 
1999 !======================================================================
2000 ! parallel processing
2001 !======================================================================
2005 
2006  Allocate( shifted_connectivity(dimension_3,15) )
2007 
2008  shifted_connectivity = connectivity
2009 
2010 ! Replace local node index by global node index before gathering the array
2011  DO ijk = ijkstart3,ijkend3
2012  DO l=1,number_of_nodes(ijk)
2013  IF(connectivity(ijk,l) <= ijkend3) THEN
2014  i = i_of(connectivity(ijk,l))
2015  j = j_of(connectivity(ijk,l))
2016  k = k_of(connectivity(ijk,l))
2017  shifted_connectivity(ijk,l) = funijk_gl(i,j,k)
2018  ENDIF
2019  ENDDO
2020  ENDDO
2021 
2022 ! Allocate, initialize and gather arrays
2023  IF(.NOT.global_var_allocated) THEN
2024 
2025  IF (mype == pe_io) THEN
2026  allocate (global_i_of(ijkmax3))
2027  allocate (global_j_of(ijkmax3))
2028  allocate (global_k_of(ijkmax3))
2029  allocate (global_connectivity(ijkmax3,15))
2030  allocate (global_number_of_nodes(ijkmax3))
2031  allocate (global_interior_cell_at(ijkmax3))
2032  allocate (global_blocked_cell_at(ijkmax3))
2033  allocate (global_standard_cell_at(ijkmax3))
2034  allocate (global_cut_cell_at(ijkmax3))
2035 
2036  ELSE
2037  allocate (global_i_of(1))
2038  allocate (global_j_of(1))
2039  allocate (global_k_of(1))
2040  allocate (global_connectivity(1,15))
2041  allocate (global_number_of_nodes(1))
2042  allocate (global_interior_cell_at(1))
2043  allocate (global_blocked_cell_at(1))
2044  allocate (global_standard_cell_at(1))
2045  allocate (global_cut_cell_at(1))
2046  ENDIF
2047 
2048  global_var_allocated = .true.
2049 
2050  ENDIF
2051 
2052 
2053  global_interior_cell_at = .false.
2054  global_blocked_cell_at = .false.
2055  global_cut_cell_at = .false.
2056  global_standard_cell_at = .true.
2057 
2058  call gather (i_of,global_i_of,root)
2059  call gather (j_of,global_j_of,root)
2060  call gather (k_of,global_k_of,root)
2061  call gather (shifted_connectivity,global_connectivity,root)
2067 
2068  deAllocate( shifted_connectivity )
2069 
2070 
2071 
2072 ! Count the number of cells
2074 
2075  IF (mype == pe_io) THEN
2076 
2077  poly_counter = 0
2078 
2079  number_of_cells = 0
2080 
2082 
2084 
2086 
2087  DO ijk = 1, ijkmax3
2088 
2089  IF (global_interior_cell_at(ijk)) THEN
2090 
2092 
2096 
2098 
2099  ENDIF
2100 
2101  END DO
2102 
2103 ! There are no new points since there a no cut cells
2105 
2106  ENDIF
2107 
2108  RETURN
2109 
2110 
2111  END SUBROUTINE setup_vtk_no_cutcell
2112 
2113 
2114 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2115 ! C
2116 ! Module name: PRINT_GRID_STATISTICS C
2117 ! Purpose: PRINT_GRID_STATISTICS ON SCREEN C
2118 ! C
2119 ! C
2120 ! Author: Jeff Dietiker Date: 21-Feb-08 C
2121 ! Reviewer: Date: C
2122 ! C
2123 ! Revision Number # Date: ##-###-## C
2124 ! Author: # C
2125 ! Purpose: # C
2126 ! C
2127 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
2128  SUBROUTINE print_grid_statistics
2130  USE param
2131  USE param1
2132  USE parallel
2133  USE constant
2134  USE run
2135  USE toleranc
2136  USE geometry
2137  USE indices
2138  USE compar
2139  USE mpi_utility
2140  USE sendrecv
2141  USE quadric
2142  USE cutcell
2143  USE fldvar
2144  USE vtk
2145  USE functions
2146 
2147  IMPLICIT NONE
2148 
2149  INTEGER :: IJK
2150 
2151  INTEGER :: IERR
2152 
2153  DOUBLE PRECISION :: MIN_VOL, MAX_VOL, GLOBAL_MIN_VOL,GLOBAL_MAX_VOL
2154  DOUBLE PRECISION :: MIN_AYZ, MAX_AYZ, GLOBAL_MIN_AYZ,GLOBAL_MAX_AYZ
2155  DOUBLE PRECISION :: MIN_AXZ, MAX_AXZ, GLOBAL_MIN_AXZ,GLOBAL_MAX_AXZ
2156  DOUBLE PRECISION :: MIN_AXY, MAX_AXY, GLOBAL_MIN_AXY,GLOBAL_MAX_AXY
2157  DOUBLE PRECISION :: MIN_CUT, MAX_CUT, GLOBAL_MIN_CUT,GLOBAL_MAX_CUT
2158  DOUBLE PRECISION :: LOCAL_MIN_Q,LOCAL_MAX_Q, GLOBAL_MIN_Q,GLOBAL_MAX_Q
2159 
2160  IF (mype == pe_io) THEN
2161 
2162  IF(.NOT.grid_info_printed_on_screen) THEN
2163  WRITE(*,5) 'GRID STATISTICS:'
2164  WRITE(*,5) 'NUMBER OF CELLS = ', number_of_cells
2165  WRITE(*,10)'NUMBER OF STANDARD CELLS = ', &
2167  WRITE(*,10)'NUMBER OF CUT CELLS = ', &
2169  WRITE(*,10)'NUMBER OF BLOCKED CELLS = ', &
2171 
2172 5 FORMAT(1x,a,i8)
2173 10 FORMAT(1x,a,i8,' (',f6.2,' % of Total)')
2174 
2175  ENDIF
2176 
2178 
2179  ENDIF
2180 
2181 
2182 !======================================================================
2183 ! Scalar Cell volumes and areas
2184 !======================================================================
2185 
2186  min_vol = large_number
2187  max_vol = - large_number
2188  min_ayz = large_number
2189  max_ayz = - large_number
2190  min_axz = large_number
2191  max_axz = - large_number
2192  min_axy = large_number
2193  max_axy = - large_number
2194 
2195  DO ijk = ijkstart3, ijkend3
2196  IF(standard_cell_at(ijk)) THEN ! STANDARD CELLS
2197 
2198  min_vol = dmin1(min_vol,vol(ijk))
2199  max_vol = dmax1(max_vol,vol(ijk))
2200  min_ayz = dmin1(min_ayz,ayz(ijk))
2201  max_ayz = dmax1(max_ayz,ayz(ijk))
2202  min_axz = dmin1(min_axz,axz(ijk))
2203  max_axz = dmax1(max_axz,axz(ijk))
2204  min_axy = dmin1(min_axy,axy(ijk))
2205  max_axy = dmax1(max_axy,axy(ijk))
2206 
2207  ENDIF
2208  END DO
2209 
2210  call global_min(min_vol, global_min_vol, pe_io, ierr )
2211  call global_max(max_vol, global_max_vol, pe_io, ierr )
2212  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2213  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2214  call global_min(min_axz, global_min_axz, pe_io, ierr )
2215  call global_max(max_axz, global_max_axz, pe_io, ierr )
2216  call global_min(min_axy, global_min_axy, pe_io, ierr )
2217  call global_max(max_axy, global_max_axy, pe_io, ierr )
2218 
2219  IF (mype == pe_io) THEN
2220  WRITE(unit_cut_cell_log,1000) '################################################################'
2221  WRITE(unit_cut_cell_log,1000) ' CELLS STATISTICS '
2222  WRITE(unit_cut_cell_log,1000) '################################################################'
2223  WRITE(unit_cut_cell_log,1000) 'SCALAR STANDARD CELLS:'
2224  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2225  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2226  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2227  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2228  ENDIF
2229 
2230 
2231  min_vol = large_number
2232  max_vol = - large_number
2233  min_ayz = large_number
2234  max_ayz = - large_number
2235  min_axz = large_number
2236  max_axz = - large_number
2237  min_axy = large_number
2238  max_axy = - large_number
2239 
2240  DO ijk = ijkstart3, ijkend3
2241  IF(cut_cell_at(ijk)) THEN ! CUT CELLS
2242 
2243  min_vol = dmin1(min_vol,vol(ijk))
2244  max_vol = dmax1(max_vol,vol(ijk))
2245  min_ayz = dmin1(min_ayz,ayz(ijk))
2246  max_ayz = dmax1(max_ayz,ayz(ijk))
2247  min_axz = dmin1(min_axz,axz(ijk))
2248  max_axz = dmax1(max_axz,axz(ijk))
2249  min_axy = dmin1(min_axy,axy(ijk))
2250  max_axy = dmax1(max_axy,axy(ijk))
2251 
2252  ENDIF
2253  END DO
2254 
2255  call global_min(min_vol, global_min_vol, pe_io, ierr )
2256  call global_max(max_vol, global_max_vol, pe_io, ierr )
2257  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2258  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2259  call global_min(min_axz, global_min_axz, pe_io, ierr )
2260  call global_max(max_axz, global_max_axz, pe_io, ierr )
2261  call global_min(min_axy, global_min_axy, pe_io, ierr )
2262  call global_max(max_axy, global_max_axy, pe_io, ierr )
2263 
2264  IF (mype == pe_io) THEN
2265  WRITE(unit_cut_cell_log,1000) 'SCALAR CUT CELLS:'
2266  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2267  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2268  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2269  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2270  WRITE(unit_cut_cell_log,1010) 'NUMBER OF SMALL SCALAR CELLS = ', number_of_small_cells
2271  WRITE(unit_cut_cell_log,1000) '################################################################'
2272  ENDIF
2273 
2274 1000 FORMAT(a,e14.7,2x,e14.7)
2275 1010 FORMAT(a,i8)
2276 
2277 !======================================================================
2278 ! U-Momentum Cell volumes and areas
2279 !======================================================================
2280 
2281  min_vol = large_number
2282  max_vol = - large_number
2283  min_ayz = large_number
2284  max_ayz = - large_number
2285  min_axz = large_number
2286  max_axz = - large_number
2287  min_axy = large_number
2288  max_axy = - large_number
2289 
2290  DO ijk = ijkstart3, ijkend3
2291  IF(standard_u_cell_at(ijk)) THEN ! STANDARD CELLS
2292 
2293  min_vol = dmin1(min_vol,vol_u(ijk))
2294  max_vol = dmax1(max_vol,vol_u(ijk))
2295  min_ayz = dmin1(min_ayz,ayz_u(ijk))
2296  max_ayz = dmax1(max_ayz,ayz_u(ijk))
2297  min_axz = dmin1(min_axz,axz_u(ijk))
2298  max_axz = dmax1(max_axz,axz_u(ijk))
2299  min_axy = dmin1(min_axy,axy_u(ijk))
2300  max_axy = dmax1(max_axy,axy_u(ijk))
2301 
2302  ENDIF
2303  END DO
2304 
2305  call global_min(min_vol, global_min_vol, pe_io, ierr )
2306  call global_max(max_vol, global_max_vol, pe_io, ierr )
2307  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2308  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2309  call global_min(min_axz, global_min_axz, pe_io, ierr )
2310  call global_max(max_axz, global_max_axz, pe_io, ierr )
2311  call global_min(min_axy, global_min_axy, pe_io, ierr )
2312  call global_max(max_axy, global_max_axy, pe_io, ierr )
2313 
2314  IF (mype == pe_io) THEN
2315  WRITE(unit_cut_cell_log,1000) 'U-MOMENTUM STANDARD CELLS:'
2316  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2317  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2318  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2319  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2320  ENDIF
2321 
2322  min_vol = large_number
2323  max_vol = - large_number
2324  min_ayz = large_number
2325  max_ayz = - large_number
2326  min_axz = large_number
2327  max_axz = - large_number
2328  min_axy = large_number
2329  max_axy = - large_number
2330  min_cut = large_number
2331  max_cut = - large_number
2332 
2333  DO ijk = ijkstart3, ijkend3
2334  IF(cut_u_cell_at(ijk).AND.(.NOT.wall_u_at(ijk))) THEN ! CUT CELLS
2335 
2336  min_vol = dmin1(min_vol,vol_u(ijk))
2337  max_vol = dmax1(max_vol,vol_u(ijk))
2338  min_ayz = dmin1(min_ayz,ayz_u(ijk))
2339  max_ayz = dmax1(max_ayz,ayz_u(ijk))
2340  min_axz = dmin1(min_axz,axz_u(ijk))
2341  max_axz = dmax1(max_axz,axz_u(ijk))
2342  min_axy = dmin1(min_axy,axy_u(ijk))
2343  max_axy = dmax1(max_axy,axy_u(ijk))
2344  min_cut = dmin1(min_cut,area_u_cut(ijk))
2345  max_cut = dmax1(max_cut,area_u_cut(ijk))
2346 
2347  ENDIF
2348  END DO
2349 
2350  call global_min(min_vol, global_min_vol, pe_io, ierr )
2351  call global_max(max_vol, global_max_vol, pe_io, ierr )
2352  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2353  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2354  call global_min(min_axz, global_min_axz, pe_io, ierr )
2355  call global_max(max_axz, global_max_axz, pe_io, ierr )
2356  call global_min(min_axy, global_min_axy, pe_io, ierr )
2357  call global_max(max_axy, global_max_axy, pe_io, ierr )
2358  call global_min(min_cut, global_min_cut, pe_io, ierr )
2359  call global_max(max_cut, global_max_cut, pe_io, ierr )
2360 
2361  IF (mype == pe_io) THEN
2362  WRITE(unit_cut_cell_log,1000) 'U-MOMENTUM CUT CELLS:'
2363  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2364  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2365  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2366  WRITE(unit_cut_cell_log,1000) 'RANGE OF CUT AREA = ', global_min_cut,global_max_cut
2367  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2368  WRITE(unit_cut_cell_log,1010) 'NUMBER OF U WALL CELLS = ', number_of_u_wall_cells
2369  WRITE(unit_cut_cell_log,1000) '################################################################'
2370  ENDIF
2371 !======================================================================
2372 ! V-Momentum Cell volumes and areas
2373 !======================================================================
2374 
2375 
2376  min_vol = large_number
2377  max_vol = - large_number
2378  min_ayz = large_number
2379  max_ayz = - large_number
2380  min_axz = large_number
2381  max_axz = - large_number
2382  min_axy = large_number
2383  max_axy = - large_number
2384 
2385  DO ijk = ijkstart3, ijkend3
2386  IF(standard_v_cell_at(ijk)) THEN ! STANDARD CELLS
2387 
2388  min_vol = dmin1(min_vol,vol_v(ijk))
2389  max_vol = dmax1(max_vol,vol_v(ijk))
2390  min_ayz = dmin1(min_ayz,ayz_v(ijk))
2391  max_ayz = dmax1(max_ayz,ayz_v(ijk))
2392  min_axz = dmin1(min_axz,axz_v(ijk))
2393  max_axz = dmax1(max_axz,axz_v(ijk))
2394  min_axy = dmin1(min_axy,axy_v(ijk))
2395  max_axy = dmax1(max_axy,axy_v(ijk))
2396 
2397  ENDIF
2398  END DO
2399 
2400  call global_min(min_vol, global_min_vol, pe_io, ierr )
2401  call global_max(max_vol, global_max_vol, pe_io, ierr )
2402  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2403  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2404  call global_min(min_axz, global_min_axz, pe_io, ierr )
2405  call global_max(max_axz, global_max_axz, pe_io, ierr )
2406  call global_min(min_axy, global_min_axy, pe_io, ierr )
2407  call global_max(max_axy, global_max_axy, pe_io, ierr )
2408 
2409  IF (mype == pe_io) THEN
2410  WRITE(unit_cut_cell_log,1000) 'V-MOMENTUM STANDARD CELLS:'
2411  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2412  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2413  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2414  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2415  ENDIF
2416 
2417  min_vol = large_number
2418  max_vol = - large_number
2419  min_ayz = large_number
2420  max_ayz = - large_number
2421  min_axz = large_number
2422  max_axz = - large_number
2423  min_axy = large_number
2424  max_axy = - large_number
2425  min_cut = large_number
2426  max_cut = - large_number
2427 
2428  DO ijk = ijkstart3, ijkend3
2429  IF(cut_v_cell_at(ijk).AND.(.NOT.wall_v_at(ijk))) THEN ! CUT CELLS
2430 
2431  min_vol = dmin1(min_vol,vol_v(ijk))
2432  max_vol = dmax1(max_vol,vol_v(ijk))
2433  min_ayz = dmin1(min_ayz,ayz_v(ijk))
2434  max_ayz = dmax1(max_ayz,ayz_v(ijk))
2435  min_axz = dmin1(min_axz,axz_v(ijk))
2436  max_axz = dmax1(max_axz,axz_v(ijk))
2437  min_axy = dmin1(min_axy,axy_v(ijk))
2438  max_axy = dmax1(max_axy,axy_v(ijk))
2439  min_cut = dmin1(min_cut,area_v_cut(ijk))
2440  max_cut = dmax1(max_cut,area_v_cut(ijk))
2441 
2442  ENDIF
2443  END DO
2444 
2445  call global_min(min_vol, global_min_vol, pe_io, ierr )
2446  call global_max(max_vol, global_max_vol, pe_io, ierr )
2447  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2448  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2449  call global_min(min_axz, global_min_axz, pe_io, ierr )
2450  call global_max(max_axz, global_max_axz, pe_io, ierr )
2451  call global_min(min_axy, global_min_axy, pe_io, ierr )
2452  call global_max(max_axy, global_max_axy, pe_io, ierr )
2453  call global_min(min_cut, global_min_cut, pe_io, ierr )
2454  call global_max(max_cut, global_max_cut, pe_io, ierr )
2455 
2456  IF (mype == pe_io) THEN
2457  WRITE(unit_cut_cell_log,1000) 'V-MOMENTUM CUT CELLS:'
2458  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2459  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2460  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2461  WRITE(unit_cut_cell_log,1000) 'RANGE OF CUT AREA = ', global_min_cut,global_max_cut
2462  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2463  WRITE(unit_cut_cell_log,1010) 'NUMBER OF V WALL CELLS = ', number_of_v_wall_cells
2464  WRITE(unit_cut_cell_log,1000) '################################################################'
2465  ENDIF
2466 
2467 !======================================================================
2468 ! W-Momentum Cell volumes and areas
2469 !======================================================================
2470 
2471 
2472  IF(do_k) THEN
2473 
2474  min_vol = large_number
2475  max_vol = - large_number
2476  min_ayz = large_number
2477  max_ayz = - large_number
2478  min_axz = large_number
2479  max_axz = - large_number
2480  min_axy = large_number
2481  max_axy = - large_number
2482 
2483  DO ijk = ijkstart3, ijkend3
2484  IF(standard_w_cell_at(ijk)) THEN ! STANDARD CELLS
2485 
2486  min_vol = dmin1(min_vol,vol_w(ijk))
2487  max_vol = dmax1(max_vol,vol_w(ijk))
2488  min_ayz = dmin1(min_ayz,ayz_w(ijk))
2489  max_ayz = dmax1(max_ayz,ayz_w(ijk))
2490  min_axz = dmin1(min_axz,axz_w(ijk))
2491  max_axz = dmax1(max_axz,axz_w(ijk))
2492  min_axy = dmin1(min_axy,axy_w(ijk))
2493  max_axy = dmax1(max_axy,axy_w(ijk))
2494 
2495  ENDIF
2496  END DO
2497 
2498  call global_min(min_vol, global_min_vol, pe_io, ierr )
2499  call global_max(max_vol, global_max_vol, pe_io, ierr )
2500  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2501  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2502  call global_min(min_axz, global_min_axz, pe_io, ierr )
2503  call global_max(max_axz, global_max_axz, pe_io, ierr )
2504  call global_min(min_axy, global_min_axy, pe_io, ierr )
2505  call global_max(max_axy, global_max_axy, pe_io, ierr )
2506 
2507  IF (mype == pe_io) THEN
2508  WRITE(unit_cut_cell_log,1000) 'W-MOMENTUM STANDARD CELLS:'
2509  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2510  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2511  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2512  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2513  ENDIF
2514 
2515  min_vol = large_number
2516  max_vol = - large_number
2517  min_ayz = large_number
2518  max_ayz = - large_number
2519  min_axz = large_number
2520  max_axz = - large_number
2521  min_axy = large_number
2522  max_axy = - large_number
2523  min_cut = large_number
2524  max_cut = - large_number
2525 
2526  DO ijk = ijkstart3, ijkend3
2527  IF(cut_w_cell_at(ijk).AND.(.NOT.wall_w_at(ijk))) THEN ! CUT CELLS
2528 
2529  min_vol = dmin1(min_vol,vol_w(ijk))
2530  max_vol = dmax1(max_vol,vol_w(ijk))
2531  min_ayz = dmin1(min_ayz,ayz_w(ijk))
2532  max_ayz = dmax1(max_ayz,ayz_w(ijk))
2533  min_axz = dmin1(min_axz,axz_w(ijk))
2534  max_axz = dmax1(max_axz,axz_w(ijk))
2535  min_axy = dmin1(min_axy,axy_w(ijk))
2536  max_axy = dmax1(max_axy,axy_w(ijk))
2537  min_cut = dmin1(min_cut,area_w_cut(ijk))
2538  max_cut = dmax1(max_cut,area_w_cut(ijk))
2539 
2540  ENDIF
2541  END DO
2542 
2543  call global_min(min_vol, global_min_vol, pe_io, ierr )
2544  call global_max(max_vol, global_max_vol, pe_io, ierr )
2545  call global_min(min_ayz, global_min_ayz, pe_io, ierr )
2546  call global_max(max_ayz, global_max_ayz, pe_io, ierr )
2547  call global_min(min_axz, global_min_axz, pe_io, ierr )
2548  call global_max(max_axz, global_max_axz, pe_io, ierr )
2549  call global_min(min_axy, global_min_axy, pe_io, ierr )
2550  call global_max(max_axy, global_max_axy, pe_io, ierr )
2551  call global_min(min_cut, global_min_cut, pe_io, ierr )
2552  call global_max(max_cut, global_max_cut, pe_io, ierr )
2553 
2554  IF (mype == pe_io) THEN
2555  WRITE(unit_cut_cell_log,1000) 'W-MOMENTUM CUT CELLS:'
2556  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXY = ', global_min_axy,global_max_axy
2557  WRITE(unit_cut_cell_log,1000) 'RANGE OF AXZ = ', global_min_axz,global_max_axz
2558  WRITE(unit_cut_cell_log,1000) 'RANGE OF AYZ = ', global_min_ayz,global_max_ayz
2559  WRITE(unit_cut_cell_log,1000) 'RANGE OF CUT AREA = ', global_min_cut,global_max_cut
2560  WRITE(unit_cut_cell_log,1000) 'RANGE OF VOLUME = ', global_min_vol,global_max_vol
2561  WRITE(unit_cut_cell_log,1010) 'NUMBER OF W WALL CELLS = ', number_of_w_wall_cells
2562  WRITE(unit_cut_cell_log,1000) '################################################################'
2563  ENDIF
2564 
2565  ENDIF
2566 
2567 
2568 
2569  local_min_q = minval(alpha_ue_c)
2570  local_max_q = maxval(alpha_ue_c)
2571  call global_min(local_min_q, global_min_q, pe_io, ierr )
2572  call global_max(local_max_q, global_max_q, pe_io, ierr )
2573  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Alpha_Ue_c = ', global_min_q, global_max_q
2574 
2575  local_min_q = minval(alpha_un_c)
2576  local_max_q = maxval(alpha_un_c)
2577  call global_min(local_min_q, global_min_q, pe_io, ierr )
2578  call global_max(local_max_q, global_max_q, pe_io, ierr )
2579  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Alpha_Un_c = ', global_min_q, global_max_q
2580 
2581  local_min_q = minval(alpha_ut_c)
2582  local_max_q = maxval(alpha_ut_c)
2583  call global_min(local_min_q, global_min_q, pe_io, ierr )
2584  call global_max(local_max_q, global_max_q, pe_io, ierr )
2585  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Alpha_Ut_c = ', global_min_q, global_max_q
2586  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2587 
2588  local_min_q = minval(theta_ue)
2589  local_max_q = maxval(theta_ue)
2590  call global_min(local_min_q, global_min_q, pe_io, ierr )
2591  call global_max(local_max_q, global_max_q, pe_io, ierr )
2592  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Theta_Ue = ', global_min_q, global_max_q
2593 
2594  local_min_q = minval(theta_un)
2595  local_max_q = maxval(theta_un)
2596  call global_min(local_min_q, global_min_q, pe_io, ierr )
2597  call global_max(local_max_q, global_max_q, pe_io, ierr )
2598  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Theta_Un = ', global_min_q, global_max_q
2599 
2600  local_min_q = minval(theta_ut)
2601  local_max_q = maxval(theta_ut)
2602  call global_min(local_min_q, global_min_q, pe_io, ierr )
2603  call global_max(local_max_q, global_max_q, pe_io, ierr )
2604  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Theta_Ut = ', global_min_q, global_max_q
2605  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2606 
2607  local_min_q = minval(theta_u_ne)
2608  local_max_q = maxval(theta_u_ne)
2609  call global_min(local_min_q, global_min_q, pe_io, ierr )
2610  call global_max(local_max_q, global_max_q, pe_io, ierr )
2611  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Theta_U_ne = ', global_min_q, global_max_q
2612 
2613  local_min_q = minval(theta_u_te)
2614  local_max_q = maxval(theta_u_te)
2615  call global_min(local_min_q, global_min_q, pe_io, ierr )
2616  call global_max(local_max_q, global_max_q, pe_io, ierr )
2617  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM Theta_U_te = ', global_min_q, global_max_q
2618  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2619 
2620  local_min_q = minval(noc_u_e)
2621  local_max_q = maxval(noc_u_e)
2622  call global_min(local_min_q, global_min_q, pe_io, ierr )
2623  call global_max(local_max_q, global_max_q, pe_io, ierr )
2624  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM NOC_U_E = ', global_min_q, global_max_q
2625 
2626  local_min_q = minval(noc_u_n)
2627  local_max_q = maxval(noc_u_n)
2628  call global_min(local_min_q, global_min_q, pe_io, ierr )
2629  call global_max(local_max_q, global_max_q, pe_io, ierr )
2630  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM NOC_U_N = ', global_min_q, global_max_q
2631 
2632  local_min_q = minval(noc_u_t)
2633  local_max_q = maxval(noc_u_t)
2634  call global_min(local_min_q, global_min_q, pe_io, ierr )
2635  call global_max(local_max_q, global_max_q, pe_io, ierr )
2636  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM NOC_U_T = ', global_min_q, global_max_q
2637  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2638 
2639  local_min_q = minval(delh_u)
2640  local_max_q = maxval(delh_u)
2641  call global_min(local_min_q, global_min_q, pe_io, ierr )
2642  call global_max(local_max_q, global_max_q, pe_io, ierr )
2643  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF U-MOMENTUM DELH_U = ', global_min_q, global_max_q
2644  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) '################################################################'
2645 
2646 
2647 
2648  local_min_q = minval(alpha_ve_c)
2649  local_max_q = maxval(alpha_ve_c)
2650  call global_min(local_min_q, global_min_q, pe_io, ierr )
2651  call global_max(local_max_q, global_max_q, pe_io, ierr )
2652  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Alpha_Ve_c = ', global_min_q, global_max_q
2653  local_min_q = minval(alpha_vn_c)
2654  local_max_q = maxval(alpha_vn_c)
2655  call global_min(local_min_q, global_min_q, pe_io, ierr )
2656  call global_max(local_max_q, global_max_q, pe_io, ierr )
2657  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Alpha_Vn_c = ', global_min_q, global_max_q
2658  local_min_q = minval(alpha_vt_c)
2659  local_max_q = maxval(alpha_vt_c)
2660  call global_min(local_min_q, global_min_q, pe_io, ierr )
2661  call global_max(local_max_q, global_max_q, pe_io, ierr )
2662  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Alpha_Vt_c = ', global_min_q, global_max_q
2663  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2664  local_min_q = minval(theta_ve)
2665  local_max_q = maxval(theta_ve)
2666  call global_min(local_min_q, global_min_q, pe_io, ierr )
2667  call global_max(local_max_q, global_max_q, pe_io, ierr )
2668  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Theta_Ve = ', global_min_q, global_max_q
2669  local_min_q = minval(theta_vn)
2670  local_max_q = maxval(theta_vn)
2671  call global_min(local_min_q, global_min_q, pe_io, ierr )
2672  call global_max(local_max_q, global_max_q, pe_io, ierr )
2673  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Theta_Vn = ', global_min_q, global_max_q
2674  local_min_q = minval(theta_vt)
2675  local_max_q = maxval(theta_vt)
2676  call global_min(local_min_q, global_min_q, pe_io, ierr )
2677  call global_max(local_max_q, global_max_q, pe_io, ierr )
2678  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Theta_Vt = ', global_min_q, global_max_q
2679  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2680  local_min_q = minval(theta_v_ne)
2681  local_max_q = maxval(theta_v_ne)
2682  call global_min(local_min_q, global_min_q, pe_io, ierr )
2683  call global_max(local_max_q, global_max_q, pe_io, ierr )
2684  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Theta_V_ne = ', global_min_q, global_max_q
2685  local_min_q = minval(theta_v_nt)
2686  local_max_q = maxval(theta_v_nt)
2687  call global_min(local_min_q, global_min_q, pe_io, ierr )
2688  call global_max(local_max_q, global_max_q, pe_io, ierr )
2689  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM Theta_V_nt = ', global_min_q, global_max_q
2690  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2691  local_min_q = minval(noc_v_e)
2692  local_max_q = maxval(noc_v_e)
2693  call global_min(local_min_q, global_min_q, pe_io, ierr )
2694  call global_max(local_max_q, global_max_q, pe_io, ierr )
2695  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM NOC_V_E = ', global_min_q, global_max_q
2696  local_min_q = minval(noc_v_n)
2697  local_max_q = maxval(noc_v_n)
2698  call global_min(local_min_q, global_min_q, pe_io, ierr )
2699  call global_max(local_max_q, global_max_q, pe_io, ierr )
2700  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM NOC_V_N = ', global_min_q, global_max_q
2701  local_min_q = minval(noc_v_t)
2702  local_max_q = maxval(noc_v_t)
2703  call global_min(local_min_q, global_min_q, pe_io, ierr )
2704  call global_max(local_max_q, global_max_q, pe_io, ierr )
2705  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM NOC_V_T = ', global_min_q, global_max_q
2706  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2707  local_min_q = minval(delh_v)
2708  local_max_q = maxval(delh_v)
2709  call global_min(local_min_q, global_min_q, pe_io, ierr )
2710  call global_max(local_max_q, global_max_q, pe_io, ierr )
2711  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF V-MOMENTUM DELH_V = ', global_min_q, global_max_q
2712  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) '################################################################'
2713 
2714 
2715  IF(do_k) THEN
2716 
2717  local_min_q = minval(alpha_we_c)
2718  local_max_q = maxval(alpha_we_c)
2719  call global_min(local_min_q, global_min_q, pe_io, ierr )
2720  call global_max(local_max_q, global_max_q, pe_io, ierr )
2721  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Alpha_We_c = ', global_min_q, global_max_q
2722  local_min_q = minval(alpha_wn_c)
2723  local_max_q = maxval(alpha_wn_c)
2724  call global_min(local_min_q, global_min_q, pe_io, ierr )
2725  call global_max(local_max_q, global_max_q, pe_io, ierr )
2726  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Alpha_Wn_c = ', global_min_q, global_max_q
2727  local_min_q = minval(alpha_wt_c)
2728  local_max_q = maxval(alpha_wt_c)
2729  call global_min(local_min_q, global_min_q, pe_io, ierr )
2730  call global_max(local_max_q, global_max_q, pe_io, ierr )
2731  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Alpha_Wt_c = ', global_min_q, global_max_q
2732  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2733  local_min_q = minval(theta_we)
2734  local_max_q = maxval(theta_we)
2735  call global_min(local_min_q, global_min_q, pe_io, ierr )
2736  call global_max(local_max_q, global_max_q, pe_io, ierr )
2737  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Theta_We = ', global_min_q, global_max_q
2738  local_min_q = minval(theta_wn)
2739  local_max_q = maxval(theta_wn)
2740  call global_min(local_min_q, global_min_q, pe_io, ierr )
2741  call global_max(local_max_q, global_max_q, pe_io, ierr )
2742  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Theta_Wn = ', global_min_q, global_max_q
2743  local_min_q = minval(theta_wt)
2744  local_max_q = maxval(theta_wt)
2745  call global_min(local_min_q, global_min_q, pe_io, ierr )
2746  call global_max(local_max_q, global_max_q, pe_io, ierr )
2747  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Theta_Wt = ', global_min_q, global_max_q
2748  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2749  local_min_q = minval(theta_w_te)
2750  local_max_q = maxval(theta_w_te)
2751  call global_min(local_max_q, global_min_q, pe_io, ierr )
2752  call global_max(local_min_q, global_max_q, pe_io, ierr )
2753  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Theta_W_te = ', global_min_q, global_max_q
2754  local_min_q = minval(theta_w_tn)
2755  local_max_q = maxval(theta_w_tn)
2756  call global_min(local_min_q, global_min_q, pe_io, ierr )
2757  call global_max(local_max_q, global_max_q, pe_io, ierr )
2758  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM Theta_W_tn = ', global_min_q, global_max_q
2759  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2760  local_min_q = minval(noc_w_e)
2761  local_max_q = maxval(noc_w_e)
2762  call global_min(local_min_q, global_min_q, pe_io, ierr )
2763  call global_max(local_max_q, global_max_q, pe_io, ierr )
2764  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM NOC_W_E = ', global_min_q, global_max_q
2765  local_min_q = minval(noc_w_n)
2766  local_max_q = maxval(noc_w_n)
2767  call global_min(local_min_q, global_min_q, pe_io, ierr )
2768  call global_max(local_max_q, global_max_q, pe_io, ierr )
2769  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM NOC_W_N = ', global_min_q, global_max_q
2770  local_min_q = minval(noc_w_t)
2771  local_max_q = maxval(noc_w_t)
2772  call global_min(local_min_q, global_min_q, pe_io, ierr )
2773  call global_max(local_max_q, global_max_q, pe_io, ierr )
2774  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM NOC_W_T = ', global_min_q, global_max_q
2775  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000)
2776  local_min_q = minval(delh_w)
2777  local_max_q = maxval(delh_w)
2778  call global_min(local_min_q, global_min_q, pe_io, ierr )
2779  call global_max(local_max_q, global_max_q, pe_io, ierr )
2780  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) ' RANGE OF W-MOMENTUM DELH_W = ', global_min_q, global_max_q
2781  IF (mype == pe_io) WRITE(unit_cut_cell_log,1000) '################################################################'
2782 
2783  ENDIF
2784 
2785  RETURN
2786 
2787  END SUBROUTINE print_grid_statistics
2788 
2789 
2790 
2791 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2792 ! C
2793 ! Module name: CLEAN_GEOMETRY C
2794 ! Purpose: Clean-up the list of point and only keep points C
2795 ! that are used in the connectivity list. C
2796 ! C
2797 ! Author: Jeff Dietiker Date: 19-Dec-14 C
2798 ! Reviewer: Date: C
2799 ! C
2800 ! Revision Number # Date: ##-###-## C
2801 ! Author: # C
2802 ! Purpose: # C
2803 ! C
2804 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
2805 
2806  SUBROUTINE clean_geometry
2808  USE param
2809  USE param1
2810  USE parallel
2811  USE constant
2812  USE run
2813  USE toleranc
2814  USE geometry
2815  USE indices
2816  USE compar
2817  USE mpi_utility
2818  USE sendrecv
2819  USE quadric
2820  USE cutcell
2821  USE fldvar
2822  USE vtk
2823  USE cdist
2824  USE functions
2825 
2826  IMPLICIT NONE
2827 
2828  INTEGER :: IJK,L
2829 
2830  INTEGER ::POINT_ID,IJKC
2831  INTEGER , ALLOCATABLE :: POINT_NEW_ID(:)
2832  INTEGER , ALLOCATABLE :: NEW_POINT_NEW_ID(:)
2833  LOGICAL , ALLOCATABLE :: KEEP_POINT(:)
2834  LOGICAL , ALLOCATABLE :: KEEP_NEW_POINT(:)
2835 
2836 
2837  IF (mype == pe_io.AND.(.NOT.bdist_io)) THEN
2838 
2839  IF(ALLOCATED(global_cleaned_connectivity)) DEALLOCATE(global_cleaned_connectivity)
2840  IF(ALLOCATED(keep_new_point)) DEALLOCATE (keep_new_point)
2841  IF(ALLOCATED(point_new_id)) DEALLOCATE (point_new_id)
2842  IF(ALLOCATED(new_point_new_id)) DEALLOCATE (new_point_new_id)
2843  IF(ALLOCATED(keep_point)) DEALLOCATE (keep_point)
2844 
2845  ALLOCATE (global_cleaned_connectivity(ijkmax3,15))
2846  ALLOCATE (keep_new_point(global_number_of_new_points))
2847 
2848  ALLOCATE (point_new_id(ijkmax3))
2849  ALLOCATE (new_point_new_id(ijkmax3))
2850  ALLOCATE (keep_point(ijkmax3))
2851 
2852 ! Step 1: Go through connectivity list and only keep points that are used.
2853 ! For background cell corners, assign KEEP_POINT = .TRUE.
2854 ! For cut cells, the new intersection points were called NEW_POINTS,
2855 ! so assign KEEP_NEW_POINT = .TRUE.
2856 ! A NEW_POINT had an IJK index larger than IJKMAX3
2857 
2858  keep_point = .false.
2859  keep_new_point = .false.
2860 
2861  DO ijk = 1,ijkmax3
2862  IF (belongs_to_vtk_subdomain(ijk)) THEN
2863  DO l=1,global_number_of_nodes(ijk)
2864  ijkc = global_connectivity(ijk,l)
2865  IF(ijkc<=ijkmax3) keep_point(ijkc) = .true.
2866  IF(ijkc>ijkmax3) keep_new_point(ijkc-ijkmax3) = .true.
2867  ENDDO
2868  ENDIF
2869  END DO
2870 
2871 
2872 ! Step 2: Clean-up list of used points and store cleaned connectivity
2873  point_new_id = -1
2874  new_point_new_id = -1
2875  point_id = 1
2876 ! This is for the background grid cell corners
2877  DO ijk = 1,ijkmax3
2878  IF(keep_point(ijk)) THEN
2879  point_new_id(ijk) = point_id
2880  point_id = point_id + 1
2881  ENDIF
2882  END DO
2883 ! This is for the cut cell new corners
2884  DO ijk = 1,global_number_of_new_points
2885  IF(keep_new_point(ijk)) THEN
2886  new_point_new_id(ijk) = point_id
2887  point_id = point_id + 1
2888  ENDIF
2889  END DO
2890 
2891 ! Update the true (clean) number of points
2892  number_of_points = point_id - 1
2893 
2894 ! Now, store a list of coordinates for all used points
2895  IF(ALLOCATED(global_coords_of_points)) DEALLOCATE(global_coords_of_points)
2896 
2898 
2899  point_id = 1
2900 ! This is for the background grid cell corners
2901  DO ijk = 1,ijkmax3
2902  IF(keep_point(ijk)) THEN
2903  global_coords_of_points(1:3,point_id) = &
2904  (/REAL(XG_E(GLOBAL_I_OF(IJK))),REAL(YG_N(GLOBAL_J_OF(IJK))),REAL(zg_t(global_k_of(ijk)))/)
2905  point_id = point_id + 1
2906  ENDIF
2907  END DO
2908 ! This is for the cut cell new corners
2909  DO ijk = 1,global_number_of_new_points
2910  IF(keep_new_point(ijk)) THEN
2911  new_point_new_id(ijk) = point_id
2912  global_coords_of_points(1:3,point_id) = &
2913  (/REAL(GLOBAL_X_NEW_POINT(IJK)),REAL(GLOBAL_Y_NEW_POINT(IJK)),REAL(global_z_new_point(ijk))/)
2914  point_id = point_id + 1
2915  ENDIF
2916  END DO
2917 
2918 
2919 ! Step 3: Shift connectivity with new point indices
2920  DO ijk = 1,ijkmax3
2921  IF (belongs_to_vtk_subdomain(ijk)) THEN
2922  DO l=1,global_number_of_nodes(ijk)
2923  IF(global_connectivity(ijk,l)<=ijkmax3) THEN
2924  global_cleaned_connectivity(ijk,l) = point_new_id(global_connectivity(ijk,l))
2925  ELSE
2926  global_cleaned_connectivity(ijk,l) = new_point_new_id(global_connectivity(ijk,l)-ijkmax3)
2927  ENDIF
2928  ENDDO
2929  ENDIF
2930  END DO
2931 
2932 
2933 
2934  ELSEIF(bdist_io) THEN
2935 
2936 
2937  IF(ALLOCATED(cleaned_connectivity)) DEALLOCATE (cleaned_connectivity)
2938  IF(ALLOCATED(keep_new_point)) DEALLOCATE (keep_new_point)
2939  IF(ALLOCATED(point_new_id)) DEALLOCATE (point_new_id)
2940  IF(ALLOCATED(new_point_new_id)) DEALLOCATE (new_point_new_id)
2941  IF(ALLOCATED(keep_point)) DEALLOCATE (keep_point)
2942 
2943  ALLOCATE (cleaned_connectivity(ijkend3,15))
2944  ALLOCATE (keep_new_point(number_of_new_points))
2945 
2946  ALLOCATE (point_new_id(ijkend3))
2947  ALLOCATE (new_point_new_id(ijkend3))
2948  ALLOCATE (keep_point(ijkend3))
2949 
2950 ! Step 1: Go through connectivity list and only keep points that are used.
2951 ! For background cell corners, assign KEEP_POINT = .TRUE.
2952 ! For cut cells, the new intersection points were called NEW_POINTS,
2953 ! so assign KEEP_NEW_POINT = .TRUE.
2954 ! A NEW_POINT had an IJK index larger than IJKMAX3
2955 
2956  keep_point = .false.
2957  keep_new_point = .false.
2958 
2959  DO ijk = 1,ijkend3
2960  IF (belongs_to_vtk_subdomain(ijk)) THEN
2961  DO l=1,number_of_nodes(ijk)
2962  ijkc = connectivity(ijk,l)
2963  IF(ijkc<=ijkend3) keep_point(ijkc) = .true.
2964  IF(ijkc>ijkend3) keep_new_point(ijkc-ijkend3) = .true.
2965  ENDDO
2966  ENDIF
2967  END DO
2968 
2969 
2970 ! Step 2: Clean-up list of used points and store cleaned connectivity
2971  point_new_id = -1
2972  new_point_new_id = -1
2973  point_id = 1
2974 ! This is for the background grid cell corners
2975  DO ijk = 1,ijkend3
2976  IF(keep_point(ijk)) THEN
2977  point_new_id(ijk) = point_id
2978  point_id = point_id + 1
2979  ENDIF
2980  END DO
2981 ! This is for the cut cell new corners
2982  DO ijk = 1,number_of_new_points
2983  IF(keep_new_point(ijk)) THEN
2984  new_point_new_id(ijk) = point_id
2985  point_id = point_id + 1
2986  ENDIF
2987  END DO
2988 
2989 ! Update the true (clean) number of points
2990  number_of_points = point_id - 1
2991 
2992 ! Now, store a list of coordinates for all used points
2993  IF(ALLOCATED(coords_of_points)) DEALLOCATE(coords_of_points)
2994 
2995  ALLOCATE(coords_of_points(number_of_points,3))
2996 
2997  point_id = 1
2998 ! This is for the background grid cell corners
2999  DO ijk = 1,ijkend3
3000  IF(keep_point(ijk)) THEN
3001  coords_of_points(point_id,1:3) = &
3002  (/REAL(XG_E(I_OF(IJK))),REAL(YG_N(J_OF(IJK))),REAL(zg_t(k_of(ijk)))/)
3003  point_id = point_id + 1
3004  ENDIF
3005  END DO
3006 ! This is for the cut cell new corners
3007  DO ijk = 1,number_of_new_points
3008  IF(keep_new_point(ijk)) THEN
3009  new_point_new_id(ijk) = point_id
3010  coords_of_points(point_id,1:3) = &
3011  (/REAL(X_NEW_POINT(IJK)),REAL(Y_NEW_POINT(IJK)),REAL(z_new_point(ijk))/)
3012  point_id = point_id + 1
3013  ENDIF
3014  END DO
3015 
3016 
3017 ! Step 3: Shift connectivity with new point indices
3018  DO ijk = 1,ijkend3
3019  IF (belongs_to_vtk_subdomain(ijk)) THEN
3020  DO l=1,number_of_nodes(ijk)
3021  IF(connectivity(ijk,l)<=ijkend3) THEN
3022  cleaned_connectivity(ijk,l) = point_new_id(connectivity(ijk,l))
3023  ELSE
3024  cleaned_connectivity(ijk,l) = new_point_new_id(connectivity(ijk,l)-ijkend3)
3025  ENDIF
3026  ENDDO
3027  ENDIF
3028  END DO
3029 
3030  ENDIF
3031 
3032  RETURN
3033 
3034  END SUBROUTINE clean_geometry
3035 
3036 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
3037 ! C
3038 ! Module name: SETUP_VTK_REGION C
3039 ! C
3040 ! Purpose: Filter the cells based on the VTK region bounds and C
3041 ! set the flag BELONGS_TO_VTK_SUBDOMAIN(IJK) to .TRUE. C
3042 ! to keep the cell. C
3043 ! C
3044 ! Author: Jeff Dietiker Date: 19-Dec-14 C
3045 ! Reviewer: Date: C
3046 ! C
3047 ! Revision Number # Date: ##-###-## C
3048 ! Author: # C
3049 ! Purpose: # C
3050 ! C
3051 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
3052  SUBROUTINE setup_vtk_region
3054  USE cdist
3055  USE compar, only: mype, pe_io, ijkend3
3056  USE cutcell
3057  USE geometry
3058  USE indices, only: i_of, j_of, k_of
3059  USE vtk
3060 
3061  IMPLICIT NONE
3062 
3063  INTEGER :: IJK,I,J,K,I_E,I_W,J_N,J_S,K_T,K_B
3064  INTEGER :: NXS,NYS,NZS,NS,I_TMP,J_TMP,K_TMP
3065  INTEGER :: I_SLICE(dim_i),J_SLICE(dim_j),K_SLICE(dim_k)
3066  DOUBLE PRECISION :: XE,XW,YS,YN,ZB,ZT
3067  DOUBLE PRECISION :: XSLICE,YSLICE,ZSLICE
3068  LOGICAL :: KEEP_XDIR,KEEP_YDIR,KEEP_ZDIR
3069 
3070 ! Get VTK region bounds
3071  xe = vtk_x_e(vtk_region)
3072  xw = vtk_x_w(vtk_region)
3073  ys = vtk_y_s(vtk_region)
3074  yn = vtk_y_n(vtk_region)
3075  zb = vtk_z_b(vtk_region)
3076  zt = vtk_z_t(vtk_region)
3077 
3078  nxs = vtk_nxs(vtk_region)
3079  nys = vtk_nys(vtk_region)
3080  nzs = vtk_nzs(vtk_region)
3081 
3082  CALL calc_cell (xmin, vtk_x_w(vtk_region), dx, imax, i_w)
3083  i_w = i_w !+ 1
3084  CALL calc_cell (xmin, vtk_x_e(vtk_region), dx, imax, i_e)
3085 
3086 
3087  CALL calc_cell (zero, vtk_y_s(vtk_region), dy, jmax, j_s)
3088  j_s = j_s !+ 1
3089  CALL calc_cell (zero, vtk_y_n(vtk_region), dy, jmax, j_n)
3090 
3091  IF (no_k) THEN
3092  k_b = 1
3093  k_t = 1
3094  ELSE
3095  CALL calc_cell (zero, vtk_z_b(vtk_region), dz, kmax, k_b)
3096  k_b = k_b !+ 1
3097  CALL calc_cell (zero, vtk_z_t(vtk_region), dz, kmax, k_t)
3098  ENDIF
3099 
3100 ! get slice(s) location
3101  DO ns = 1,nxs
3102  xslice = xw + (xe-xw)/(nxs-1)*(ns-1)
3103  CALL calc_cell (xmin, xslice, dx, imax, i_tmp)
3104  i_slice(ns) = max(min(i_tmp,imax1),imin1)
3105  ENDDO
3106 
3107  DO ns = 1,nys
3108  yslice = ys + (yn-ys)/(nys-1)*(ns-1)
3109  CALL calc_cell (zero, yslice, dy, jmax, j_tmp)
3110  j_slice(ns) = max(min(j_tmp,jmax1),jmin1)
3111  ENDDO
3112 
3113  DO ns = 1,nzs
3114  zslice = zb + (zt-zb)/(nzs-1)*(ns-1)
3115  CALL calc_cell (zero, zslice, dz, kmax, k_tmp)
3116  k_slice(ns) = max(min(k_tmp,kmax1),kmin1)
3117  ENDDO
3118 
3119  IF (mype == pe_io.AND.(.NOT.bdist_io)) THEN
3120 
3121  IF(ALLOCATED(belongs_to_vtk_subdomain)) DEALLOCATE(belongs_to_vtk_subdomain)
3122 
3123  ALLOCATE (belongs_to_vtk_subdomain(ijkmax3))
3124 
3125 ! Filter the cells based on the VTK region bounds and set the
3126 ! flag BELONGS_TO_VTK_SUBDOMAIN(IJK) to .TRUE. to keep the cell.
3127 
3128  belongs_to_vtk_subdomain = .false.
3130 
3131  DO ijk = 1,ijkmax3
3132  IF (global_interior_cell_at(ijk)) THEN
3133  IF (.NOT.global_blocked_cell_at(ijk)) THEN
3134  i = global_i_of(ijk)
3135  j = global_j_of(ijk)
3136  k = global_k_of(ijk)
3137 
3138  IF(vtk_cutcell_only(vtk_region)) THEN
3139  IF(i==imin1.OR.i==imax1.OR. &
3140  j==jmin1.OR.j==jmax1.OR. &
3141  k==kmin1.OR.k==kmax1.OR. &
3142  global_cut_cell_at(ijk)) THEN
3143 
3144  belongs_to_vtk_subdomain(ijk) = .true.
3146  ENDIF
3147  cycle
3148  ENDIF
3149 
3150 
3151 ! X-direction
3152  keep_xdir=.false.
3153  IF(nxs==0) THEN
3154  IF(i_w<=i.AND.i<=i_e) keep_xdir=.true.
3155  ELSE
3156  DO ns = 1,nxs
3157  IF(i==i_slice(ns)) keep_xdir=.true.
3158  ENDDO
3159  ENDIF
3160 
3161 ! Y-direction
3162  keep_ydir=.false.
3163  IF(nys==0) THEN
3164  IF(j_s<=j.AND.j<=j_n) keep_ydir=.true.
3165  ELSE
3166  DO ns = 1,nys
3167  IF(j==j_slice(ns)) keep_ydir=.true.
3168  ENDDO
3169  ENDIF
3170 
3171 ! Z-direction
3172  keep_zdir=.false.
3173  IF(nzs==0) THEN
3174  IF(k_b<=k.AND.k<=k_t) keep_zdir=.true.
3175  ELSE
3176  DO ns = 1,nzs
3177  IF(k==k_slice(ns)) keep_zdir=.true.
3178  ENDDO
3179  ENDIF
3180 
3181 ! Now combine
3182  IF(keep_xdir.AND.keep_ydir.AND.keep_zdir) THEN
3183  belongs_to_vtk_subdomain(ijk) = .true.
3185  ENDIF
3186  ENDIF
3187  ENDIF
3188  END DO
3189 
3190  ELSE ! BDIST_IO
3191 
3192  IF(ALLOCATED(belongs_to_vtk_subdomain)) DEALLOCATE(belongs_to_vtk_subdomain)
3193 
3194  ALLOCATE (belongs_to_vtk_subdomain(ijkend3))
3195 
3196 ! Filter the cells based on the VTK region bounds and set the
3197 ! flag BELONGS_TO_VTK_SUBDOMAIN(IJK) to .TRUE. to keep the cell.
3198 
3199  belongs_to_vtk_subdomain = .false.
3201 
3202  DO ijk = 1,ijkend3
3203  IF (interior_cell_at(ijk)) THEN
3204  IF (.NOT.blocked_cell_at(ijk)) THEN
3205  i = i_of(ijk)
3206  j = j_of(ijk)
3207  k = k_of(ijk)
3208 
3209 ! X-direction
3210  keep_xdir=.false.
3211  IF(nxs==0) THEN
3212  IF(i_w<=i.AND.i<=i_e) keep_xdir=.true.
3213  ELSE
3214  DO ns = 1,nxs
3215  IF(i==i_slice(ns)) keep_xdir=.true.
3216  ENDDO
3217  ENDIF
3218 
3219 ! Y-direction
3220  keep_ydir=.false.
3221  IF(nys==0) THEN
3222  IF(j_s<=j.AND.j<=j_n) keep_ydir=.true.
3223  ELSE
3224  DO ns = 1,nys
3225  IF(j==j_slice(ns)) keep_ydir=.true.
3226  ENDDO
3227  ENDIF
3228 
3229 ! Z-direction
3230  keep_zdir=.false.
3231  IF(nzs==0) THEN
3232  IF(k_b<=k.AND.k<=k_t) keep_zdir=.true.
3233  ELSE
3234  DO ns = 1,nzs
3235  IF(k==k_slice(ns)) keep_zdir=.true.
3236  ENDDO
3237  ENDIF
3238 
3239 ! Now combine
3240  IF(keep_xdir.AND.keep_ydir.AND.keep_zdir) THEN
3241  belongs_to_vtk_subdomain(ijk) = .true.
3243  ENDIF
3244  ENDIF
3245  ENDIF
3246  END DO
3247 
3248  ENDIF
3249 
3250  RETURN
3251 
3252  END SUBROUTINE setup_vtk_region
double precision, dimension(:), allocatable theta_wn
Definition: cutcell_mod.f:302
character(len=255) vtu_filename
Definition: vtk_mod.f:28
integer dimension_i
Definition: param_mod.f:7
logical, dimension(dimension_vtk) vtk_v_g
Definition: vtk_mod.f:162
double precision, dimension(:), allocatable theta_un
Definition: cutcell_mod.f:222
logical, dimension(:), allocatable global_blocked_cell_at
Definition: vtk_mod.f:56
double precision, dimension(:,:), allocatable v_s
Definition: fldvar_mod.f:105
logical re_indexing
Definition: cutcell_mod.f:16
logical, dimension(dimension_vtk) vtk_defined
Definition: vtk_mod.f:108
subroutine calc_cell(RMIN, REACTOR_LOC, D_DIR, N_DIR, CELL_LOC)
Definition: calc_cell.f:14
double precision, dimension(:), allocatable vol_w
Definition: geometry_mod.f:242
double precision, dimension(:), allocatable alpha_ut_c
Definition: cutcell_mod.f:231
double precision, dimension(:), allocatable theta_w_tn
Definition: cutcell_mod.f:290
integer imax2
Definition: geometry_mod.f:61
double precision, dimension(:), allocatable noc_v_t
Definition: cutcell_mod.f:270
integer, dimension(:), allocatable i_of
Definition: indices_mod.f:45
character(len=255) vtu_dir
Definition: vtk_mod.f:27
logical, dimension(:), allocatable wall_u_at
Definition: cutcell_mod.f:126
logical, dimension(dimension_vtk, dim_m) vtk_w_s
Definition: vtk_mod.f:171
integer number_of_cut_cells
Definition: vtk_mod.f:18
logical, dimension(:), allocatable cut_u_cell_at
Definition: cutcell_mod.f:356
logical, dimension(dimension_vtk, dim_m) vtk_u_s
Definition: vtk_mod.f:159
logical, dimension(:), allocatable standard_u_cell_at
Definition: cutcell_mod.f:369
subroutine setup_vtk_no_cutcell
Definition: vtk_out.f:1887
double precision, dimension(:,:), allocatable debug_cg
Definition: cutcell_mod.f:444
integer ijkend3
Definition: compar_mod.f:80
integer kend1
Definition: compar_mod.f:80
logical time_dependent_filename
Definition: vtk_mod.f:25
double precision, dimension(:), allocatable alpha_vt_c
Definition: cutcell_mod.f:269
double precision, dimension(:), allocatable yg_n
Definition: cutcell_mod.f:45
real, dimension(:,:), allocatable global_coords_of_points
Definition: vtk_mod.f:53
subroutine open_pvd_file
Definition: vtk_out.f:1299
double precision, dimension(:), allocatable ep_g
Definition: fldvar_mod.f:15
integer istart1
Definition: compar_mod.f:80
logical bdist_io
Definition: cdist_mod.f:4
double precision, dimension(:), allocatable global_f_at
Definition: vtk_mod.f:60
double precision, dimension(:), allocatable theta_u_ne
Definition: cutcell_mod.f:213
logical, dimension(:), allocatable global_interior_cell_at
Definition: vtk_mod.f:55
integer, dimension(:,:), allocatable cleaned_connectivity
Definition: vtk_mod.f:43
double precision, dimension(:), allocatable k_turb_g
Definition: fldvar_mod.f:161
double precision, dimension(:), allocatable global_x_new_point
Definition: vtk_mod.f:62
integer iend1
Definition: compar_mod.f:80
subroutine update_and_close_pvd_file
Definition: vtk_out.f:1403
Definition: pgcor_mod.f:1
integer, dimension(:,:), allocatable global_connectivity
Definition: vtk_mod.f:49
logical, dimension(dimension_vtk) vtk_t_g
Definition: vtk_mod.f:186
logical, dimension(dimension_vtk) vtk_dbg_file
Definition: vtk_mod.f:93
integer dimension_3
Definition: param_mod.f:11
integer number_of_vtk_cells
Definition: vtk_mod.f:21
double precision, dimension(:), allocatable xg_e
Definition: cutcell_mod.f:44
subroutine gatherv_1d(lbuf, sendcnt, gbuf, rcount, disp, mroot, idebug)
subroutine reorder_polygon(NP, COORD, NORMAL, IERROR)
subroutine write_vtp_file(LCV, MODE)
Definition: vtp_mod.f:573
double precision, dimension(:), allocatable delh_u
Definition: cutcell_mod.f:200
integer pvtu_unit
Definition: vtk_mod.f:36
double precision, dimension(:), allocatable axy
Definition: geometry_mod.f:210
integer dimension_k
Definition: param_mod.f:9
integer, dimension(:), allocatable global_number_of_nodes
Definition: vtk_mod.f:51
double precision, dimension(dimension_vtk) vtk_x_e
Definition: vtk_mod.f:114
double precision fac_dim_max_cut_cell
Definition: cutcell_mod.f:28
double precision, dimension(:), allocatable noc_w_n
Definition: cutcell_mod.f:306
double precision, dimension(0:15) z_node
Definition: cutcell_mod.f:76
double precision, dimension(:,:), allocatable w_s
Definition: fldvar_mod.f:117
logical reset_frame_at_time_zero
Definition: vtk_mod.f:26
Definition: rxns_mod.f:1
subroutine gather_data
Definition: vtk_out.f:1700
character(len=60) run_name
Definition: run_mod.f:24
double precision, dimension(:), allocatable noc_v_n
Definition: cutcell_mod.f:264
subroutine send_receive_1d_logical(L1D, NLAYERS)
double precision, dimension(:), allocatable alpha_ve_c
Definition: cutcell_mod.f:260
logical, dimension(dimension_vtk) vtk_p_star
Definition: vtk_mod.f:153
double precision, dimension(0:15) y_node
Definition: cutcell_mod.f:75
double precision, dimension(:), allocatable t_g
Definition: fldvar_mod.f:63
logical, dimension(:), allocatable wall_v_at
Definition: cutcell_mod.f:127
subroutine open_vtu_file_bin(MODE)
Definition: vtk_out.f:355
double precision, dimension(:), allocatable global_y_new_point
Definition: vtk_mod.f:63
integer number_of_w_wall_cells
Definition: cutcell_mod.f:403
logical, dimension(:), allocatable standard_w_cell_at
Definition: cutcell_mod.f:371
integer, dimension(:), allocatable global_j_of
Definition: vtk_mod.f:47
double precision, dimension(:,:), allocatable scalar
Definition: fldvar_mod.f:155
double precision, dimension(0:dim_j) dy
Definition: geometry_mod.f:70
double precision, dimension(:), allocatable theta_wt
Definition: cutcell_mod.f:293
integer, dimension(dimension_vtk) vtk_nxs
Definition: vtk_mod.f:129
subroutine write_geometry_in_vtu_bin(PASS)
Definition: vtk_out.f:552
subroutine allgather_1i(lbuf, gbuf, idebug)
double precision, dimension(:), allocatable global_z_new_point
Definition: vtk_mod.f:64
double precision, dimension(:), allocatable p_g
Definition: fldvar_mod.f:26
Definition: vtk_mod.f:1
double precision, dimension(:), allocatable noc_u_n
Definition: cutcell_mod.f:226
double precision, dimension(:,:), allocatable normal_s
Definition: cutcell_mod.f:120
double precision, dimension(:), allocatable ayz_u
Definition: geometry_mod.f:218
logical, dimension(dimension_vtk) pvd_file_initialized
Definition: vtk_mod.f:105
logical full_log
Definition: output_mod.f:31
double precision, dimension(dimension_vtk) vtk_y_n
Definition: vtk_mod.f:120
integer mpierr
Definition: compar_mod.f:27
character(len=512) buffer
Definition: vtk_mod.f:30
double precision, dimension(0:dim_k) dz
Definition: geometry_mod.f:72
double precision, dimension(:), allocatable ayz
Definition: geometry_mod.f:206
integer number_of_new_points
Definition: cutcell_mod.f:81
character(len=32), dimension(dim_n_g) species_alias_g
Definition: rxns_mod.f:48
double precision, dimension(:), allocatable delh_v
Definition: cutcell_mod.f:238
logical, dimension(:), allocatable global_snap
Definition: vtk_mod.f:59
double precision, dimension(:), allocatable a
Definition: scalars_mod.f:29
logical, dimension(dimension_vtk) vtk_u_g
Definition: vtk_mod.f:156
double precision, dimension(:), allocatable alpha_ue_c
Definition: cutcell_mod.f:219
double precision, dimension(:,:), allocatable u_s
Definition: fldvar_mod.f:93
integer kstart1
Definition: compar_mod.f:80
logical, dimension(dimension_vtk) vtk_dwall
Definition: vtk_mod.f:223
double precision, dimension(:), allocatable ayz_v
Definition: geometry_mod.f:227
integer numpes
Definition: compar_mod.f:24
integer imax
Definition: geometry_mod.f:47
double precision, dimension(:), allocatable axz_u
Definition: geometry_mod.f:220
integer ijkmax3
Definition: geometry_mod.f:82
double precision, dimension(:), allocatable theta_ut
Definition: cutcell_mod.f:228
character(len=1), dimension(dimension_vtk) vtk_data
Definition: vtk_mod.f:99
double precision, dimension(:), allocatable alpha_wt_c
Definition: cutcell_mod.f:308
integer, dimension(:), allocatable k_of
Definition: indices_mod.f:47
integer, dimension(dimension_vtk) frame
Definition: vtk_mod.f:102
integer pe_io
Definition: compar_mod.f:30
double precision, dimension(:), allocatable f_at
Definition: cutcell_mod.f:486
double precision, dimension(:), allocatable alpha_wn_c
Definition: cutcell_mod.f:305
logical, dimension(dimension_vtk, vtk_nrrmax) vtk_rrate
Definition: vtk_mod.f:204
subroutine get_global_cell_node_coordinates(IJK, TYPE_OF_CELL)
integer poly_counter
Definition: vtk_mod.f:41
integer pvd_unit
Definition: vtk_mod.f:35
double precision, dimension(:), allocatable partition
Definition: cutcell_mod.f:117
logical, dimension(dimension_vtk, dim_m) vtk_t_s
Definition: vtk_mod.f:189
logical, dimension(dimension_vtk, dim_m) vtk_vel_s
Definition: vtk_mod.f:177
logical use_rrates
Definition: rxns_mod.f:21
logical, dimension(dimension_vtk) vtk_ep_g
Definition: vtk_mod.f:147
integer kmax1
Definition: geometry_mod.f:58
integer, dimension(:), allocatable global_i_of
Definition: vtk_mod.f:46
logical, dimension(dimension_vtk, dim_scalar) vtk_scalar
Definition: vtk_mod.f:201
double precision, dimension(:), allocatable noc_w_e
Definition: cutcell_mod.f:300
double precision, dimension(:), allocatable area_v_cut
Definition: cutcell_mod.f:133
subroutine mfix_exit(myID, normal_termination)
Definition: exit.f:5
logical use_stl
Definition: cutcell_mod.f:428
integer mmax
Definition: physprop_mod.f:19
logical, dimension(dimension_vtk) vtk_lambda_2
Definition: vtk_mod.f:214
double precision, dimension(:,:,:), allocatable x_s
Definition: fldvar_mod.f:78
integer, dimension(:), allocatable j_of
Definition: indices_mod.f:46
integer imax1
Definition: geometry_mod.f:54
logical, dimension(:), allocatable wall_w_at
Definition: cutcell_mod.f:128
integer vtk_region
Definition: vtk_mod.f:86
double precision, dimension(:), allocatable theta_v_ne
Definition: cutcell_mod.f:248
double precision, dimension(:), allocatable axy_v
Definition: geometry_mod.f:231
integer jmax2
Definition: geometry_mod.f:63
double precision, dimension(:,:), allocatable t_s
Definition: fldvar_mod.f:66
logical, dimension(dimension_vtk) vtk_w_g
Definition: vtk_mod.f:168
double precision, dimension(:,:), allocatable x_g
Definition: fldvar_mod.f:75
Definition: stl_mod.f:1
logical, dimension(dimension_vtk) vtk_k_turb_g
Definition: vtk_mod.f:207
integer root
Definition: compar_mod.f:41
subroutine write_vtu_file(LCV, MODE)
Definition: vtk_out.f:46
Definition: exit.f:2
Definition: cdist_mod.f:2
integer number_of_blocked_cells
Definition: vtk_mod.f:19
double precision, dimension(:), allocatable alpha_un_c
Definition: cutcell_mod.f:225
double precision, dimension(dimension_vtk) vtk_y_s
Definition: vtk_mod.f:117
subroutine get_cell_node_coordinates(IJK, TYPE_OF_CELL)
integer vtu_offset
Definition: vtk_mod.f:76
double precision, dimension(:), allocatable noc_u_t
Definition: cutcell_mod.f:232
character(len=16) run_type
Definition: run_mod.f:33
subroutine print_grid_statistics
Definition: vtk_out.f:2129
logical, dimension(dimension_vtk, dim_n_g) vtk_x_g
Definition: vtk_mod.f:192
double precision, dimension(:,:), allocatable theta_m
Definition: fldvar_mod.f:149
double precision, dimension(:), allocatable theta_vn
Definition: cutcell_mod.f:251
double precision, dimension(:), allocatable alpha_we_c
Definition: cutcell_mod.f:299
subroutine close_vtu_file_bin(MODE)
Definition: vtk_out.f:1212
double precision, dimension(:), allocatable v_g
Definition: fldvar_mod.f:99
logical, dimension(dimension_vtk, dim_m) vtk_theta_m
Definition: vtk_mod.f:198
double precision tol_f
Definition: quadric_mod.f:76
integer number_of_u_wall_cells
Definition: cutcell_mod.f:401
integer, dimension(:,:), allocatable global_cleaned_connectivity
Definition: vtk_mod.f:50
logical, dimension(:), allocatable cut_w_cell_at
Definition: cutcell_mod.f:358
double precision, dimension(0:dim_i) dx
Definition: geometry_mod.f:68
subroutine clean_geometry
Definition: vtk_out.f:2807
logical, dimension(dimension_vtk, dim_m, dim_n_s) vtk_x_s
Definition: vtk_mod.f:195
logical do_j
Definition: geometry_mod.f:26
integer n_usr_def
Definition: cutcell_mod.f:424
double precision, dimension(:), allocatable noc_u_e
Definition: cutcell_mod.f:220
logical, dimension(dimension_vtk, 15) vtk_debug
Definition: vtk_mod.f:238
double precision, dimension(:), allocatable theta_ue
Definition: cutcell_mod.f:210
double precision, dimension(:), allocatable z_new_point
Definition: cutcell_mod.f:84
character(len=1), parameter end_rec
Definition: vtk_mod.f:32
integer kmax2
Definition: geometry_mod.f:65
double precision, dimension(:), allocatable y_new_point
Definition: cutcell_mod.f:83
logical global_var_allocated
Definition: vtk_mod.f:69
double precision, dimension(:), allocatable w_g
Definition: fldvar_mod.f:111
integer number_of_cells
Definition: vtk_mod.f:17
subroutine write_scalar_in_vtu_bin(VAR_NAME, VAR, PASS)
Definition: vtk_out.f:895
double precision, dimension(:), allocatable x_new_point
Definition: cutcell_mod.f:82
Definition: pscor_mod.f:1
integer, dimension(:), allocatable number_of_nodes
Definition: cutcell_mod.f:105
integer jmax1
Definition: geometry_mod.f:56
Definition: run_mod.f:13
logical, dimension(dimension_vtk) vtk_bc_id
Definition: vtk_mod.f:220
integer nrr
Definition: rxns_mod.f:10
double precision, dimension(:), allocatable theta_v_nt
Definition: cutcell_mod.f:254
double precision, dimension(:), allocatable axz
Definition: geometry_mod.f:208
double precision, dimension(:), allocatable ayz_w
Definition: geometry_mod.f:236
double precision, parameter large_number
Definition: param1_mod.f:23
Definition: vtp_mod.f:1
Definition: param_mod.f:2
double precision, dimension(:), allocatable delh_w
Definition: cutcell_mod.f:277
logical, dimension(dimension_vtk) vtk_vel_g
Definition: vtk_mod.f:174
integer kmax
Definition: geometry_mod.f:51
logical no_k
Definition: geometry_mod.f:28
logical, dimension(:), allocatable cut_v_cell_at
Definition: cutcell_mod.f:357
logical, dimension(dimension_vtk) vtk_p_g
Definition: vtk_mod.f:150
integer number_of_points
Definition: vtk_mod.f:41
double precision, dimension(:), allocatable vorticity
Definition: cutcell_mod.f:406
logical, dimension(dimension_vtk, dim_m) vtk_rop_s
Definition: vtk_mod.f:180
logical, dimension(:), allocatable interior_cell_at
Definition: cutcell_mod.f:40
integer, dimension(0:dim_m) nmax
Definition: physprop_mod.f:119
logical, dimension(dimension_vtk) vtk_partition
Definition: vtk_mod.f:217
integer jmin1
Definition: geometry_mod.f:42
subroutine write_dbg_vtu_and_vtp_files
Definition: vtk_out.f:16
logical do_k
Definition: geometry_mod.f:30
logical, dimension(:), allocatable cut_cell_at
Definition: cutcell_mod.f:355
logical, dimension(dimension_vtk, dim_m) vtk_v_s
Definition: vtk_mod.f:165
double precision, dimension(:,:), allocatable p_s
Definition: fldvar_mod.f:123
logical, dimension(:), allocatable global_standard_cell_at
Definition: vtk_mod.f:57
logical k_epsilon
Definition: run_mod.f:97
integer mype
Definition: compar_mod.f:24
double precision, dimension(:,:), allocatable reactionrates
Definition: rxns_mod.f:7
double precision, dimension(:), allocatable axz_w
Definition: geometry_mod.f:238
integer dimension_max_cut_cell
Definition: cutcell_mod.f:25
integer ijkstart3
Definition: compar_mod.f:80
double precision, dimension(:), allocatable noc_w_t
Definition: cutcell_mod.f:309
logical use_msh
Definition: cutcell_mod.f:430
integer number_of_small_cells
Definition: cutcell_mod.f:399
subroutine calc_vorticity
Definition: calc_vort_out.f:15
logical, dimension(dimension_vtk) vtk_cutcell_only
Definition: vtk_mod.f:141
double precision, dimension(:), allocatable axy_w
Definition: geometry_mod.f:240
integer global_number_of_new_points
Definition: vtk_mod.f:66
double precision, dimension(dimension_vtk) vtk_z_b
Definition: vtk_mod.f:123
double precision, dimension(:), allocatable theta_w_te
Definition: cutcell_mod.f:287
integer n_polygon
Definition: polygon_mod.f:6
integer, dimension(dimension_vtk) vtk_nys
Definition: vtk_mod.f:132
integer, dimension(:), allocatable bc_id
Definition: cutcell_mod.f:433
double precision, dimension(:), allocatable u_g
Definition: fldvar_mod.f:87
integer nscalar
Definition: scalars_mod.f:7
double precision xmin
Definition: geometry_mod.f:75
integer number_of_v_wall_cells
Definition: cutcell_mod.f:402
double precision, dimension(:), allocatable vol_u
Definition: geometry_mod.f:224
integer jmax
Definition: geometry_mod.f:49
logical, dimension(:), allocatable standard_cell_at
Definition: cutcell_mod.f:368
logical, dimension(:), allocatable standard_v_cell_at
Definition: cutcell_mod.f:370
double precision, dimension(:), allocatable area_w_cut
Definition: cutcell_mod.f:134
logical do_i
Definition: geometry_mod.f:22
logical, dimension(dimension_vtk) vtk_e_turb_g
Definition: vtk_mod.f:208
logical, dimension(dimension_vtk) vtk_vorticity
Definition: vtk_mod.f:211
double precision, dimension(:,:), allocatable rop_s
Definition: fldvar_mod.f:51
integer, dimension(:), allocatable global_k_of
Definition: vtk_mod.f:48
double precision, dimension(:), allocatable theta_vt
Definition: cutcell_mod.f:266
character(len=255), dimension(dimension_vtk) vtk_filebase
Definition: vtk_mod.f:144
double precision, dimension(:), allocatable area_u_cut
Definition: cutcell_mod.f:132
subroutine setup_vtk_region
Definition: vtk_out.f:3053
double precision, dimension(:), allocatable theta_we
Definition: cutcell_mod.f:296
double precision, dimension(:), allocatable dwall
Definition: cutcell_mod.f:488
integer number_of_standard_cells
Definition: vtk_mod.f:20
integer vtu_frame_unit
Definition: vtk_mod.f:36
logical grid_info_printed_on_screen
Definition: vtk_mod.f:71
double precision, dimension(dimension_vtk) vtk_z_t
Definition: vtk_mod.f:126
double precision, dimension(:), allocatable lambda2
Definition: cutcell_mod.f:407
subroutine eval_f(METHOD, x1, x2, x3, Q, f, CLIP_FLAG)
logical, dimension(:), allocatable belongs_to_vtk_subdomain
Definition: vtk_mod.f:78
integer, parameter dimension_vtk
Definition: vtk_mod.f:81
real, dimension(:,:), allocatable coords_of_points
Definition: vtk_mod.f:44
logical, dimension(:), allocatable blocked_cell_at
Definition: cutcell_mod.f:361
subroutine write_cut_surface_vtk
Definition: vtk_out.f:1470
integer imin1
Definition: geometry_mod.f:40
double precision time
Definition: run_mod.f:45
subroutine unshift_dp_array(ARRAY_1, ARRAY_2)
character(len=18), dimension(dim_n_all) species_name
Definition: rxns_mod.f:44
character(len=32), dimension(dim_m, dim_n_s) species_alias_s
Definition: rxns_mod.f:52
double precision, dimension(:), allocatable e_turb_g
Definition: fldvar_mod.f:162
double precision, dimension(:), allocatable vol
Definition: geometry_mod.f:212
double precision, dimension(:), allocatable zg_t
Definition: cutcell_mod.f:46
double precision, dimension(:), allocatable axy_u
Definition: geometry_mod.f:222
subroutine write_vector_in_vtu_bin(VAR_NAME, VARX, VARY, VARZ, PASS)
Definition: vtk_out.f:1048
double precision, dimension(:), allocatable theta_u_te
Definition: cutcell_mod.f:216
double precision, dimension(dimension_vtk) vtk_x_w
Definition: vtk_mod.f:111
integer, dimension(0:15) ijk_of_node
Definition: cutcell_mod.f:78
logical, dimension(:), allocatable global_cut_cell_at
Definition: vtk_mod.f:58
integer kmin1
Definition: geometry_mod.f:44
double precision, dimension(:), allocatable noc_v_e
Definition: cutcell_mod.f:261
double precision, dimension(:), allocatable theta_ve
Definition: cutcell_mod.f:257
double precision, dimension(:), allocatable axz_v
Definition: geometry_mod.f:229
logical, dimension(dimension_vtk) vtk_ijk
Definition: vtk_mod.f:232
double precision, dimension(:), allocatable x
Definition: geometry_mod.f:129
double precision, parameter zero
Definition: param1_mod.f:27
double precision, dimension(0:15) x_node
Definition: cutcell_mod.f:74
integer jend1
Definition: compar_mod.f:80
integer jstart1
Definition: compar_mod.f:80
double precision, dimension(:), allocatable alpha_vn_c
Definition: cutcell_mod.f:263
integer, dimension(dimension_vtk) vtk_nzs
Definition: vtk_mod.f:135
integer, dimension(:,:), allocatable connectivity
Definition: cutcell_mod.f:111
double precision, dimension(:), allocatable vol_v
Definition: geometry_mod.f:233
integer dimension_j
Definition: param_mod.f:8
integer unit_cut_cell_log
Definition: cutcell_mod.f:8
character(len=255) vtu_frame_filename
Definition: vtk_mod.f:29
logical, dimension(:), allocatable snap
Definition: cutcell_mod.f:413
integer vtu_unit
Definition: vtk_mod.f:35