Additional Features¶
Nek5000 includes several features which make case setup, running, and data processing much easier. These capabilities are outlined here.
History Points¶
Assuming a case named foo
, a list of monitor points can be defined in the file foo.his
to evaluate velocity, temperature, pressure and passive scalars.
Values for each scalar will be spectrally interpolated to each point and appended to this file each time the subroutine hpts()
is called.
The values printed to this file are controlled by the writeToFieldFile
parameters for each quantity in the .par
file.
Depending on the number of monitoring points, you may need to increase parameter lhis
in SIZE.
Note that the monitoring points are assumed to be distributed among MPI ranks, so it may be necessary to increase lhis
above the actual number of history points requested.
Usage example:
setup an ASCII file called
foo.his
, e.g.:3 !number of monitoring points 1.1 -1.2 1.0 . . . x y z
add
call hpts()
touserchk()
Grid-to-Grid Interpolation¶
Nek5000 includes the capability to transfer a solution from one mesh to an entirely different mesh.
This allows a user to restart from an existing field file with a new mesh.
This is accomplished by calling the generic field reader which will spectrally interpolate a result file from a previous case.
If using the latest master branch from github, this can be done by specifying the int
restart option in the .par
file.
For example, to interpolate a result from case foo
:
[GENERAL]
startFrom = foo0.f00001 int
When invoked from the .par
file, the int
option is compatible with all other restart options, except for X
.
To use this feature in V19, add a call to the gfldr
subroutine in userchk
in the user routines file.
c-----------------------------------------------------------------------
subroutine userchk()
c implicit none
include 'SIZE'
include 'TOTAL'
if(istep.eq.0) call gfldr("foo0.f00001")
return
end
c-----------------------------------------------------------------------
Note that foo0.f00001
must include the coordinates, i.e.it must have been created from a run with writeToFieldFile = yes
in the [MESH]
section of the .par
file and that selection of specific fields to read is not currently supported.
That means all fields included in foo.f00001
will be overwritten.
Restart Options¶
By default, Nek5000 will read all available variables from the restart file.
Restart options can be added after the filename in the .par
file to control which variables are loaded, to manually set the time, and in the latest github version, to specify an interpolated restart (see Grid-to-Grid Interpolation).
To control which variables are read from a restart file, add an additional string composed of the following:
Variable | Characters |
---|---|
Coordinates | X |
Velocity | U |
Pressure | P |
Temperature | T |
Passive scalar ‘i’ | Si |
reset time | time=0.0 |
Example: | The following directive will load only velocity and passive scalar 2, and set the physical time to 5.0. |
---|
[GENERAL]
startFrom = foo0.f00001 US2 time=5.0
Additionally, Nek5000 supports loading multiple files simultaneously by providing a list of comma separated filenames. With the default behavior, each subsequent file will overwrite everything loaded with the previous file. However, by making use of the restart options, a Frankenstein’s Monster type of restart can be created combining fields from multiple solution files.
Example: | The following directive will load the mesh coordinates from the first file, mshfoo0.f00001 , then interpolate velocity, pressure and temperature from the second file, foo0.f00001 onto the new coordinates. |
---|
[GENERAL]
startFrom = mshfoo0.f00001 X, foo0.f00001 UPT int
Note
If a restart file contains coordinates, Nek5000 will overwrite the coordinates generated from the .re2
file. This behavior may or may not be desirable, use the restart options to control it!
Averaging¶
When running a high fidelity case with DNS or LES turbulence models, it is often necessary to time-average the solution fields to extract meaningful quantities.
This may sometimes even be useful for a URANS case as well.
Nek5000 includes a subroutine for calculating a running time-average of all the primitive variables, i.e. \(u\), \(v\), \(w\), \(p\), and \(T\), as well as the second order terms \(u^2\), \(v^2\), \(w^2\), \(uv\), \(uw\), and \(vw\).
Which can be used to reconstruct the Reynolds stresses.
To activate time-averaging, simply call avg_all
in userchk
.
c-----------------------------------------------------------------------
subroutine userchk()
c implicit none
include 'SIZE'
include 'TOTAL'
call avg_all
return
end
c-----------------------------------------------------------------------
Adding this call to userchk
will output three additional files, avgfoo
, rmsfoo
, and rm2foo
, where “foo” is your case name.
As the case is running, the running averages are stored in memory in the arrays documented here.
Warning
Averaging files are written in double precision by default and can very quickly consume a large amount of disk space!
These files will be written at the same interval as the standard restart output.
When the files are written, the averaging restarts.
The average files thus only contain averages over the window specified by writeInterval
in the .par
file.
The complete list of variables, including which file they are written to and the scalar position they occupy in that file are specified in the table below. Additionally, the width of the time-window is recorded as the physical time in each average file.
Variable | filename | scalar |
---|---|---|
\(\overline{u}\) | avgfoo0.f00000 | u-velocity |
\(\overline{v}\) | avgfoo0.f00000 | v-velocity |
\(\overline{w}\) | avgfoo0.f00000 | w-velocity |
\(\overline{p}\) | avgfoo0.f00000 | pressure |
\(\overline{T}\) | avgfoo0.f00000 | temperature |
\(\overline{\phi_i}\) | avgfoo0.f00000 | scalar i |
\(\overline{u^2}\) | rmsfoo0.f00000 | u-velocity |
\(\overline{v^2}\) | rmsfoo0.f00000 | v-velocity |
\(\overline{w^2}\) | rmsfoo0.f00000 | w-velocity |
\(\overline{p^2}\) | rmsfoo0.f00000 | pressure |
\(\overline{T^2}\) | rmsfoo0.f00000 | temperature |
\(\overline{\phi_i^2}\) | rmsfoo0.f00000 | scalar i |
\(\overline{uv}\) | rm2foo0.f00000 | u-velocity |
\(\overline{vw}\) | rm2foo0.f00000 | v-velocity |
\(\overline{uw}\) | rm2foo0.f00000 | w-velocity |
\(\overline{p^2}\) | rm2foo0.f00000 | pressure |
\(\overline{T^2}\) | rm2foo0.f00000 | temperature |
Note
avg_all
does NOT output enough information to reconstruct the turbulent heat fluxes by default.
Currently, custom user code is necessary to accomplish this.
The averaging files can then be reloaded into Nek5000 as a standard restart file for post processing. The files contain enough information to reconstruct Reynolds stresses considering that for a sufficiently large time-window at statistically steady state:
Post Processing¶
TODO…
Objects¶
TODO…