module_io.f90
Go to the documentation of this file.
1 !*************************************************************************
2 ! This file is part of HERCULES.
3 !
4 ! HERCULES is free software; you can redistribute it and/or modify
5 ! it under the terms of the GNU General Public License as published by
6 ! the Free Software Foundation; either version 3 of the License, or
7 ! (at your option) any later version.
8 !
9 ! HERCULES is distributed in the hope that it will be useful,
10 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
11 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 ! GNU General Public License for more details.
13 !
14 ! You should have received a copy of the GNU General Public License
15 ! along with this program. If not, see <http://www.gnu.org/licenses/>.
16 !
17 ! Copyright 2015, Ping He, MEAS, NCSU
18 !*************************************************************************
19 
20 
21 
22 module module_io
23 
25  implicit none
26 
27 contains
28 
29 
30 
31  ! this function prints cpu time, CFL number etc to screen
32  subroutine screen_cpu_time(u,v,w)
33 
34  use mpi
35 
36  implicit none
37 
38  real(mytype),dimension(:,:,:),intent(in) :: u, v, w
39 
40  real(mytype) :: cpu_t,cfl
41  integer :: i,j,k
42 
43  ! calculate cfl number
44  cfl=0.d0
45  do k=kstr3,kend3,1
46  do j=jstr3,jend3,1
47  do i=istr3,iend3,1
48 
49  if ( ( dt*abs(u(i,j,k))/dx + &
50  dt*abs(v(i,j,k))/dy + &
51  dt*abs(w(i,j,k))/dz(k) ) .gt. cfl) then
52 
53  cfl=dt*abs(u(i,j,k))/dx + &
54  dt*abs(v(i,j,k))/dy + &
55  dt*abs(w(i,j,k))/dz(k)
56 
57  endif
58 
59  enddo
60  enddo
61  enddo
62 
63  ! exchange CFL between processors
64  ! and find the max CFL number in the domain
65  call mpi_reduce(cfl,cfl_max,1,real_type,mpi_max,0,mpi_comm_world,ierr)
66 
67  ! calculate the simulation time
69 
70  ! print the CFL max to screen
71  if (myid .eq. 0) then
72  call cpu_time(cpu_t)
73  print *, " "
74  write (*,184) time_sim,cpu_t,cfl_max
75  184 format(' Time = ',es12.5,'; CPU Time = ',es11.4,'; CFL max = ',f5.2)
76  endif
77 
78  return
79 
80  end subroutine
81 
82 
83 
84  ! this function calculates the divergence of the flow fields and prints to screen
85  subroutine screen_div_error(u,v,w,isfinal,ishistout)
86 
87  use mpi
88  use module_tools
89 
90  implicit none
91 
92  real(mytype),dimension(:,:,:),intent(in) :: u,v,w
93  integer,intent(in) :: isfinal,ishistout
94 
95  integer :: i,j,k
96  real(mytype) :: dudx,dvdy,dwdz
97  real(mytype) :: div_avg, div_max
98  real(mytype) :: div_avg_all,div_max_all
99  real(mytype) :: div_hist_max_all
100 
101  ! calculate the divergence
102  div_avg=0.d0
103  div_max=0.d0
104  do k=kstr3,kend3,1
105  do j=jstr3,jend3,1
106  do i=istr3,iend3,1
107 
108  ! check to use 2nd or 4th schemes
109  if ( cds .eq. 1) then
110 
111  dudx = ( u(i,j,k)-u(i-1,j,k) ) / dx
112  dvdy = ( v(i,j,k)-v(i,j-1,k) ) / dy
113 
114  elseif ( cds .eq. 2) then
115 
116  dudx = ( -u(i+1,j,k) + &
117  27.d0*u(i,j,k) - &
118  27.d0*u(i-1,j,k) + &
119  u(i-2,j,k) ) / 24.d0 / dx
120 
121  dvdy = ( -v(i,j+1,k) + &
122  27.d0*v(i,j,k) - &
123  27.d0*v(i,j-1,k) + &
124  v(i,j-2,k) ) / 24.d0 / dy
125 
126  endif
127 
128  dwdz = ( w(i,j,k)-w(i,j,k-1) ) / dz(k)
129  div_avg=div_avg+abs(dudx+dvdy+dwdz)
130 
131  if ( abs(dudx+dvdy+dwdz) .gt. div_max ) then
132  div_max=(dudx+dvdy+dwdz)
133  endif
134 
135  if ( u(i,j,k) .ne. u(i,j,k) .or. &
136  v(i,j,k) .ne. v(i,j,k) .or. &
137  w(i,j,k) .ne. w(i,j,k) ) then
138  print *, 'NaN error! ',i,' ',j,' ',k
139  stop
140  endif
141 
142  if (isfinal .eq. 1) then
143  if ( abs(dudx+dvdy+dwdz) .gt. div_hist_max ) then
144  div_hist_max=(dudx+dvdy+dwdz)
145  endif
146  endif
147 
148  enddo
149  enddo
150  enddo
151  div_avg=div_avg/nx/ny/nz
152 
153  ! exchange div_avg div_max between processors
154  ! sum up all the div_avg between all processors
155  call mpi_reduce(div_avg,div_avg_all,1,real_type,mpi_sum,0, &
156  mpi_comm_world,ierr)
157  ! find the max div_max between all processors
158  call mpi_reduce(div_max,div_max_all,1,real_type,mpi_max,0, &
159  mpi_comm_world,ierr)
160  ! find the max div_hist_max between all processors
161  call mpi_reduce(div_hist_max,div_hist_max_all,1,real_type,mpi_max,0, &
162  mpi_comm_world,ierr)
163 
164  ! print div to screen
165  if (myid .eq. 0) then
166 
167  if ( isfinal .eq. 0) then
168 
169  write (*,104) div_avg_all,div_max_all
170  104 format(' DIV-ERR INIT : DIV_AVG = ',es11.4,'; DIV_MAX = ',es11.4)
171 
172  elseif ( isfinal .eq. 1 .and. ishistout .eq. 0) then
173 
174  write (*,105) div_avg_all,div_max_all
175  105 format(' DIV-ERR FINAL: DIV_AVG = ',es11.4,'; DIV_MAX = ',es11.4)
176 
177  elseif ( isfinal .eq. 1 .and. ishistout .eq. 1) then
178 
179  write (*,115) div_avg_all,div_max_all
180  115 format(' DIV-ERR FINAL: DIV_AVG = ',es11.4,'; DIV_MAX = ',es11.4)
181  write (*,106) div_hist_max_all
182  106 format(' DIV-ERR HIST_MAX = ', es11.4)
183 
184  ! print ri if is_ri_var=1
185  if (is_ri_var .eq. 1) then
186  write (*,258) got
187  258 format(' Ri_b = ', es12.5)
188  endif
189 
190  ! print bulk velocity if dp_opt=2
191  if (dp_opt .eq. 2) then
192  write (*,259) ubar+u_mrf
193  259 format(' Ubar = ', es12.5)
194  endif
195 
196  else
197  print *,'div out flag error!'
198  stop
199  endif
200 
201  endif
202 
203 
204  return
205 
206  end subroutine
207 
208 
209 
210  ! output time-series
211  subroutine output_time_series(u,v,w,p,t)
213  use mpi
214 
215  implicit none
216 
217  real(mytype),dimension(:,:,:),intent(in) :: u,v,w,t,p
218 
219  integer :: k,itmsr,jtmsr
220 
221  itmsr=ghst+1
222  jtmsr=ghst+1
223 
224  open(99,file='results/tmsr_u.dat',status='old', &
225  position='append',action='write')
226  do k=kstr3,kend3,1
227 
228  if (k .eq. kstr3) then
229  ! note that here we need to add the speed of moving reference frame
230  write(99,'(2f14.8)',advance='no') time_sim,u(itmsr,jtmsr,k)+u_mrf
231  elseif (k .eq. kend3) then
232  write(99,'(1f14.8)') u(itmsr,jtmsr,k)+u_mrf
233  else
234  write(99,'(1f14.8)',advance='no') u(itmsr,jtmsr,k)+u_mrf
235  endif
236 
237  enddo
238  close(99)
239 
240  open(98,file='results/tmsr_v.dat',status='old', &
241  position='append',action='write')
242  do k=kstr3,kend3,1
243 
244  if (k .eq. kstr3) then
245  write(98,'(2f14.8)',advance='no') time_sim,v(itmsr,jtmsr,k)
246  elseif (k .eq. kend3) then
247  write(98,'(1f14.8)') v(itmsr,jtmsr,k)
248  else
249  write(98,'(1f14.8)',advance='no') v(itmsr,jtmsr,k)
250  endif
251 
252  enddo
253  close(98)
254 
255 
256  open(97,file='results/tmsr_w.dat',status='old', &
257  position='append',action='write')
258  do k=kstr3,kend3,1
259 
260  if (k .eq. kstr3) then
261  write(97,'(2f14.8)',advance='no') time_sim,w(itmsr,jtmsr,k)
262  elseif (k .eq. kend3) then
263  write(97,'(1f14.8)') w(itmsr,jtmsr,k)
264  else
265  write(97,'(1f14.8)',advance='no') w(itmsr,jtmsr,k)
266  endif
267 
268  enddo
269  close(97)
270 
271 
272  open(96,file='results/tmsr_t.dat',status='old', &
273  position='append',action='write')
274  do k=kstr3,kend3,1
275 
276  if (k .eq. kstr3) then
277  write(96,'(2f14.8)',advance='no') time_sim,t(itmsr,jtmsr,k)
278  elseif (k .eq. kend3) then
279  write(96,'(1f14.8)') t(itmsr,jtmsr,k)
280  else
281  write(96,'(1f14.8)',advance='no') t(itmsr,jtmsr,k)
282  endif
283 
284  enddo
285  close(96)
286 
287  open(95,file='results/tmsr_p.dat',status='old', &
288  position='append',action='write')
289  do k=kstr3,kend3,1
290 
291  if (k .eq. kstr3) then
292  write(95,'(2f14.8)',advance='no') time_sim,p(itmsr,jtmsr,k)
293  elseif (k .eq. kend3) then
294  write(95,'(1f14.8)') p(itmsr,jtmsr,k)
295  else
296  write(95,'(1f14.8)',advance='no') p(itmsr,jtmsr,k)
297  endif
298 
299  enddo
300  close(95)
301 
302 
303  return
304 
305  end subroutine
306 
307 
308 
309  ! note that we store all the backup fields in staggered grids
310  ! instead of cell center
311  subroutine output_backup(u,v,w,p,t,ff,gg,hh,ss)
313  use decomp_2d_io
314 
315  implicit none
316 
317  real(mytype),dimension(:,:,:),intent(in) :: u,v,w,p,t,ff,gg,hh,ss
318 
319  if (myid .eq. 0) then
320  print *, 'Writing Backup Files...'
321  endif
322 
323  ! note that here we need to add the speed of moving reference frame
324  call decomp_2d_write_one(3,u(istr3:iend3,jstr3:jend3,kstr3:kend3)+u_mrf, &
325  'results/bak_u.dat')
326 
327  call decomp_2d_write_one(3,v(istr3:iend3,jstr3:jend3,kstr3:kend3), &
328  'results/bak_v.dat')
329 
330  call decomp_2d_write_one(3,w(istr3:iend3,jstr3:jend3,kstr3:kend3), &
331  'results/bak_w.dat')
332 
333  call decomp_2d_write_one(3,p(istr3:iend3,jstr3:jend3,kstr3:kend3), &
334  'results/bak_p.dat')
335 
336  call decomp_2d_write_one(3,t(istr3:iend3,jstr3:jend3,kstr3:kend3), &
337  'results/bak_t.dat')
338 
339  ! for RK3, we don't need to write f,g,h,s
340  ! note here the ff gg hh ss are based on the relative velocity (u-u_mrf)
341  ! if running a restart simulation using AB2 or AB3, one should not change u_mrf
342  if (dts .ne. 2) then
343 
344  call decomp_2d_write_one(3,ff(istr3:iend3,jstr3:jend3,kstr3:kend3), &
345  'results/bak_ff.dat')
346 
347  call decomp_2d_write_one(3,gg(istr3:iend3,jstr3:jend3,kstr3:kend3), &
348  'results/bak_gg.dat')
349 
350  call decomp_2d_write_one(3,hh(istr3:iend3,jstr3:jend3,kstr3:kend3), &
351  'results/bak_hh.dat')
352 
353  call decomp_2d_write_one(3,ss(istr3:iend3,jstr3:jend3,kstr3:kend3), &
354  'results/bak_ss.dat')
355 
356  endif
357 
358  return
359 
360  end subroutine
361 
362 
363 
364  ! note that we store all the instantaneous fields in staggered grids
365  ! instead of cell center
366  subroutine output_inst_fields(u,v,w,p,t,time_step)
368  use decomp_2d_io
369 
370  implicit none
371 
372  real(mytype),dimension(:,:,:),intent(in) :: u, v, w, p, t
373  integer, intent(in) :: time_step
374 
375  character(len=8) :: fmt1
376  character(len=6) :: xx1
377 
378  fmt1= '(i3.3)'
379  write (xx1,fmt1) time_step/iinstfl
380 
381  if (myid .eq. 0) then
382  print *, 'Writing Instantaneous Fields...'
383  endif
384 
385  ! note that here we need to add the speed of moving reference frame
386  call decomp_2d_write_one(3,u(istr3:iend3,jstr3:jend3,kstr3:kend3)+u_mrf, &
387  'results/inst_u'//trim(xx1)//'.dat')
388 
389  call decomp_2d_write_one(3,v(istr3:iend3,jstr3:jend3,kstr3:kend3), &
390  'results/inst_v'//trim(xx1)//'.dat')
391 
392  call decomp_2d_write_one(3,w(istr3:iend3,jstr3:jend3,kstr3:kend3), &
393  'results/inst_w'//trim(xx1)//'.dat')
394 
395  call decomp_2d_write_one(3,p(istr3:iend3,jstr3:jend3,kstr3:kend3), &
396  'results/inst_p'//trim(xx1)//'.dat')
397 
398  call decomp_2d_write_one(3,t(istr3:iend3,jstr3:jend3,kstr3:kend3), &
399  'results/inst_t'//trim(xx1)//'.dat')
400 
401  return
402 
403  end subroutine
404 
405 
406 
407  ! write 2d slice
408  subroutine output_2d_slices(u,v,w,p,t,time_step)
410  use decomp_2d_io
411 
412  implicit none
413 
414  real(mytype),dimension(:,:,:),intent(in) :: u, v, w, p, t
415  integer, intent(in) :: time_step
416 
417  character(len=8) :: fmt0
418  character(len=8) :: fmt1
419  character(len=6) :: ff1
420  character(len=6) :: ii1
421 
422  integer :: size_xy=0,size_xz=0,size_yz=0,n,test_size,indx
423 
424  fmt0= '(i6.6)'
425  fmt1= '(i4.4)'
426  write (ff1,fmt0) time_step/intv_2d
427 
428  ! output xy slices
429  if (isxy2d .eq. 1) then
430 
431  if (myid .eq. 0) then
432  print *, 'Writing X-Y Slices...'
433  endif
434 
435  ! we need to check how many slices to extract
436  test_size=size(xy2d_id)
437  indx=1
438  do n=1,test_size,1
439  if (xy2d_id(n) .eq. 0) then
440  size_xy=indx-1
441  exit
442  endif
443  indx=indx+1
444  enddo
445 
446  do n=1,size_xy,1
447 
448  write (ii1,fmt1) xy2d_id(n)
449  ! note that here we need to add the speed of moving reference frame
450  call decomp_2d_write_plane(3,u(istr3:iend3,jstr3:jend3,kstr3:kend3)+u_mrf, &
451  3,xy2d_id(n),'results/slices/xy2d_u_k'//trim(ii1)//'_t'//trim(ff1) )
452  call decomp_2d_write_plane(3,v(istr3:iend3,jstr3:jend3,kstr3:kend3), &
453  3,xy2d_id(n),'results/slices/xy2d_v_k'//trim(ii1)//'_t'//trim(ff1) )
454  call decomp_2d_write_plane(3,w(istr3:iend3,jstr3:jend3,kstr3:kend3), &
455  3,xy2d_id(n),'results/slices/xy2d_w_k'//trim(ii1)//'_t'//trim(ff1) )
456  call decomp_2d_write_plane(3,p(istr3:iend3,jstr3:jend3,kstr3:kend3), &
457  3,xy2d_id(n),'results/slices/xy2d_p_k'//trim(ii1)//'_t'//trim(ff1) )
458  call decomp_2d_write_plane(3,t(istr3:iend3,jstr3:jend3,kstr3:kend3), &
459  3,xy2d_id(n),'results/slices/xy2d_t_k'//trim(ii1)//'_t'//trim(ff1) )
460 
461  enddo
462 
463  endif
464 
465  ! output xz slices
466  if (isxz2d .eq. 1) then
467 
468  if (myid .eq. 0) then
469  print *, 'Writing X-Z Slices...'
470  endif
471 
472  ! we need to check how many slices to extract
473  test_size=size(xz2d_id)
474  indx=1
475  do n=1,test_size,1
476  if (xz2d_id(n) .eq. 0) then
477  size_xz=indx-1
478  exit
479  endif
480  indx=indx+1
481  enddo
482 
483  do n=1,size_xz,1
484 
485  write (ii1,fmt1) xz2d_id(n)
486  ! note that here we need to add the speed of moving reference frame
487  call decomp_2d_write_plane(3,u(istr3:iend3,jstr3:jend3,kstr3:kend3)+u_mrf, &
488  2,xz2d_id(n),'results/slices/xz2d_u_j'//trim(ii1)//'_t'//trim(ff1) )
489  call decomp_2d_write_plane(3,v(istr3:iend3,jstr3:jend3,kstr3:kend3), &
490  2,xz2d_id(n),'results/slices/xz2d_v_j'//trim(ii1)//'_t'//trim(ff1) )
491  call decomp_2d_write_plane(3,w(istr3:iend3,jstr3:jend3,kstr3:kend3), &
492  2,xz2d_id(n),'results/slices/xz2d_w_j'//trim(ii1)//'_t'//trim(ff1) )
493  call decomp_2d_write_plane(3,p(istr3:iend3,jstr3:jend3,kstr3:kend3), &
494  2,xz2d_id(n),'results/slices/xz2d_p_j'//trim(ii1)//'_t'//trim(ff1) )
495  call decomp_2d_write_plane(3,t(istr3:iend3,jstr3:jend3,kstr3:kend3), &
496  2,xz2d_id(n),'results/slices/xz2d_t_j'//trim(ii1)//'_t'//trim(ff1) )
497 
498  enddo
499 
500  endif
501 
502  ! output yz slices
503  if (isyz2d .eq. 1) then
504 
505  if (myid .eq. 0) then
506  print *, 'Writing Y-Z Slices...'
507  endif
508 
509  ! we need to check how many slices to extract
510  test_size=size(yz2d_id)
511  indx=1
512  do n=1,test_size,1
513  if (yz2d_id(n) .eq. 0) then
514  size_yz=indx-1
515  exit
516  endif
517  indx=indx+1
518  enddo
519 
520  do n=1,size_yz,1
521 
522  write (ii1,fmt1) yz2d_id(n)
523  ! note that here we need to add the speed of moving reference frame
524  call decomp_2d_write_plane(3,u(istr3:iend3,jstr3:jend3,kstr3:kend3)+u_mrf, &
525  1,yz2d_id(n),'results/slices/yz2d_u_i'//trim(ii1)//'_t'//trim(ff1) )
526  call decomp_2d_write_plane(3,v(istr3:iend3,jstr3:jend3,kstr3:kend3), &
527  1,yz2d_id(n),'results/slices/yz2d_v_i'//trim(ii1)//'_t'//trim(ff1) )
528  call decomp_2d_write_plane(3,w(istr3:iend3,jstr3:jend3,kstr3:kend3), &
529  1,yz2d_id(n),'results/slices/yz2d_w_i'//trim(ii1)//'_t'//trim(ff1) )
530  call decomp_2d_write_plane(3,p(istr3:iend3,jstr3:jend3,kstr3:kend3), &
531  1,yz2d_id(n),'results/slices/yz2d_p_i'//trim(ii1)//'_t'//trim(ff1) )
532  call decomp_2d_write_plane(3,t(istr3:iend3,jstr3:jend3,kstr3:kend3), &
533  1,yz2d_id(n),'results/slices/yz2d_t_i'//trim(ii1)//'_t'//trim(ff1) )
534 
535  enddo
536 
537  endif
538 
539 
540  return
541 
542  end subroutine
543 
544 
545 
546  subroutine output_mean_fields(mean_u,mean_v,mean_w,mean_t,mean_uu, &
547  mean_vv,mean_ww,mean_uw,mean_vw,mean_tt,mean_tw,time_step)
549  use decomp_2d_io
550 
551  implicit none
552 
553  real(mytype),dimension(:,:,:),intent(in) :: &
554  mean_u,mean_v,mean_w,mean_t,mean_uu,mean_vv,mean_ww,mean_uw,mean_vw, &
555  mean_tt,mean_tw
556  integer, intent(in) :: time_step
557 
558  if (myid .eq. 0) then
559  print *, 'Writing Mean Fields...'
560  endif
561 
562  call decomp_2d_write_one(3,mean_u(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
563  'results/mean_u.dat')
564 
565  call decomp_2d_write_one(3,mean_v(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
566  'results/mean_v.dat')
567 
568  call decomp_2d_write_one(3,mean_w(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
569  'results/mean_w.dat')
570 
571  call decomp_2d_write_one(3,mean_t(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
572  'results/mean_t.dat')
573 
574  call decomp_2d_write_one(3,mean_uu(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
575  'results/mean_uu.dat')
576 
577  call decomp_2d_write_one(3,mean_vv(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
578  'results/mean_vv.dat')
579 
580  call decomp_2d_write_one(3,mean_ww(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
581  'results/mean_ww.dat')
582 
583  call decomp_2d_write_one(3,mean_uw(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
584  'results/mean_uw.dat')
585 
586  call decomp_2d_write_one(3,mean_vw(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
587  'results/mean_vw.dat')
588 
589  call decomp_2d_write_one(3,mean_tt(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
590  'results/mean_tt.dat')
591 
592  call decomp_2d_write_one(3,mean_tw(istr3:iend3,jstr3:jend3,kstr3:kend3)/time_step, &
593  'results/mean_tw.dat')
594 
595 
596  return
597 
598  end subroutine
599 
600 
601 
602 end module
603 
604 
605 
subroutine output_time_series(u, v, w, p, t)
Definition: module_io.f90:212
integer, save, protected is_ri_var
integer, save, protected myid
integer, save, protected isxy2d
real(mytype), save, protected ubar
integer, save, protected ny
integer, save, protected isyz2d
integer, dimension(50), save, protected yz2d_id
subroutine screen_div_error(u, v, w, isfinal, ishistout)
Definition: module_io.f90:86
real(mytype), save, protected dt
integer, save, protected dts
integer, save, protected jend3
integer, save, protected nz
integer, save, protected intv_2d
subroutine output_2d_slices(u, v, w, p, t, time_step)
Definition: module_io.f90:409
integer, save, protected iinstfl
real(mytype), dimension(:), allocatable, save, protected dz
integer, dimension(50), save, protected xy2d_id
integer, save, protected istr3
integer, save, protected isxz2d
real(mytype), save, protected dx
integer, save, protected jstr3
subroutine output_backup(u, v, w, p, t, ff, gg, hh, ss)
Definition: module_io.f90:312
real(mytype), save, protected u_mrf
integer, parameter ghst
integer, save, protected cds
real(mytype), save time_sim
real(mytype), save, protected dy
integer, save, protected nx
real(mytype), save, protected got
subroutine output_inst_fields(u, v, w, p, t, time_step)
Definition: module_io.f90:367
real(mytype), save div_hist_max
integer, save, protected iend3
subroutine output_mean_fields(mean_u, mean_v, mean_w, mean_t, mean_uu, mean_vv, mean_ww, mean_uw, mean_vw, mean_tt, mean_tw, time_step)
Definition: module_io.f90:548
subroutine screen_cpu_time(u, v, w)
Definition: module_io.f90:33
integer, dimension(50), save, protected xz2d_id
integer, save, protected dp_opt
real(mytype), save cfl_max
integer, save, protected kstr3
integer, save, protected kend3