KTH framework for Nek5000 toolboxes; testing version  0.0.1
prepost.f
Go to the documentation of this file.
1 c-----------------------------------------------------------------------
2  subroutine set_outfld
3 
4 c Check if we are going to checkpoint at this timestep
5 c and set ifoutfld accordingly
6 
7  include 'SIZE'
8  include 'TOTAL'
9  include 'CTIMER'
10 
11  common /rdump/ ntdump
12 
13  ifoutfld = .false.
14 
15  if (iostep.le.0 .and. timeio.le.0) return
16 
17 c if (istep.ge.nsteps) lastep=1
18 
19  if (iostep.gt.0) then
20  if(mod(istep,iostep).eq.0) ifoutfld=.true.
21  else if (timeioe.ne.0.0) then
22  if (dnekclock_sync()-etimes .ge. (ntdump + 1)*timeio) then
23  ntdump=ntdump+1
24  ifoutfld=.true.
25  endif
26  else if (timeio.ne.0.0) then
27  if (time.ge.(ntdump + 1)*timeio) then
28  ntdump=ntdump+1
29  ifoutfld=.true.
30  endif
31  endif
32 
33  if (ioinfodmp.ne.0 .or. lastep.eq.1) ifoutfld=.true.
34 
35  return
36  end
37 c-----------------------------------------------------------------------
38  subroutine check_ioinfo
39 
40 c Check for io request in file 'ioinfo'
41 
42  include 'SIZE'
43  include 'TSTEP'
44  include 'INPUT'
45 
46  parameter(lxyz=lx1*ly1*lz1)
47  parameter(lpsc9=ldimt1+9)
48  common /ctmp1/ tdump(lxyz,lpsc9)
49  real*4 tdump
50  real tdmp(4)
51  equivalence(tdump,tdmp)
52 
53  integer maxstep
54  save maxstep
55  data maxstep /999999999/
56 
57  character*132 fname
58  character*1 fname1(132)
59  equivalence(fname,fname1)
60 
61  ioinfodmp=0
62  if (nid.eq.0 .and. (mod(istep,10).eq.0 .or. istep.lt.200)) then
63  call blank(fname1,size(fname1))
64  len = ltrunc(path,132)
65  call chcopy(fname1,path,len)
66  call chcopy(fname1(len+1),'ioinfo',6)
67  open(unit=87,file=fname,status='old',err=88)
68  read(87,*,end=87,err=87) idummy
69  if (ioinfodmp.eq.0) ioinfodmp=idummy
70  if (idummy.ne.0) then ! overwrite last i/o request
71  rewind(87)
72  write(87,86)
73  86 format(' 0')
74  endif
75  87 continue
76  close(unit=87)
77  88 continue
78  if (ioinfodmp.ne.0) write(6,*) 'Output:',ioinfodmp
79  endif
80 
81  tdmp(1)=ioinfodmp
82  call gop(tdmp,tdmp(3),'+ ',1)
83  ioinfodmp=tdmp(1)
84  if (ioinfodmp.lt.0) maxstep=abs(ioinfodmp)
85  if (istep.ge.maxstep.or.ioinfodmp.eq.-2) lastep=1
86 
87  return
88  end
89 c-----------------------------------------------------------------------
90  subroutine prepost(ifdoin,prefin)
91 
92 c Store results for later postprocessing
93 c
94 c Recent updates:
95 c
96 c p65 now indicates the number of parallel i/o files; iff p66 >= 6
97 c
98 c we now check whether we are going to checkpoint in set_outfld
99 c
100  include 'SIZE'
101  include 'TOTAL'
102  include 'CTIMER'
103 
104  common /scrcg/ pm1(lx1,ly1,lz1,lelv)
105 
106  character*3 prefin,prefix
107 
108  logical ifdoin
109 
110  if (ioinfodmp.eq.-2) return
111 
112 #ifdef TIMER
113  etime1=dnekclock_sync()
114 #endif
115 
116  prefix = prefin
117  if (prefix.eq.'his') prefix = ' '
118 
119  if (ifdoin) then
120  icalld=icalld+1
121  nprep=icalld
122 
123  call prepost_map(0) ! map pr and axisymm. arrays
124  call outfld(prefix)
125  call prepost_map(1) ! map back axisymm. arrays
126 
127 #ifdef TIMER
128  tprep=tprep+dnekclock_sync()-etime1
129 #endif
130  endif
131 
132  return
133  end
134 c-----------------------------------------------------------------------
135  subroutine prepost_map(isave) ! isave=0-->fwd, isave=1-->bkwd
136 
137 c Store results for later postprocessing
138 
139  include 'SIZE'
140  include 'TOTAL'
141 C
142 C Work arrays and temporary arrays
143 C
144  common /scruz/ vxax(lx1,ly1,lelv)
145  $ , vyax(lx1,ly1,lelv)
146  $ , prax(lx2,ly2,lelv)
147  $ , yax(lx1,ly1,lelt)
148  common /scrmg/ tax(lx1,ly1,lelt,ldimt)
149  common /scrcg/ pm1(lx1,ly1,lz1,lelv)
150 C
151 c
152  common /prepst/ pa(lx1,ly2,lz2),pb(lx1,ly1,lz2)
153  integer e
154 
155  if (isave.eq.0) then ! map to GLL grid
156 
157  if (ifaxis) then
158  ntotm1 = lx1*ly1*nelt
159  call copy (yax,ym1,ntotm1)
160  do 5 e=1,nelt
161  if (ifrzer(e)) then
162  call mxm (ym1(1,1,1,e),lx1,iatjl1,ly1,pb,ly1)
163  call copy (ym1(1,1,1,e),pb,lx1*ly1)
164  endif
165  5 continue
166  if (ifflow) then
167  ntotm1 = lx1*ly1*nelt
168  ntotm2 = lx2*ly2*nelt
169  call copy (vxax,vx,ntotm1)
170  call copy (vyax,vy,ntotm1)
171  call copy (prax,pr,ntotm2)
172  do 10 e=1,nelt
173  if (ifrzer(e)) then
174  call mxm (vx(1,1,1,e),lx1,iatjl1,ly1,pb,ly1)
175  call copy (vx(1,1,1,e),pb,lx1*ly1)
176  call mxm (vy(1,1,1,e),lx1,iatjl1,ly1,pb,ly1)
177  call copy (vy(1,1,1,e),pb,lx1*ly1)
178  call mxm (pr(1,1,1,e),lx2,iatjl2,ly2,pb,ly2)
179  call copy (pr(1,1,1,e),pb,lx2*ly2)
180  endif
181  10 continue
182  endif
183  if (ifheat) then
184  ntotm1 = lx1*ly1*nelt
185  do 15 ifldt=1,npscal+1
186  call copy (tax(1,1,1,ifldt),t(1,1,1,1,ifldt),ntotm1)
187  15 continue
188  do 30 e=1,nelt
189  if (ifrzer(e)) then
190  do 25 ifldt=1,npscal+1
191  call mxm (t(1,1,1,e,ifldt),lx1,iatjl1,ly1,
192  $ pb,ly1)
193  call copy (t(1,1,1,e,ifldt),pb,lx1*ly1)
194  25 continue
195  endif
196  30 continue
197  endif
198  endif
199 C Map the pressure onto the velocity mesh
200 C
201  ntott = lx1*ly1*lz1*nelt
202  ntot1 = lx1*ly1*lz1*nelt
203  nyz2 = ly2*lz2
204  nxy1 = lx1*ly1
205  nxyz = lx1*ly1*lz1
206  nxyz2 = lx2*ly2*lz2
207 C
208 
209  call rzero(pm1,ntott)
210  if (ifsplit) then
211  call copy(pm1,pr,ntot1)
212  elseif (if_full_pres) then
213  call rzero(pm1,ntot1)
214  do e=1,nelt
215  call copy(pm1(1,1,1,e),pr(1,1,1,e),nxyz2)
216  enddo
217  else
218  do 1000 e=1,nelt
219  call mxm (ixm21,lx1,pr(1,1,1,e),lx2,pa(1,1,1),nyz2)
220  do 100 iz=1,lz2
221  call mxm (pa(1,1,iz),lx1,iytm21,ly2,pb(1,1,iz),ly1)
222  100 continue
223  call mxm (pb(1,1,1),nxy1,iztm21,lz2,pm1(1,1,1,e),lz1)
224  1000 continue
225  endif
226 
227  else ! map back
228 
229  if (ifaxis) then
230  ntot1 = lx1*ly1*nelt
231  call copy (ym1,yax,ntot1)
232  if (ifflow) then
233  ntot1 = lx1*ly1*nelt
234  ntot2 = lx2*ly2*nelt
235  call copy (vx,vxax,ntot1)
236  call copy (vy,vyax,ntot1)
237  call copy (pr,prax,ntot2)
238  endif
239  if (ifheat) then
240  ntot1 = lx1*ly1*nelt
241  do 3000 ifldt=1,npscal+1
242  call copy (t(1,1,1,1,ifldt),tax(1,1,1,ifldt),ntot1)
243  3000 continue
244  endif
245  endif
246 
247  endif
248  return
249  end
250 c-----------------------------------------------------------------------
251  subroutine outfld(prefix)
252 
253 c output .fld file
254 
255  include 'SIZE'
256  include 'TOTAL'
257  include 'RESTART'
258 C
259 C Work arrays and temporary arrays
260 C
261  common /scrcg/ pm1(lx1,ly1,lz1,lelv)
262 c
263 c note, this usage of CTMP1 will be less than elsewhere if NELT ~> 3.
264  parameter(lxyz=lx1*ly1*lz1)
265  parameter(lpsc9=ldimt1+9)
266  common /ctmp1/ tdump(lxyz,lpsc9)
267  real*4 tdump
268  real tdmp(4)
269  equivalence(tdump,tdmp)
270 
271  real*4 test_pattern
272 
273  character*3 prefix
274  character*1 fhdfle1(132)
275  character*132 fhdfle
276  equivalence(fhdfle,fhdfle1)
277 
278  character*1 excode(30)
279  character*10 frmat
280 
281  common /nopenf/ nopen(99)
282 
283  common /rdump/ ntdump
284  data ndumps / 0 /
285 
286  logical ifxyo_s
287 
288  if(nio.eq.0) then
289  WRITE(6,1001) istep,time
290  1001 FORMAT(/,i9,1pe12.4,' Write checkpoint')
291  endif
292  call nekgsync()
293 
294  p66 = param(66)
295  if (abs(p66).eq.6) then
296  call mfo_outfld(prefix)
297  return
298  endif
299 
300  ifxyo_s = ifxyo ! Save ifxyo
301 
302  iprefix = i_find_prefix(prefix,99)
303 
304  ierr = 0
305  if (nid.eq.0) then
306 
307 c Open new file for each dump on /cfs
308  nopen(iprefix)=nopen(iprefix)+1
309 
310  if (prefix.eq.' '.and.nopen(iprefix).eq.1) ifxyo = .true. ! 1st file
311 
312  if (prefix.eq.'rst'.and.max_rst.gt.0)
313  $ nopen(iprefix) = mod1(nopen(iprefix),max_rst) ! restart
314 
315  call file2(nopen(iprefix),prefix)
316  if (p66.lt.1.0) then
317  open(unit=24,file=fldfle,form='formatted',status='unknown')
318  else
319  call byte_open (fldfle,ierr)
320 c write header as character string
321  call blank(fhdfle,132)
322  endif
323  endif
324  call bcast(ifxyo,lsize)
325  if(p66.ge.1.0)
326  $ call err_chk(ierr,'Error opening file in outfld. Abort. $')
327 
328 C Figure out what goes in EXCODE
329  CALL blank(excode,30)
330  ndumps=ndumps+1
331  i=1
332  if (mod(p66,1.0).eq.0.0) then !old header format
333  IF(ifxyo) then
334  excode(1)='X'
335  excode(2)=' '
336  excode(3)='Y'
337  excode(4)=' '
338  i = 5
339  IF(if3d) THEN
340  excode(i) ='Z'
341  excode(i+1)=' '
342  i = i + 2
343  ENDIF
344  ENDIF
345  IF(ifvo) then
346  excode(i) ='U'
347  excode(i+1)=' '
348  i = i + 2
349  ENDIF
350  IF(ifpo) THEN
351  excode(i)='P'
352  excode(i+1)=' '
353  i = i + 2
354  ENDIF
355  IF(ifto) THEN
356  excode(i)='T '
357  excode(i+1)=' '
358  i = i + 1
359  ENDIF
360  do iip=1,ldimt1
361  if (ifpsco(iip)) then
362  write(excode(iip+i) ,'(i1)') iip
363  write(excode(iip+i+1),'(a1)') ' '
364  i = i + 1
365  endif
366  enddo
367  else
368  !new header format
369  IF (ifxyo) THEN
370  excode(i)='X'
371  i = i + 1
372  ENDIF
373  IF (ifvo) THEN
374  excode(i)='U'
375  i = i + 1
376  ENDIF
377  IF (ifpo) THEN
378  excode(i)='P'
379  i = i + 1
380  ENDIF
381  IF (ifto) THEN
382  excode(i)='T'
383  i = i + 1
384  ENDIF
385  IF (ldimt.GT.1) THEN
386  npscalo = 0
387  do k = 1,ldimt-1
388  if(ifpsco(k)) npscalo = npscalo + 1
389  enddo
390  IF (npscalo.GT.0) THEN
391  excode(i) = 'S'
392  WRITE(excode(i+1),'(I1)') npscalo/10
393  WRITE(excode(i+2),'(I1)') npscalo-(npscalo/10)*10
394  ENDIF
395  ENDIF
396  endif
397 
398 
399 C Dump header
400  ierr = 0
401  if (nid.eq.0) call dump_header(excode,p66,ierr)
402  call err_chk(ierr,'Error dumping header in outfld. Abort. $')
403 
404  call get_id(id)
405 
406  nxyz = lx1*ly1*lz1
407 
408  ierr = 0
409  do ieg=1,nelgt
410 
411  jnid = gllnid(ieg)
412  ie = gllel(ieg)
413 
414  if (nid.eq.0) then
415  if (jnid.eq.0) then
416  call fill_tmp(tdump,id,ie)
417  else
418  mtype=2000+ie
419  len=4*id*nxyz
420  dum1=0.
421  call csend (mtype,dum1,wdsize,jnid,nullpid)
422  call crecv2 (mtype,tdump,len,jnid)
423  endif
424  if(ierr.eq.0) call out_tmp(id,p66,ierr)
425  elseif (nid.eq.jnid) then
426  call fill_tmp(tdump,id,ie)
427  dum1=0.
428  mtype=2000+ie
429  len=4*id*nxyz
430  call crecv2 (mtype,dum1,wdsize,node0)
431  call csend (mtype,tdump,len,node0,nullpid)
432  endif
433  enddo
434  call err_chk(ierr,'Error writing file in outfld. Abort. $')
435 
436  ifxyo = ifxyo_s ! restore ifxyo
437 
438  if (nid.eq.0) call close_fld(p66,ierr)
439  call err_chk(ierr,'Error closing file in outfld. Abort. $')
440 
441  return
442  end
443 c-----------------------------------------------------------------------
444  subroutine file2(nopen,PREFIX)
445 C----------------------------------------------------------------------
446 C
447 C Defines machine specific input and output file names.
448 C
449 C----------------------------------------------------------------------
450  include 'SIZE'
451  include 'INPUT'
452  include 'TSTEP'
453  include 'PARALLEL'
454 C
455  CHARACTER*132 NAME
456  CHARACTER*1 SESS1(132),PATH1(132),NAM1(132)
457  equivalence(session,sess1)
458  equivalence(path,path1)
459  equivalence(name,nam1)
460  CHARACTER*1 DMP(4),FLD(4),REA(4),HIS(4),SCH(4) ,ORE(4), NRE(4)
461  CHARACTER*4 DMP4 ,FLD4 ,REA4 ,HIS4 ,SCH4 ,ORE4 , NRE4
462  equivalence(dmp,dmp4), (fld,fld4), (rea,rea4), (his,his4)
463  $ , (sch,sch4), (ore,ore4), (nre,nre4)
464  CHARACTER*1 NUMRL(0:9)
465  DATA dmp4,fld4,rea4 /'.dmp','.fld','.rea'/
466  DATA his4,sch4 /'.his','.sch'/
467  DATA ore4,nre4 /'.ore','.nre'/
468  DATA numrl /'0','1','2','3','4','5','6','7','8','9'/
469  CHARACTER*78 STRING
470 c
471  character*1 prefix(3)
472 C
473  call blank(name ,132)
474  call blank(fldfle,132)
475 C
476  ls=ltrunc(session,132)
477  lpp=ltrunc(path,132)
478  lsp=ls+lpp
479  l = 0
480 
481 c Construct file names containing full path to host:
482 c DO 100 I=1,LPP
483 c l = l+1
484 c NAM1(l)=PATH1(I)
485 c 100 CONTINUE
486 C
487  if (prefix(1).ne.' '.and.prefix(2).ne.' '.and.
488  $ prefix(3).ne.' ') then
489  do i=1,3
490  l = l+1
491  nam1(l)=prefix(i)
492  enddo
493  endif
494 C
495  DO 200 i=1,ls
496  l = l+1
497  nam1(l)=sess1(i)
498  200 CONTINUE
499 C
500 C .fld file
501  DO 300 i=1,4
502  l = l+1
503  nam1(l)=fld(i)
504  300 CONTINUE
505  if (nopen.lt.100) then
506 C less than 100 dumps....
507  iten=nopen/10
508  l = l+1
509  nam1(l)=numrl(iten)
510  ione=mod(nopen,10)
511  l = l+1
512  nam1(l)=numrl(ione)
513  elseif (nopen.lt.1000) then
514 C less than 1000 dumps....
515  ihun=nopen/100
516  l = l+1
517  nam1(l)=numrl(ihun)
518  iten=mod(nopen,100)/10
519  l = l+1
520  nam1(l)=numrl(iten)
521  ione=mod(nopen,10)
522  l = l+1
523  nam1(l)=numrl(ione)
524  elseif (nopen.lt.10000) then
525 C less than 10000 dumps....
526  itho=nopen/1000
527  l = l+1
528  nam1(l)=numrl(itho)
529  ihun=mod(nopen,1000)/100
530  l = l+1
531  nam1(l)=numrl(ihun)
532  iten=mod(nopen,100)/10
533  l = l+1
534  nam1(l)=numrl(iten)
535  ione=mod(nopen,10)
536  l = l+1
537  nam1(l)=numrl(ione)
538  endif
539  fldfle=name
540 C
541 C Write the name of the .fld file to the logfile.
542 C
543  if (nio.eq.0) then
544  call chcopy(string,fldfle,78)
545  write(6,1000) istep,time,string
546  1000 format(/,i9,1pe12.4,' OPEN: ',a78)
547  endif
548 
549  return
550  end
551 c=======================================================================
552  subroutine rzero4(a,n)
553  real*4 a(1)
554  DO 100 i = 1, n
555  100 a(i ) = 0.0
556  return
557  end
558 c=======================================================================
559  subroutine copyx4(a,b,n)
560  real*4 a(1)
561  REAL B(1)
562  DO 100 i = 1, n
563  100 a(i) = b(i)
564  return
565  end
566 c=======================================================================
567  subroutine copy4r(a,b,n)
568  real a(1)
569  real*4 b(1)
570  do i = 1, n
571  a(i) = b(i)
572  enddo
573  return
574  end
575 c=======================================================================
576  function i_find_prefix(prefix,imax)
577 c
578  character*3 prefix
579  character*3 prefixes(99)
580  save prefixes
581  data prefixes /99*'...'/
582 c
583  integer nprefix
584  save nprefix
585  data nprefix /0/
586 c
587 c Scan existing list of prefixes for a match to "prefix"
588 c
589  do i=1,nprefix
590  if (prefix.eq.prefixes(i)) then
591  i_find_prefix = i
592  return
593  endif
594  enddo
595 c
596 c If we're here, we didn't find a match.. bump list and return
597 c
598  nprefix = nprefix + 1
599  prefixes(nprefix) = prefix
600  i_find_prefix = nprefix
601 c
602 c Array bounds check on prefix list
603 c
604  if (nprefix.gt.99.or.nprefix.gt.imax) then
605  write(6,*) 'Hey! nprefix too big! ABORT in i_find_prefix'
606  $ ,nprefix,imax
607  call exitt
608  endif
609 c
610  return
611  end
612 c-----------------------------------------------------------------------
613  subroutine dump_header(excodein,p66,ierr)
614 
615  include 'SIZE'
616  include 'TOTAL'
617 
618  character*30 excodein
619 
620  character*30 excode
621  character*1 excode1(30)
622  equivalence(excode,excode1)
623 
624  real*4 test_pattern
625 
626  character*1 fhdfle1(132)
627  character*132 fhdfle
628  equivalence(fhdfle,fhdfle1)
629 
630  write(excode,'(A30)') excodein
631 
632  ikstep = istep
633  do ik=1,10
634  if (ikstep.gt.9999) ikstep = ikstep/10
635  enddo
636 
637  call blank(fhdfle,132)
638 
639 c write(6,111) ! print on screen
640 c $ nelgt,lx1,ly1,lz1,time,istep,excode
641 c
642  if (mod(p66,1.0).eq.0.0) then ! old header format
643  if (p66.lt.1.0) then !ASCII
644  if(nelgt.lt.10000) then
645  WRITE(24,'(4i4,1pe14.7,I5,1X,30A1,1X,A12)')
646  $ nelgt,lx1,ly1,lz1,time,ikstep,(excode1(i),i=1,30),
647  $ 'NELT,NX,NY,N'
648  else
649  WRITE(24,'(i10,3i4,1pe18.9,I9,1X,30A1,1X,A12)')
650  $ nelgt,lx1,ly1,lz1,time,ikstep,(excode1(i),i=1,30),
651  $ 'NELT,NX,NY,N'
652  endif
653  else !Binary
654  if (nelgt.lt.10000) then
655  WRITE(fhdfle,'(4I4,1pe14.7,I5,1X,30A1,1X,A12)')
656  $ nelgt,lx1,ly1,lz1,time,ikstep,(excode1(i),i=1,30),
657  $ ' 4 NELT,NX,NY,N'
658  else
659  write(fhdfle,'(i10,3i4,1P1e18.9,i9,1x,30a1)')
660  $ nelgt,lx1,ly1,lz1,time,istep,(excode1(i),i=1,30)
661  endif
662  call byte_write(fhdfle,20,ierr)
663  endif
664  else ! new header format
665  if (p66.eq.0.1) then
666  write(24,111)
667  $ nelgt,lx1,ly1,lz1,time,istep,excode
668  else
669  write(fhdfle,111)
670  $ nelgt,lx1,ly1,lz1,time,istep,excode
671  call byte_write(fhdfle,20,ierr)
672  endif
673  111 FORMAT(i10,1x,i2,1x,i2,1x,i2,1x,1p1e18.9,1x,i9,1x,a)
674  endif
675 
676  if(ierr.ne.0) return
677 
678  cdrror=0.0
679  if (p66.LT.1.0) then ! formatted i/o
680  WRITE(24,'(6G11.4)')(cdrror,i=1,nelgt) ! dummy
681  else
682 C write byte-ordering test pattern to byte file...
683  test_pattern = 6.54321
684  call byte_write(test_pattern,1,ierr)
685  endif
686 
687  return
688  end
689 c-----------------------------------------------------------------------
690  subroutine fill_tmp(tdump,id,ie)
691 C
692  include 'SIZE'
693  include 'TOTAL'
694 c
695  common /scrcg/ pm1(lx1,ly1,lz1,lelv)
696 C
697 C Fill work array
698 C
699  parameter(lxyz=lx1*ly1*lz1)
700  parameter(lpsc9=ldimt1+9)
701  real*4 tdump(lxyz,lpsc9)
702 C
703  nxyz = lx1*ly1*lz1
704 c
705  id=0
706  IF(ifxyo)then
707  id=id+1
708  CALL copyx4(tdump(1,id),xm1(1,1,1,ie),nxyz)
709  id=id+1
710  CALL copyx4(tdump(1,id),ym1(1,1,1,ie),nxyz)
711  IF(if3d) then
712  id=id+1
713  CALL copyx4(tdump(1,id),zm1(1,1,1,ie),nxyz)
714  ENDIF
715  ENDIF
716 c
717  IF(ifvo)then
718  IF (ie.LE.nelt) then
719  id=id+1
720  CALL copyx4(tdump(1,id),vx(1,1,1,ie),nxyz)
721  id=id+1
722  CALL copyx4(tdump(1,id),vy(1,1,1,ie),nxyz)
723  IF(if3d)then
724  id=id+1
725  CALL copyx4(tdump(1,id),vz(1,1,1,ie),nxyz)
726  ENDIF
727  ELSE
728  id=id+1
729  CALL rzero4(tdump(1,id),nxyz)
730  id=id+1
731  CALL rzero4(tdump(1,id),nxyz)
732  IF(if3d)then
733  id=id+1
734  CALL rzero4(tdump(1,id),nxyz)
735  ENDIF
736  ENDIF
737  ENDIF
738  IF(ifpo)then
739  IF (ie.LE.nelt) then
740  id=id+1
741  CALL copyx4(tdump(1,id),pm1(1,1,1,ie),nxyz)
742  ELSE
743  id=id+1
744  CALL rzero4(tdump(1,id),nxyz)
745  ENDIF
746  ENDIF
747  IF(ifto)then
748  id=id+1
749  CALL copyx4(tdump(1,id),t(1,1,1,ie,1),nxyz)
750  ENDIF
751 C PASSIVE SCALARS
752  do iip=1,ldimt1
753  if (ifpsco(iip)) then
754  id=id+1
755  call copyx4(tdump(1,id),t(1,1,1,ie,iip+1),nxyz)
756  endif
757  enddo
758 c
759  return
760  end
761 c-----------------------------------------------------------------------
762  subroutine get_id(id) ! Count amount of data to be shipped
763 
764  include 'SIZE'
765  include 'TOTAL'
766 
767  id=0
768 
769  if (ifxyo) id=id+ldim
770  if (ifvo) id=id+ldim
771  if (ifpo) id=id+1
772  if (ifto) id=id+1
773 
774  do iip=1,ldimt1
775  if (ifpsco(iip)) id=id+1 ! Passive scalars
776  enddo
777 
778  return
779  end
780 c-----------------------------------------------------------------------
781  subroutine close_fld(p66,ierr)
782 
783  include 'SIZE'
784  include 'TOTAL'
785 
786  if (nid.eq.0) then
787  if (p66.lt.1) then
788  close(unit=24)
789  else
790  call byte_close(ierr)
791  endif
792  endif
793 
794  return
795  end
796 c-----------------------------------------------------------------------
797  subroutine out_tmp(id,p66,ierr)
798 
799  include 'SIZE'
800  include 'TOTAL'
801 
802  parameter(lxyz=lx1*ly1*lz1)
803  parameter(lpsc9=ldimt1+9)
804 
805  common /ctmp1/ tdump(lxyz,lpsc9)
806  real*4 tdump
807 
808  character*11 frmat
809 
810  nxyz = lx1*ly1*lz1
811 
812  call blank(frmat,11)
813  if (id.le.9) then
814  WRITE(frmat,1801) id
815  1801 FORMAT('(1p',i1,'e14.6)')
816  else
817  WRITE(frmat,1802) id
818  1802 FORMAT('(1p',i2,'e14.6)')
819  endif
820 
821  if (p66.lt.1.0) then
822 C formatted i/o
823  WRITE(24,frmat)
824  $ ((tdump(i,ii),ii=1,id),i=1,nxyz)
825  else
826 C C binary i/o
827  do ii=1,id
828  call byte_write(tdump(1,ii),nxyz,ierr)
829  if(ierr.ne.0) goto 101
830  enddo
831  endif
832  101 continue
833 
834  return
835  end
836 c-----------------------------------------------------------------------
837  subroutine mfo_outfld(prefix) ! muti-file output
838 
839  include 'SIZE'
840  include 'TOTAL'
841  include 'RESTART'
842  common /scrcg/ pm1(lx1,ly1,lz1,lelv) ! mapped pressure
843 
844  integer*8 offs0,offs,nbyte,stride,strideB,nxyzo8
845  character*3 prefix
846  logical ifxyo_s
847 
848  common /scruz/ ur1(lxo*lxo*lxo*lelt)
849  & , ur2(lxo*lxo*lxo*lelt)
850  & , ur3(lxo*lxo*lxo*lelt)
851 
852  tiostart=dnekclock_sync()
853 
854  call io_init
855 
856  ifxyo_s = ifxyo
857  ifxyo_ = ifxyo
858  nout = nelt
859  nxo = lx1
860  nyo = ly1
861  nzo = lz1
862  if (ifreguo) then ! dump on regular (uniform) mesh
863  if (nrg.gt.lxo) then
864  if (nid.eq.0) write(6,*)
865  & 'WARNING: nrg too large, reset to lxo!'
866  nrg = lxo
867  endif
868  nxo = nrg
869  nyo = nrg
870  nzo = 1
871  if(if3d) nzo = nrg
872  endif
873  offs0 = iheadersize + 4 + isize*nelgt
874 
875  ierr=0
876  if (nid.eq.pid0) then
877  call mfo_open_files(prefix,ierr) ! open files on i/o nodes
878  endif
879  call err_chk(ierr,'Error opening file in mfo_open_files. $')
880  call bcast(ifxyo_,lsize)
881  ifxyo = ifxyo_
882  call mfo_write_hdr ! create element mapping +
883 
884 c call exitti('this is wdsizo A:$',wdsizo)
885  ! write hdr
886  nxyzo8 = nxo*nyo*nzo
887  strideb = nelb * nxyzo8*wdsizo
888  stride = nelgt* nxyzo8*wdsizo
889 
890  ioflds = 0
891  ! dump all fields based on the t-mesh to avoid different
892  ! topologies in the post-processor
893  if (ifxyo) then
894  offs = offs0 + ldim*strideb
895  call byte_set_view(offs,ifh_mbyte)
896  if (ifreguo) then
897  call map2reg(ur1,nrg,xm1,nout)
898  call map2reg(ur2,nrg,ym1,nout)
899  if (if3d) call map2reg(ur3,nrg,zm1,nout)
900  call mfo_outv(ur1,ur2,ur3,nout,nxo,nyo,nzo)
901  else
902  call mfo_outv(xm1,ym1,zm1,nout,nxo,nyo,nzo)
903  endif
904  ioflds = ioflds + ldim
905  endif
906  if (ifvo ) then
907  offs = offs0 + ioflds*stride + ldim*strideb
908  call byte_set_view(offs,ifh_mbyte)
909  if (ifreguo) then
910  call map2reg(ur1,nrg,vx,nout)
911  call map2reg(ur2,nrg,vy,nout)
912  if (if3d) call map2reg(ur3,nrg,vz,nout)
913  call mfo_outv(ur1,ur2,ur3,nout,nxo,nyo,nzo)
914  else
915  call mfo_outv(vx,vy,vz,nout,nxo,nyo,nzo) ! B-field handled thru outpost
916  endif
917  ioflds = ioflds + ldim
918  endif
919  if (ifpo ) then
920  offs = offs0 + ioflds*stride + strideb
921  call byte_set_view(offs,ifh_mbyte)
922  if (ifreguo) then
923  call map2reg(ur1,nrg,pm1,nout)
924  call mfo_outs(ur1,nout,nxo,nyo,nzo)
925  else
926  call mfo_outs(pm1,nout,nxo,nyo,nzo)
927  endif
928  ioflds = ioflds + 1
929  endif
930  if (ifto ) then
931  offs = offs0 + ioflds*stride + strideb
932  call byte_set_view(offs,ifh_mbyte)
933  if (ifreguo) then
934  call map2reg(ur1,nrg,t,nout)
935  call mfo_outs(ur1,nout,nxo,nyo,nzo)
936  else
937  call mfo_outs(t,nout,nxo,nyo,nzo)
938  endif
939  ioflds = ioflds + 1
940  endif
941  do k=1,ldimt-1
942  if(ifpsco(k)) then
943  offs = offs0 + ioflds*stride + strideb
944  call byte_set_view(offs,ifh_mbyte)
945  if (ifreguo) then
946  call map2reg(ur1,nrg,t(1,1,1,1,k+1),nout)
947  call mfo_outs(ur1,nout,nxo,nyo,nzo)
948  else
949  call mfo_outs(t(1,1,1,1,k+1),nout,nxo,nyo,nzo)
950  endif
951  ioflds = ioflds + 1
952  endif
953  enddo
954  dnbyte = 1.*ioflds*nout*wdsizo*nxo*nyo*nzo
955 
956  if (if3d) then
957  offs0 = offs0 + ioflds*stride
958  strideb = nelb *2*4 ! min/max single precision
959  stride = nelgt*2*4
960  ioflds = 0
961  ! add meta data to the end of the file
962  if (ifxyo) then
963  offs = offs0 + ldim*strideb
964  call byte_set_view(offs,ifh_mbyte)
965  call mfo_mdatav(xm1,ym1,zm1,nout)
966  ioflds = ioflds + ldim
967  endif
968  if (ifvo ) then
969  offs = offs0 + ioflds*stride + ldim*strideb
970  call byte_set_view(offs,ifh_mbyte)
971  call mfo_mdatav(vx,vy,vz,nout)
972  ioflds = ioflds + ldim
973  endif
974  if (ifpo ) then
975  offs = offs0 + ioflds*stride + strideb
976  call byte_set_view(offs,ifh_mbyte)
977  call mfo_mdatas(pm1,nout)
978  ioflds = ioflds + 1
979  endif
980  if (ifto ) then
981  offs = offs0 + ioflds*stride + strideb
982  call byte_set_view(offs,ifh_mbyte)
983  call mfo_mdatas(t,nout)
984  ioflds = ioflds + 1
985  endif
986  do k=1,ldimt-1
987  offs = offs0 + ioflds*stride + strideb
988  call byte_set_view(offs,ifh_mbyte)
989  if(ifpsco(k)) call mfo_mdatas(t(1,1,1,1,k+1),nout)
990  ioflds = ioflds + 1
991  enddo
992  dnbyte = dnbyte + 2.*ioflds*nout*wdsizo
993  endif
994 
995  ierr = 0
996  if (nid.eq.pid0) then
997  if(ifmpiio) then
998  call byte_close_mpi(ifh_mbyte,ierr)
999  else
1000  call byte_close(ierr)
1001  endif
1002  endif
1003  call err_chk(ierr,'Error closing file in mfo_outfld. Abort. $')
1004 
1005  tio = dnekclock_sync()-tiostart
1006  if (tio.le.0) tio=1.
1007 
1008  dnbyte = glsum(dnbyte,1)
1009  dnbyte = dnbyte + iheadersize + 4. + isize*nelgt
1010  dnbyte = dnbyte/1024/1024
1011  if(nio.eq.0) write(6,7) istep,time,dnbyte,dnbyte/tio,
1012  & nfileo
1013  7 format(/,i9,1pe12.4,' done :: Write checkpoint',/,
1014  & 30x,'file size = ',3pg12.2,'MB',/,
1015  & 30x,'avg data-throughput = ',0pf7.1,'MB/s',/,
1016  & 30x,'io-nodes = ',i5,/)
1017 
1018  ifxyo = ifxyo_s ! restore old value
1019 
1020  return
1021  end
1022 c-----------------------------------------------------------------------
1023  subroutine io_init ! determine which nodes will output
1024  character*132 hname
1025 
1026  include 'SIZE'
1027  include 'INPUT'
1028  include 'PARALLEL'
1029  include 'RESTART'
1030 
1031  ifdiro = .false.
1032 
1033  ifmpiio = .false.
1034  if(abs(param(65)).eq.1 .and. abs(param(66)).eq.6) ifmpiio=.true.
1035 #ifdef NOMPIIO
1036  ifmpiio = .false.
1037 #endif
1038 
1039  wdsizo = 4
1040  if (param(63).gt.0) wdsizo = 8 ! 64-bit .fld file
1041  nrg = lxo
1042 
1043  if(ifmpiio) then
1044  nfileo = np
1045  nproc_o = 1
1046  fid0 = 0
1047  pid0 = nid
1048  pid1 = 0
1049  else
1050  if(param(65).lt.0) ifdiro = .true. ! p65 < 0 --> multi subdirectories
1051  nfileo = abs(param(65))
1052  if(nfileo.eq.0) nfileo = 1
1053  if(np.lt.nfileo) nfileo=np
1054  nproc_o = np / nfileo ! # processors pointing to pid0
1055  fid0 = nid/nproc_o ! file id
1056  pid0 = nproc_o*fid0 ! my parent i/o node
1057  pid1 = min(np-1,pid0+nproc_o-1) ! range of sending procs
1058  endif
1059 
1060  ! how many elements are present up to rank nid
1061  nn = nelt
1062  nelb = igl_running_sum(nn)
1063  nelb = nelb - nelt
1064 
1065  pid00 = glmin(pid0,1)
1066 
1067  return
1068  end
1069 c-----------------------------------------------------------------------
1070  subroutine mfo_open_files(prefix,ierr) ! open files
1071 
1072  include 'SIZE'
1073  include 'INPUT'
1074  include 'PARALLEL'
1075  include 'RESTART'
1076 
1077  character*1 prefix(3)
1078  character*3 prefx
1079 
1080  character*132 fname
1081  character*1 fnam1(132)
1082  equivalence(fnam1,fname)
1083 
1084  character*6 six,str
1085  save six
1086  data six / "??????" /
1087 
1088 
1089  character*1 slash,dot
1090  save slash,dot
1091  data slash,dot / '/' , '.' /
1092 
1093  integer nopen(99,2)
1094  save nopen
1095  data nopen / 198*0 /
1096 
1097  call blank(fname,132) ! zero out for byte_open()
1098 
1099  iprefix = i_find_prefix(prefix,99)
1100  if (ifreguo) then
1101  nopen(iprefix,2) = nopen(iprefix,2)+1
1102  nfld = nopen(iprefix,2)
1103  else
1104  nopen(iprefix,1) = nopen(iprefix,1)+1
1105  nfld = nopen(iprefix,1)
1106  endif
1107 
1108  call chcopy(prefx,prefix,3) ! check for full-restart request
1109  if (prefx.eq.'rst'.and.max_rst.gt.0) nfld = mod1(nfld,max_rst)
1110 
1111  call restart_nfld( nfld, prefix ) ! Check for Restart option.
1112  if (prefx.eq.' '.and.nfld.eq.1) ifxyo_ = .true. ! 1st file
1113 
1114  if(ifmpiio) then
1115  rfileo = 1
1116  else
1117  rfileo = nfileo
1118  endif
1119  ndigit = log10(rfileo) + 1
1120 
1121  lenp = ltrunc(path,132)
1122  call chcopy(fnam1(1),path,lenp)
1123  k = 1 + lenp
1124 
1125  if (ifdiro) then ! Add directory
1126  call chcopy(fnam1(k),'A',1)
1127  k = k + 1
1128  call chcopy(fnam1(k),six,ndigit) ! put ???? in string
1129  k = k + ndigit
1130  call chcopy(fnam1(k),slash,1)
1131  k = k + 1
1132  endif
1133 
1134  if (prefix(1).ne.' '.and.prefix(2).ne.' '.and. ! Add prefix
1135  $ prefix(3).ne.' ') then
1136  call chcopy(fnam1(k),prefix,3)
1137  k = k + 3
1138  endif
1139 
1140  len=ltrunc(session,132) ! Add SESSION
1141  call chcopy(fnam1(k),session,len)
1142  k = k+len
1143 
1144  if (ifreguo) then
1145  len=4
1146  call chcopy(fnam1(k),'_reg',len)
1147  k = k+len
1148  endif
1149 
1150  call chcopy(fnam1(k),six,ndigit) ! Add file-id holder
1151  k = k + ndigit
1152 
1153  call chcopy(fnam1(k ),dot,1) ! Add .f appendix
1154  call chcopy(fnam1(k+1),'f',1)
1155  k = k + 2
1156 
1157  write(str,4) nfld ! Add nfld number
1158  4 format(i5.5)
1159  call chcopy(fnam1(k),str,5)
1160  k = k + 5
1161 
1162  call addfid(fname,fid0)
1163 
1164  if(ifmpiio) then
1165  if(nio.eq.0) write(6,*) ' FILE:',fname
1166  call byte_open_mpi(fname,ifh_mbyte,.false.,ierr)
1167  else
1168  if(nid.eq.pid0) write(6,*) ' FILE:',fname
1169  call byte_open(fname,ierr)
1170  endif
1171 
1172  return
1173  end
1174 c-----------------------------------------------------------------------
1175 
1176  subroutine restart_nfld( nfld, prefix )
1177  include 'SIZE' ! For nio
1178  character*3 prefix
1179 c
1180 c Check for Restart option and return proper nfld value.
1181 c Also, convenient spot to explain restart strategy.
1182 c
1183 c
1184 c The approach is as follows:
1185 c
1186 c Prefix rs4 would indicate 4 files in the restart cycle.
1187 c
1188 c This would be normal usage for velocity only, with
1189 c checkpoints taking place in synch with standard io.
1190 c
1191 c The resultant restart sequence might look like:
1192 c
1193 c blah.fld09 Step 0
1194 c rs4blah.fld01 1
1195 c rs4blah.fld02 2
1196 c
1197 c which implies that fld09 would be used as the i.c.
1198 c in the restart, rs4blah.fld01 would overwrite the
1199 c solution at Step 1, and rs4blah.fld02 would overwrite
1200 c Step 2. Net result is that Steps 0-2 of the restart
1201 c session have solutions identical to those computed in
1202 c the prior run. (It's important that both runs use
1203 c the same dt in this case.)
1204 c
1205 c
1206 c Another equally possible restart sequence would be:
1207 c
1208 c
1209 c blah.fld10 Step 0
1210 c rs4blah.fld03 1
1211 c rs4blah.fld04 2
1212 c
1213 c Why the 3 & 4 ? If one were to use only 1 & 2, there
1214 c is a risk that the system crashes while writing, say,
1215 c rs4blah.fld01, in which case the restart is compromised --
1216 c very frustrating at the end of a run that has been queued
1217 c for a week. By providing a toggled sequence in pairs such as
1218 c
1219 c (1,2), (3,4), (1,2), ...
1220 c
1221 c ensures that one always has at least one complete restart
1222 c sequence. In the example above, the following files would
1223 c be written, in order:
1224 c
1225 c :
1226 c :
1227 c blah.fld09
1228 c rs4blah.fld01
1229 c rs4blah.fld02
1230 c blah.fld10
1231 c rs4blah.fld03
1232 c rs4blah.fld04
1233 c blah.fld11
1234 c rs4blah.fld01 (overwriting existing rs4blah.fld01)
1235 c rs4blah.fld02 ( " " " .fld02)
1236 c blah.fld12
1237 c rs4blah.fld03 ( etc. )
1238 c rs4blah.fld04
1239 c :
1240 c :
1241 c
1242 c
1243 c Other strategies are possible, according to taste.
1244 c
1245 c Here is a data-intensive one:
1246 c
1247 c MHD + double-precision restart, but single-precision std files
1248 c
1249 c In this case, single-precision files are kept as the running
1250 c file sequence (i.e., for later post-processing) but dbl-prec.
1251 c is required for restart. A total of 12 temporary restart files
1252 c must be saved: (3 for velocity, 3 for B-field) x 2 for redundancy.
1253 c
1254 c This is expressed, using hexadecimal notation (123456789abc...),
1255 c as prefix='rsc'.
1256 c
1257 c
1258  character*16 kst
1259  save kst
1260  data kst / '0123456789abcdef' /
1261  character*1 ks1(0:15),kin
1262  equivalence(ks1,kst)
1263 
1264 c
1265 c
1266  if (indx1(prefix,'rs',2).eq.1) then
1267  read(prefix,3) kin
1268  3 format(2x,a1)
1269  do kfld=1,15
1270  if (ks1(kfld).eq.kin) goto 10
1271  enddo
1272  10 if (kfld.eq.16) kfld=4 ! std. default
1273  nfln = mod1(nfld,kfld) ! Restart A (1,2) and B (3,4)
1274  if (nio.eq.0) write(6,*) nfln,nfld,kfld,' kfld'
1275  nfld = nfln
1276  endif
1277 
1278  return
1279  end
1280 c-----------------------------------------------------------------------
1281  subroutine full_restart_save(iosave)
1282 
1283  integer iosave,save_size,nfld_save
1284  logical if_full_pres_tmp
1285 
1286  include 'SIZE'
1287  include 'INPUT'
1288  include 'TSTEP'
1289 
1290  if (param(27).lt. 0) then
1291  nfld_save=abs(param(27)) ! For full restart
1292  else
1293  nfld_save=3
1294  endif
1295  save_size=8 ! For full restart
1296 
1297  dtmp = param(63)
1298  if_full_pres_tmp = if_full_pres
1299 
1300  param(63) = 1 ! Enforce 64-bit output
1301  if_full_pres = .true. !Preserve mesh 2 pressure
1302 
1303  if (lastep.ne.1) call restart_save(iosave,nfld_save)
1304 
1305  param(63) = dtmp
1306  if_full_pres = if_full_pres_tmp
1307 
1308  return
1309  end
1310 c-----------------------------------------------------------------------
1311  subroutine restart_save(iosave,nfldi)
1312 
1313  integer iosave,nfldi
1314 
1315 
1316 c Save current fields for later restart.
1317 c
1318 c Input arguments:
1319 c
1320 c .iosave plays the usual triggering role, like iostep
1321 c
1322 c .nfldi is the number of rs files to save before overwriting
1323 c
1324 
1325  include 'SIZE'
1326  include 'TOTAL'
1327  include 'RESTART'
1328 
1329  character*3 prefix
1330 
1331  character*17 kst
1332  save kst
1333  data kst / '0123456789abcdefx' /
1334  character*1 ks1(0:16)
1335  equivalence(ks1,kst)
1336 
1337  logical if_full_pres_tmp
1338 
1339  iosav = iosave
1340 
1341  if (iosav.eq.0) iosav = iostep
1342  if (iosav.eq.0) return
1343 
1344  iotest = 0
1345 c if (iosav.eq.iostep) iotest = 1 ! currently spoiled because of
1346 c ! incompatible format of .fld
1347 c ! and multi-file i/o; the latter
1348 c ! is the only form used for restart
1349 
1350  nfld = nfldi*2
1351  nfld2 = nfld/2
1352  mfld = min(17,nfld)
1353  if (ifmhd) nfld2 = nfld/4
1354 
1355  i2 = iosav/2
1356  m1 = istep+iosav-iotest
1357  mt = mod(istep+iosav-iotest,iosav)
1358  prefix = ' '
1359 
1360  if (istep.gt.iosav/2 .and.
1361  $ mod(istep+iosav-iotest,iosav).lt.nfld2) then ! save
1362  write(prefix,'(A)') 'rs_'
1363 c write(prefix,3) ks1(mfld)
1364 c 3 format('rs',a1)
1365 
1366  p66 = param(66)
1367  param(66) = 6
1368  if (ifmhd) call outpost2(bx,by,bz,pm,t,0,prefix) ! first B
1369  call prepost (.true.,prefix)
1370  param(66) = p66
1371 
1372  endif
1373 
1374  return
1375  end
1376 c-----------------------------------------------------------------------
1377  subroutine outpost(v1,v2,v3,vp,vt,name3)
1378 
1379  include 'SIZE'
1380  include 'INPUT'
1381 
1382  real v1(1),v2(1),v3(1),vp(1),vt(1)
1383  character*3 name3
1384 
1385 
1386  itmp=0
1387  if (ifto) itmp=1
1388  call outpost2(v1,v2,v3,vp,vt,itmp,name3)
1389 
1390  return
1391  end
1392 c-----------------------------------------------------------------------
1393  subroutine outpost2(v1,v2,v3,vp,vt,nfldt,name3)
1394 
1395  include 'SIZE'
1396  include 'SOLN'
1397  include 'INPUT'
1398 
1399  parameter(ltot1=lx1*ly1*lz1*lelt)
1400  parameter(ltot2=lx2*ly2*lz2*lelv)
1401  common /outtmp/ w1(ltot1),w2(ltot1),w3(ltot1),wp(ltot2)
1402  & ,wt(ltot1,ldimt)
1403 c
1404  real v1(1),v2(1),v3(1),vp(1),vt(ltot1,1)
1405  character*3 name3
1406  logical if_save(ldimt)
1407 c
1408  ntot1 = lx1*ly1*lz1*nelt
1409  ntot1t = lx1*ly1*lz1*nelt
1410  ntot2 = lx2*ly2*lz2*nelt
1411 
1412  if(nfldt.gt.ldimt) then
1413  write(6,*) 'ABORT: outpost data too large (nfldt>ldimt)!'
1414  call exitt
1415  endif
1416 
1417 c store solution
1418  call copy(w1,vx,ntot1)
1419  call copy(w2,vy,ntot1)
1420  call copy(w3,vz,ntot1)
1421  call copy(wp,pr,ntot2)
1422  do i = 1,nfldt
1423  call copy(wt(1,i),t(1,1,1,1,i),ntot1t)
1424  enddo
1425 
1426 c swap with data to dump
1427  call copy(vx,v1,ntot1)
1428  call copy(vy,v2,ntot1)
1429  call copy(vz,v3,ntot1)
1430  call copy(pr,vp,ntot2)
1431  do i = 1,nfldt
1432  call copy(t(1,1,1,1,i),vt(1,i),ntot1t)
1433  enddo
1434 
1435 c dump data
1436  if_save(1) = ifto
1437  ifto = .false.
1438  if(nfldt.gt.0) ifto = .true.
1439  do i = 1,ldimt-1
1440  if_save(i+1) = ifpsco(i)
1441  ifpsco(i) = .false.
1442  if(i+1.le.nfldt) ifpsco(i) = .true.
1443  enddo
1444 
1445  call prepost(.true.,name3)
1446 
1447  ifto = if_save(1)
1448  do i = 1,ldimt-1
1449  ifpsco(i) = if_save(i+1)
1450  enddo
1451 
1452 c restore solution data
1453  call copy(vx,w1,ntot1)
1454  call copy(vy,w2,ntot1)
1455  call copy(vz,w3,ntot1)
1456  call copy(pr,wp,ntot2)
1457  do i = 1,nfldt
1458  call copy(t(1,1,1,1,i),wt(1,i),ntot1t)
1459  enddo
1460 
1461  return
1462  end
1463 c-----------------------------------------------------------------------
1464  subroutine mfo_mdatav(u,v,w,nel)
1465 
1466  include 'SIZE'
1467  include 'INPUT'
1468  include 'PARALLEL'
1469  include 'RESTART'
1470 
1471  real u(lx1*ly1*lz1,1),v(lx1*ly1*lz1,1),w(lx1*ly1*lz1,1)
1472 
1473  real*4 buffer(1+6*lelt)
1474 
1475  integer e
1476 
1477  call nekgsync() ! clear outstanding message queues.
1478 
1479  nxyz = lx1*ly1*lz1
1480  n = 2*ldim
1481  len = 4 + 4*(n*lelt) ! recv buffer size
1482  leo = 4 + 4*(n*nelt)
1483  ierr = 0
1484 
1485  ! Am I an I/O node?
1486  if (nid.eq.pid0) then
1487  j = 1
1488  do e=1,nel
1489  buffer(j+0) = vlmin(u(1,e),nxyz)
1490  buffer(j+1) = vlmax(u(1,e),nxyz)
1491  buffer(j+2) = vlmin(v(1,e),nxyz)
1492  buffer(j+3) = vlmax(v(1,e),nxyz)
1493  j = j + 4
1494  if(if3d) then
1495  buffer(j+0) = vlmin(w(1,e),nxyz)
1496  buffer(j+1) = vlmax(w(1,e),nxyz)
1497  j = j + 2
1498  endif
1499  enddo
1500 
1501  ! write out my data
1502  nout = n*nel
1503  if(ierr.eq.0) then
1504  if(ifmpiio) then
1505  call byte_write_mpi(buffer,nout,-1,ifh_mbyte,ierr)
1506  else
1507  call byte_write(buffer,nout,ierr)
1508  endif
1509  endif
1510 
1511  ! write out the data of my childs
1512  idum = 1
1513  do k=pid0+1,pid1
1514  mtype = k
1515  call csend(mtype,idum,4,k,0) ! handshake
1516  call crecv(mtype,buffer,len)
1517  inelp = buffer(1)
1518  nout = n*inelp
1519  if(ierr.eq.0) then
1520  if(ifmpiio) then
1521  call byte_write_mpi(buffer(2),nout,-1,ifh_mbyte,ierr)
1522  else
1523  call byte_write(buffer(2),nout,ierr)
1524  endif
1525  endif
1526  enddo
1527  else
1528  j = 1
1529  buffer(j) = nel
1530  j = j + 1
1531  do e=1,nel
1532  buffer(j+0) = vlmin(u(1,e),nxyz)
1533  buffer(j+1) = vlmax(u(1,e),nxyz)
1534  buffer(j+2) = vlmin(v(1,e),nxyz)
1535  buffer(j+3) = vlmax(v(1,e),nxyz)
1536  j = j + 4
1537  if(n.eq.6) then
1538  buffer(j+0) = vlmin(w(1,e),nxyz)
1539  buffer(j+1) = vlmax(w(1,e),nxyz)
1540  j = j + 2
1541  endif
1542  enddo
1543 
1544  ! send my data to my pararent I/O node
1545  mtype = nid
1546  call crecv(mtype,idum,4) ! hand-shake
1547  call csend(mtype,buffer,leo,pid0,0) ! u4 :=: u8
1548  endif
1549 
1550  call err_chk(ierr,'Error writing data to .f00 in mfo_mdatav. $')
1551 
1552  return
1553  end
1554 c-----------------------------------------------------------------------
1555  subroutine mfo_mdatas(u,nel)
1556 
1557  include 'SIZE'
1558  include 'INPUT'
1559  include 'PARALLEL'
1560  include 'RESTART'
1561 
1562  real u(lx1*ly1*lz1,1)
1563 
1564  real*4 buffer(1+2*lelt)
1565 
1566  integer e
1567 
1568  call nekgsync() ! clear outstanding message queues.
1569 
1570  nxyz = lx1*ly1*lz1
1571  n = 2
1572  len = 4 + 4*(n*lelt) ! recv buffer size
1573  leo = 4 + 4*(n*nelt)
1574  ierr = 0
1575 
1576  ! Am I an I/O node?
1577  if (nid.eq.pid0) then
1578  j = 1
1579  do e=1,nel
1580  buffer(j+0) = vlmin(u(1,e),nxyz)
1581  buffer(j+1) = vlmax(u(1,e),nxyz)
1582  j = j + 2
1583  enddo
1584 
1585  ! write out my data
1586  nout = n*nel
1587  if(ierr.eq.0) then
1588  if(ifmpiio) then
1589  call byte_write_mpi(buffer,nout,-1,ifh_mbyte,ierr)
1590  else
1591  call byte_write(buffer,nout,ierr)
1592  endif
1593  endif
1594 
1595  ! write out the data of my childs
1596  idum = 1
1597  do k=pid0+1,pid1
1598  mtype = k
1599  call csend(mtype,idum,4,k,0) ! handshake
1600  call crecv(mtype,buffer,len)
1601  inelp = buffer(1)
1602  nout = n*inelp
1603  if(ierr.eq.0) then
1604  if(ifmpiio) then
1605  call byte_write_mpi(buffer(2),nout,-1,ifh_mbyte,ierr)
1606  else
1607  call byte_write(buffer(2),nout,ierr)
1608  endif
1609  endif
1610  enddo
1611  else
1612  j = 1
1613  buffer(j) = nel
1614  j = j + 1
1615  do e=1,nel
1616  buffer(j+0) = vlmin(u(1,e),nxyz)
1617  buffer(j+1) = vlmax(u(1,e),nxyz)
1618  j = j + 2
1619  enddo
1620 
1621  ! send my data to my pararent I/O node
1622  mtype = nid
1623  call crecv(mtype,idum,4) ! hand-shake
1624  call csend(mtype,buffer,leo,pid0,0) ! u4 :=: u8
1625  endif
1626 
1627  call err_chk(ierr,'Error writing data to .f00 in mfo_mdatas. $')
1628 
1629  return
1630  end
1631 c-----------------------------------------------------------------------
1632  subroutine mfo_outs(u,nel,mx,my,mz) ! output a scalar field
1633 
1634  include 'SIZE'
1635  include 'INPUT'
1636  include 'PARALLEL'
1637  include 'RESTART'
1638 
1639  real u(mx,my,mz,1)
1640 
1641  common /scrns/ u4(2+lxo*lxo*lxo*2*lelt)
1642  real*4 u4
1643  real*8 u8(1+lxo*lxo*lxo*1*lelt)
1644  equivalence(u4,u8)
1645 
1646  integer e
1647 
1648  umax = glmax(u,nel*mx*my*mz)
1649  umin = glmin(u,nel*mx*my*mz)
1650  if(nid.eq.0) write(6,'(A,2g13.5)') ' min/max:', umin,umax
1651 
1652  call nekgsync() ! clear outstanding message queues.
1653  if(mx.gt.lxo .or. my.gt.lxo .or. mz.gt.lxo) then
1654  if(nid.eq.0) write(6,*) 'ABORT: lxo too small'
1655  call exitt
1656  endif
1657 
1658  nxyz = mx*my*mz
1659  len = 8 + 8*(lelt*nxyz) ! recv buffer size
1660  leo = 8 + wdsizo*(nel*nxyz)
1661  ntot = nxyz*nel
1662 
1663  idum = 1
1664  ierr = 0
1665 
1666  if (nid.eq.pid0) then
1667 
1668  if (wdsizo.eq.4) then ! 32-bit output
1669  call copyx4 (u4,u,ntot)
1670  else
1671  call copy (u8,u,ntot)
1672  endif
1673  nout = wdsizo/4 * ntot
1674  if(ierr.eq.0) then
1675  if(ifmpiio) then
1676  call byte_write_mpi(u4,nout,-1,ifh_mbyte,ierr)
1677  else
1678  call byte_write(u4,nout,ierr) ! u4 :=: u8
1679  endif
1680  endif
1681 
1682  ! write out the data of my childs
1683  idum = 1
1684  do k=pid0+1,pid1
1685  mtype = k
1686  call csend(mtype,idum,4,k,0) ! handshake
1687  call crecv(mtype,u4,len)
1688  nout = wdsizo/4 * nxyz * u8(1)
1689  if (wdsizo.eq.4.and.ierr.eq.0) then
1690  if(ifmpiio) then
1691  call byte_write_mpi(u4(3),nout,-1,ifh_mbyte,ierr)
1692  else
1693  call byte_write(u4(3),nout,ierr)
1694  endif
1695  elseif(ierr.eq.0) then
1696  if(ifmpiio) then
1697  call byte_write_mpi(u8(2),nout,-1,ifh_mbyte,ierr)
1698  else
1699  call byte_write(u8(2),nout,ierr)
1700  endif
1701  endif
1702  enddo
1703 
1704  else
1705 
1706  u8(1)= nel
1707  if (wdsizo.eq.4) then ! 32-bit output
1708  call copyx4 (u4(3),u,ntot)
1709  else
1710  call copy (u8(2),u,ntot)
1711  endif
1712 
1713  mtype = nid
1714  call crecv(mtype,idum,4) ! hand-shake
1715  call csend(mtype,u4,leo,pid0,0) ! u4 :=: u8
1716 
1717  endif
1718 
1719  call err_chk(ierr,'Error writing data to .f00 in mfo_outs. $')
1720 
1721  return
1722  end
1723 c-----------------------------------------------------------------------
1724 
1725  subroutine mfo_outv(u,v,w,nel,mx,my,mz) ! output a vector field
1726 
1727  include 'SIZE'
1728  include 'INPUT'
1729  include 'PARALLEL'
1730  include 'RESTART'
1731 
1732  real u(mx*my*mz,1),v(mx*my*mz,1),w(mx*my*mz,1)
1733 
1734  common /scrns/ u4(2+lxo*lxo*lxo*6*lelt)
1735  real*4 u4
1736  real*8 u8(1+lxo*lxo*lxo*3*lelt)
1737  equivalence(u4,u8)
1738 
1739  integer e
1740 
1741  umax = glmax(u,nel*mx*my*mz)
1742  vmax = glmax(v,nel*mx*my*mz)
1743  wmax = glmax(w,nel*mx*my*mz)
1744  umin = glmin(u,nel*mx*my*mz)
1745  vmin = glmin(v,nel*mx*my*mz)
1746  wmin = glmin(w,nel*mx*my*mz)
1747  if(nid.eq.0) write(6,'(A,6g13.5)') ' min/max:',
1748  $ umin,umax, vmin,vmax, wmin,wmax
1749 
1750  call nekgsync() ! clear outstanding message queues.
1751  if(mx.gt.lxo .or. my.gt.lxo .or. mz.gt.lxo) then
1752  if(nid.eq.0) write(6,*) 'ABORT: lxo too small'
1753  call exitt
1754  endif
1755 
1756  nxyz = mx*my*mz
1757  len = 8 + 8*(lelt*nxyz*ldim) ! recv buffer size (u4)
1758  leo = 8 + wdsizo*(nel*nxyz*ldim)
1759  idum = 1
1760  ierr = 0
1761 
1762  if (nid.eq.pid0) then
1763  j = 0
1764  if (wdsizo.eq.4) then ! 32-bit output
1765  do iel = 1,nel
1766  call copyx4 (u4(j+1),u(1,iel),nxyz)
1767  j = j + nxyz
1768  call copyx4 (u4(j+1),v(1,iel),nxyz)
1769  j = j + nxyz
1770  if(if3d) then
1771  call copyx4 (u4(j+1),w(1,iel),nxyz)
1772  j = j + nxyz
1773  endif
1774  enddo
1775  else
1776  do iel = 1,nel
1777  call copy (u8(j+1),u(1,iel),nxyz)
1778  j = j + nxyz
1779  call copy (u8(j+1),v(1,iel),nxyz)
1780  j = j + nxyz
1781  if(if3d) then
1782  call copy (u8(j+1),w(1,iel),nxyz)
1783  j = j + nxyz
1784  endif
1785  enddo
1786  endif
1787  nout = wdsizo/4 * ldim*nel * nxyz
1788  if(ierr.eq.0) then
1789  if(ifmpiio) then
1790  call byte_write_mpi(u4,nout,-1,ifh_mbyte,ierr)
1791  else
1792  call byte_write(u4,nout,ierr) ! u4 :=: u8
1793  endif
1794  endif
1795 
1796  ! write out the data of my childs
1797  do k=pid0+1,pid1
1798  mtype = k
1799  call csend(mtype,idum,4,k,0) ! handshake
1800  call crecv(mtype,u4,len)
1801  nout = wdsizo/4 * ldim*nxyz * u8(1)
1802 
1803  if (wdsizo.eq.4.and.ierr.eq.0) then
1804  if(ifmpiio) then
1805  call byte_write_mpi(u4(3),nout,-1,ifh_mbyte,ierr)
1806  else
1807  call byte_write(u4(3),nout,ierr)
1808  endif
1809  elseif(ierr.eq.0) then
1810  if(ifmpiio) then
1811  call byte_write_mpi(u8(2),nout,-1,ifh_mbyte,ierr)
1812  else
1813  call byte_write(u8(2),nout,ierr)
1814  endif
1815  endif
1816  enddo
1817  else
1818 
1819  u8(1) = nel
1820  if (wdsizo.eq.4) then ! 32-bit output
1821  j = 2
1822  do iel = 1,nel
1823  call copyx4 (u4(j+1),u(1,iel),nxyz)
1824  j = j + nxyz
1825  call copyx4 (u4(j+1),v(1,iel),nxyz)
1826  j = j + nxyz
1827  if(if3d) then
1828  call copyx4 (u4(j+1),w(1,iel),nxyz)
1829  j = j + nxyz
1830  endif
1831  enddo
1832  else
1833  j = 1
1834  do iel = 1,nel
1835  call copy (u8(j+1),u(1,iel),nxyz)
1836  j = j + nxyz
1837  call copy (u8(j+1),v(1,iel),nxyz)
1838  j = j + nxyz
1839  if(if3d) then
1840  call copy (u8(j+1),w(1,iel),nxyz)
1841  j = j + nxyz
1842  endif
1843  enddo
1844  endif
1845 
1846  mtype = nid
1847  call crecv(mtype,idum,4) ! hand-shake
1848  call csend(mtype,u4,leo,pid0,0) ! u4 :=: u8
1849 
1850  endif
1851 
1852  call err_chk(ierr,'Error writing data to .f00 in mfo_outv. $')
1853  return
1854  end
1855 c-----------------------------------------------------------------------
1856  subroutine mfo_write_hdr ! write hdr, byte key, els.
1857 
1858  include 'SIZE'
1859  include 'SOLN'
1860  include 'INPUT'
1861  include 'PARALLEL'
1862  include 'RESTART'
1863  include 'TSTEP'
1864  real*4 test_pattern
1865  common /ctmp0/ lglist(0:lelt)
1866 
1867  character*132 hdr
1868  integer*8 ioff
1869  logical if_press_mesh
1870 
1871  call nekgsync()
1872  idum = 1
1873 
1874  if(ifmpiio) then
1875  nfileoo = 1 ! all data into one file
1876  nelo = nelgt
1877  else
1878  nfileoo = nfileo
1879  if(nid.eq.pid0) then ! how many elements to dump
1880  nelo = nelt
1881  do j = pid0+1,pid1
1882  mtype = j
1883  call csend(mtype,idum,4,j,0) ! handshake
1884  call crecv(mtype,inelp,4)
1885  nelo = nelo + inelp
1886  enddo
1887  else
1888  mtype = nid
1889  call crecv(mtype,idum,4) ! hand-shake
1890  call csend(mtype,nelt,4,pid0,0) ! u4 :=: u8
1891  endif
1892  endif
1893 
1894  ierr = 0
1895  if(nid.eq.pid0) then
1896 
1897  call blank(hdr,132) ! write header
1898  call blank(rdcode1,10)
1899  i = 1
1900  IF (ifxyo) THEN
1901  rdcode1(i)='X'
1902  i = i + 1
1903  ENDIF
1904  IF (ifvo) THEN
1905  rdcode1(i)='U'
1906  i = i + 1
1907  ENDIF
1908  IF (ifpo) THEN
1909  rdcode1(i)='P'
1910  i = i + 1
1911  ENDIF
1912  IF (ifto) THEN
1913  rdcode1(i)='T'
1914  i = i + 1
1915  ENDIF
1916  IF (ldimt.GT.1) THEN
1917  npscalo = 0
1918  do k = 1,ldimt-1
1919  if(ifpsco(k)) npscalo = npscalo + 1
1920  enddo
1921  IF (npscalo.GT.0) THEN
1922  rdcode1(i) = 'S'
1923  WRITE(rdcode1(i+1),'(I1)') npscalo/10
1924  WRITE(rdcode1(i+2),'(I1)') npscalo-(npscalo/10)*10
1925  ENDIF
1926  ENDIF
1927 
1928 c check pressure format
1929  if_press_mesh = .false.
1930  if (.not.ifsplit.and.if_full_pres) if_press_mesh = .true.
1931 
1932  write(hdr,1) wdsizo,nxo,nyo,nzo,nelo,nelgt,time,istep,fid0,nfileoo
1933  $ ,(rdcode1(i),i=1,10),p0th,if_press_mesh
1934  1 format('#std',1x,i1,1x,i2,1x,i2,1x,i2,1x,i10,1x,i10,1x,e20.13,
1935  & 1x,i9,1x,i6,1x,i6,1x,10a,1pe15.7,1x,l1)
1936 
1937  test_pattern = 6.54321 ! write test pattern for byte swap
1938 
1939  if(ifmpiio) then
1940  ! only rank0 (pid00) will write hdr + test_pattern
1941  call byte_write_mpi(hdr,iheadersize/4,pid00,ifh_mbyte,ierr)
1942  call byte_write_mpi(test_pattern,1,pid00,ifh_mbyte,ierr)
1943  else
1944  call byte_write(hdr,iheadersize/4,ierr)
1945  call byte_write(test_pattern,1,ierr)
1946  endif
1947 
1948  endif
1949 
1950  call err_chk(ierr,'Error writing header in mfo_write_hdr. $')
1951 
1952  ! write global element numbering for this group
1953  if(nid.eq.pid0) then
1954  if(ifmpiio) then
1955  ioff = iheadersize + 4 + nelb*isize
1956  call byte_set_view (ioff,ifh_mbyte)
1957  call byte_write_mpi(lglel,nelt,-1,ifh_mbyte,ierr)
1958  else
1959  call byte_write(lglel,nelt,ierr)
1960  endif
1961 
1962  do j = pid0+1,pid1
1963  mtype = j
1964  call csend(mtype,idum,4,j,0) ! handshake
1965  len = 4*(lelt+1)
1966  call crecv(mtype,lglist,len)
1967  if(ierr.eq.0) then
1968  if(ifmpiio) then
1969  call byte_write_mpi(lglist(1),lglist(0),-1,ifh_mbyte,ierr)
1970  else
1971  call byte_write(lglist(1),lglist(0),ierr)
1972  endif
1973  endif
1974  enddo
1975  else
1976  mtype = nid
1977  call crecv(mtype,idum,4) ! hand-shake
1978 
1979  lglist(0) = nelt
1980  call icopy(lglist(1),lglel,nelt)
1981 
1982  len = 4*(nelt+1)
1983  call csend(mtype,lglist,len,pid0,0)
1984  endif
1985 
1986  call err_chk(ierr,'Error writing global nums in mfo_write_hdr$')
1987  return
1988  end
1989 c-----------------------------------------------------------------------
#define byte_open
Definition: byte.c:35
void exitt()
Definition: comm_mpi.f:604
#define byte_write
Definition: byte.c:39
#define byte_close
Definition: byte.c:36
subroutine byte_open_mpi(fnamei, mpi_fh, ifro, ierr)
Definition: byte_mpi.f:13
subroutine byte_write_mpi(buf, icount, iorank, mpi_fh, ierr)
Definition: byte_mpi.f:65
subroutine byte_close_mpi(mpi_fh, ierr)
Definition: byte_mpi.f:84
subroutine byte_set_view(ioff_in, mpi_fh)
Definition: byte_mpi.f:97
subroutine crecv2(mtype, buf, lenm, jnid)
Definition: comm_mpi.f:333
function igl_running_sum(in)
Definition: comm_mpi.f:700
subroutine crecv(mtype, buf, lenm)
Definition: comm_mpi.f:313
subroutine gop(x, w, op, n)
Definition: comm_mpi.f:201
subroutine nekgsync()
Definition: comm_mpi.f:502
subroutine csend(mtype, buf, len, jnid, jpid)
Definition: comm_mpi.f:303
subroutine bcast(buf, len)
Definition: comm_mpi.f:431
subroutine err_chk(ierr, string)
Definition: comm_mpi.f:558
real *8 function dnekclock_sync()
Definition: comm_mpi.f:401
integer function gllel(ieg)
Definition: dprocmap.f:183
integer function gllnid(ieg)
Definition: dprocmap.f:161
subroutine addfid(fname, fid)
Definition: ic.f:2512
integer function indx1(S1, S2, L2)
Definition: ic.f:1357
function glmin(a, n)
Definition: math.f:973
function mod1(i, n)
Definition: math.f:509
subroutine icopy(a, b, n)
Definition: math.f:289
subroutine ione(a, n)
Definition: math.f:223
real function vlmax(vec, n)
Definition: math.f:396
subroutine copy(a, b, n)
Definition: math.f:260
subroutine blank(A, N)
Definition: math.f:19
function ltrunc(string, l)
Definition: math.f:494
function glsum(x, n)
Definition: math.f:861
subroutine chcopy(a, b, n)
Definition: math.f:281
real function vlmin(vec, n)
Definition: math.f:357
subroutine rzero(a, n)
Definition: math.f:208
function glmax(a, n)
Definition: math.f:960
subroutine mxm(a, n1, b, n2, c, n3)
Definition: mxm_wrapper.f:2
subroutine map2reg(ur, n, u, nel)
Definition: postpro.f:544
subroutine mfo_mdatas(u, nel)
Definition: prepost.f:1556
subroutine restart_save(iosave, nfldi)
Definition: prepost.f:1312
subroutine dump_header(excodein, p66, ierr)
Definition: prepost.f:614
subroutine mfo_outfld(prefix)
Definition: prepost.f:838
subroutine close_fld(p66, ierr)
Definition: prepost.f:782
subroutine set_outfld
Definition: prepost.f:3
subroutine out_tmp(id, p66, ierr)
Definition: prepost.f:798
subroutine mfo_outv(u, v, w, nel, mx, my, mz)
Definition: prepost.f:1726
subroutine mfo_open_files(prefix, ierr)
Definition: prepost.f:1071
subroutine prepost(ifdoin, prefin)
Definition: prepost.f:91
subroutine file2(nopen, PREFIX)
Definition: prepost.f:445
subroutine full_restart_save(iosave)
Definition: prepost.f:1282
subroutine copy4r(a, b, n)
Definition: prepost.f:568
subroutine mfo_mdatav(u, v, w, nel)
Definition: prepost.f:1465
subroutine prepost_map(isave)
Definition: prepost.f:136
subroutine mfo_outs(u, nel, mx, my, mz)
Definition: prepost.f:1633
subroutine outfld(prefix)
Definition: prepost.f:252
subroutine rzero4(a, n)
Definition: prepost.f:553
subroutine copyx4(a, b, n)
Definition: prepost.f:560
subroutine get_id(id)
Definition: prepost.f:763
subroutine outpost(v1, v2, v3, vp, vt, name3)
Definition: prepost.f:1378
subroutine outpost2(v1, v2, v3, vp, vt, nfldt, name3)
Definition: prepost.f:1394
subroutine fill_tmp(tdump, id, ie)
Definition: prepost.f:691
subroutine io_init
Definition: prepost.f:1024
subroutine restart_nfld(nfld, prefix)
Definition: prepost.f:1177
subroutine check_ioinfo
Definition: prepost.f:39
subroutine mfo_write_hdr
Definition: prepost.f:1857
function i_find_prefix(prefix, imax)
Definition: prepost.f:577