Postprocessing/Output in Nek: Extracting values at certain spatial values and time values

[sohaib malik]

Hi Neks,

I am a new user and gradually progressing with case setup courtesy helpful advice from the group members. Please apologize for a long email, being a novice these basic questions are very important and will be very helpful. Just wanted to cover post processing in one email which can benefit other new users. I have tried to read output flow files in VisIt and do some post processing. Now, I am trying to do output post processing inside Nek i.e. Extracting values at certain spatial values and time values. On the lines of advice in some posts, I am trying to use write (6,*) command, outpost, hpts(), loadf_fld(), cp_fld(), interp_nfld(), inside usercheck(), etc. But facing some issues being a novice regarding how they work, and simple examples of the codes of them to start. 

1. I tried to run e3q case to see how postprocessing works,  

err0.f**1 file is generated when Nek is in postprocessing mode timesteps=0 in par file, and read in VisIt, what is to plot. Can I have the values extracted arranged in columns w.r.t. say domain like (x,y, or z) one column and values like (temperature,velocity,etc.) in other column. Can these be read and a x-y plot can be generated through any command. I saw 'grep' and 'gnuplot' command but unable to understand, how it works. Some useful lines of codes to start may be helpful.

I included part of code from e3q usrfile and logfile to request what highlighted lines means. 

2. I modified hillp tutorial to remove the mesh modification in userdat2() and make it a simple pipe. I want to include some operations of postprocessing with comands like outpost(). Afterwards, I included those lines in a case of simple pipe flow modifying hillp to simple pipe to understand postprocessing as mentioned in Sr. No.2 below, and I find some errors using outpost(), write, etc. commands that are also highlighted.

3. I have read in one of the posts about using grep command to read some values in a file and then plot using gnuplot. Can anyone please provide example of using that command.

4. I want to get values of say velocity at different (y) values for a certain (x), to say plot parabolic profile of velocity in a pipe. how to get those velocity values, outpost or write them to file, and generate x-y-plot. I have read some posts using some loop and getting values in the form vx(i,1,1,1).

5. Use of hpts(). I wrote .his file with 6 points and included hpts() in usrcheck(). I did find some lines in logfile showing 'dumping history points' at every time step. Now, the question is how to use them, as in VisIT I can load the fieldfiles and plot pseudo color, etc. But, how to use these history point values at certain spatial values w.r.t. time using any command, to write them to file, and plot x-y plots, etc.

6. Similarly another command load_fld(), I can understand that it will be used to load a flow field file, then question will be similar to point 4, how to extract say the domain x, y pts (I think interpolate command) is for this, and also velocity values and then outpost, write them, generating x-y plots, etc.

<e3q usrcheck() in .usr file>

if (istep.eq.0) then        !  Reset velocity & pressure to eliminate

         if(nid.eq.0) write(6,*) 'setting vx,vy,pr ', istep, time, visc 

         call copy (vx,ue,n)      !  start-up contributions to

      endif

      if (istep.eq.0) call outpost(ue,ve,vx,pe,t,'   ')

      call sub3   (ud,ue,vx,n)

      if (istep.eq.nsteps) call outpost(ud,vd,vx,pd,t,'err')

      umx = glamax(vx,n)

      vmx = glamax(vy,n)

      pdx = glamax(pd,n2)

      if (nid.eq.0) then

         write(6,11) istep,time,udx,umx,uex,u0,'  X err'

         write(6,11) istep,time,vdx,vmx,vex,v0,'  Y err'

         write(6,12) istep,time,pdx,pmx,pex,'                P err'

   11    format(i5,1p5e14.6,a7)

   12    format(i5,1p4e14.6,a21)

      endif

<logfile output for nsteps=0, postprocessing mode>

nsteps=0 -> skip time loop

 running solver in post processing mode

 call userchk

 setting vx,vy,pr            0   0.0000000000000000        5.0000000000000001E-003

        0  0.0000E+00 Write checkpoint

       FILE:/home/mtdl/Nek5000/run/e3q-try/e3q0.f00001                                                                                          

 min/max:   0.0000       6.2832       0.0000       6.2832       0.0000       1.0000    

 min/max:  -1.0000       3.0000      -1.3771       1.9771      -1.0000       3.0000    

 min/max:  -3.6249       1.3410    

       

        0  0.0000E+00 Write checkpoint

       FILE:/home/mtdl/Nek5000/run/e3q-try/erre3q0.f00001                                                                                       

  done :: userchk

<my usercheck subroutine on simple pipe case to try postprocessing like e3q>

      subroutine userchk

      include 'SIZE'

      include 'TOTAL'

      include 'RESTART'

      integer n,n2

      n = nx1*ny1*nz1*nelv

      n2=nx2*ny2*nz2*nelv

      if (istep.eq.0) call outpost(vx,vy,vz,pr,t,'out')

       

      call hpts()

         

      if (mod(istep,2).eq.0) 

      call outpost(vx,vy,vz,pr,t,'values')

      endif

      

      if (nio.eq.0) write(6,1) vx,vy,vz

      1 format(i9,1p5e12.4,'vx,vy,vz')

      return

      end

<terminal output showing errors>

         endif                                                          

            1

Error: Expecting END SUBROUTINE statement at (1)

/home/mtdl/Nek5000/run/pipetry3/pipetry3.f:135.72:

      if (mod(istep,2).eq.0)                                            

                                                                        1

Error: Cannot assign to a named constant at (1)

makefile:129: recipe for target 'usrfile' failed

make: *** [usrfile] Error 1

      1 format(i9,1p5e12.4,'vx,vy,vz')                                  

       1

Error: Invalid character in name at (1)

<log file showing just success of using call hpts() subroutine>

call userchk

 dump history points

 reading history points

 found            6  points

 done :: userchk

 Final time step =    1.4910192403760902E-003

             Solving for fluid

      26467  Project PRES           8.4160E-08   8.4324E-08   1.0019E+00   8   8

      26467  PRES gmres         1   1.1452E-08   1.6197E-08   1.0000E-05   1.4453E-04   2.7767E-04    F

      26467  Hmholtz VELX       7   9.6702E-09   1.3325E-02   1.0000E-08

      26467  Hmholtz VELY       1   5.4904E-10   4.4163E-09   1.0000E-08

      26467  PRES gmres         1   2.0396E-15   4.0804E-15   1.0000E-05   1.4081E-04   2.8700E-04    F

      26467  Hmholtz VELX      10   1.8927E-09   9.9940E-01   1.0000E-08

      26467  Hmholtz VELY       1   1.6374E-16   1.2760E-15   1.0000E-08

      26467  Volflow X              2.0000E+02   1.3333E-02   5.0974E-05   2.9999E+00   3.0000E+00

             L1/L2 DIV(V)          -5.6816E-17   6.4474E-10

             L1/L2 QTL              0.0000E+00   0.0000E+00

             L1/L2 DIV(V)-QTL      -5.6816E-17   6.4474E-10

      26467  Fluid done  2.0000E+02  3.0880E-03

 dump history points

    26467  2.0000E+02 Write checkpoint

       FILE:/home/mtdl/Nek5000/run/pipetry3/pipetry30.f00010                                                                                    

<>

Kind help is highly appreciated. I hope this help in doing postprocessing will be a big step ahead in my Nek5000 journey. 

Sohaib

[stek...@gmail.com]

I can't really help you with most of your questions, but the compilation error you show is most likely because you're missing a 'then' in your if statement:

      if (mod(istep,2).eq.0) then

      call outpost(vx,vy,vz,pr,t,'values')

      endif

For plotting custom csv-files there are many data processing tools that you can use, e.g. gnuplot, matlab, mathematica, ...