PyFR can be obtained here.
PyFR 1.11.0 has a hard dependency on Python 3.6+ and the following Python packages:
Note that due to a bug in numpy PyFR is not compatible with 32-bit Python distributions.
The CUDA backend targets NVIDIA GPUs with a compute capability of 3.0 or greater. The backend requires:
The HIP backend targets AMD GPUs which are supported by the ROCm stack. The backend requires:
The OpenCL backend targets a range of accelerators including GPUs from AMD and NVIDIA. The backend requires:
The OpenMP backend targets multi-core CPUs. The backend requires:
PyFR 1.11.0 can be installed using pip and virtualenv, as shown in the quick-start guides below.
Alternatively, PyFR 1.11.0 can be installed from source. To install the software from source, use the provided setup.py installer or add the root PyFR directory to PYTHONPATH using:
user@computer ~/PyFR$ export PYTHONPATH=.:$PYTHONPATH
When installing from source, we strongly recommend using pip and virtualenv to manage the Python dependencies.
We recommend using the package manager homebrew. Open the terminal and install the dependencies with the following commands:
brew install python3 open-mpi metis
pip3 install virtualenv
For visualisation of results, either install Paraview from the command line:
brew cask install paraview
or dowload the app from the Paraview website. Then create a virtual environment and activate it:
virtualenv --python=python3 ENV3
source ENV3/bin/activate
Finally install PyFR with pip in the virtual environment:
pip install pyfr
This concludes the installation. In order to run PyFR with the OpenMP backend (see Running PyFR), use the following settings in the configuration file (.ini):
[backend-openmp]
cc = gcc-8
cblas = /usr/lib/libblas.dylib
cblas-type = parallel
Note the version of the compiler which must support the openmp flag. This has been tested on macOS 10.14.
Open the terminal and install the dependencies with the following commands:
sudo apt install python3 python3-pip libopenmpi-dev openmpi-bin
sudo apt install metis libmetis-dev libblas3
pip3 install virtualenv
For visualisation of results, either install Paraview from the command line:
sudo apt install paraview
or dowload the app from the Paraview website. Then create a virtual environment and activate it:
python3 -m virtualenv ENV3
source ENV3/bin/activate
Finally install PyFR with pip in the virtual environment:
pip install pyfr
This concludes the installation. In order to run PyFR with the OpenMP backend (see Running PyFR), use the following settings in the configuration file (.ini):
[backend-openmp]
cc = gcc
cblas = /usr/lib/x86_64-linux-gnu/blas/libblas.so.3
cblas-type = parallel
This has been tested on Ubuntu 18.04.
PyFR 1.11.0 uses three distinct file formats:
The following commands are available from the pyfr program:
pyfr import — convert a Gmsh .msh file into a PyFR .pyfrm file.
Example:
pyfr import mesh.msh mesh.pyfrm
pyfr partition — partition an existing mesh and associated solution files.
Example:
pyfr partition 2 mesh.pyfrm solution.pyfrs .
pyfr run — start a new PyFR simulation. Example:
pyfr run mesh.pyfrm configuration.ini
pyfr restart — restart a PyFR simulation from an existing solution file. Example:
pyfr restart mesh.pyfrm solution.pyfrs
pyfr export — convert a PyFR .pyfrs file into an unstructured VTK .vtu or .pvtu file. Example:
pyfr export mesh.pyfrm solution.pyfrs solution.vtu
pyfr can be run in parallel. To do so prefix pyfr with mpiexec -n <cores/devices>. Note that the mesh must be pre-partitioned, and the number of cores or devices must be equal to the number of partitions.
The .ini configuration file parameterises the simulation. It is written in the INI format. Parameters are grouped into sections. The roles of each section and their associated parameters are described below.
Parameterises the backend with
precision — number precision:
single | double
rank-allocator — MPI rank allocator:
linear | random
Example:
[backend]
precision = double
rank-allocator = linear
Parameterises the CUDA backend with
device-id — method for selecting which device(s) to run on:
int | round-robin | local-rank
mpi-type — type of MPI library that is being used:
standard | cuda-aware
block-1d — block size for one dimensional pointwise kernels:
int
block-2d — block size for two dimensional pointwise kernels:
int
Example:
[backend-cuda]
device-id = round-robin
mpi-type = standard
block-1d = 64
block-2d = 128
Parameterises the HIP backend with
device-id — method for selecting which device(s) to run on:
int | local-rank
mpi-type — type of MPI library that is being used:
standard | hip-aware
block-1d — block size for one dimensional pointwise kernels:
int
block-2d — block size for two dimensional pointwise kernels:
int
Example:
[backend-hip]
device-id = local-rank
gimmik-max-nnz = 512
mpi-type = standard
block-1d = 64
block-2d = 128
Parameterises the OpenCL backend with
platform-id — for selecting platform id:
int | string
device-type — for selecting what type of device(s) to run on:
all | cpu | gpu | accelerator
device-id — for selecting which device(s) to run on:
int | string | local-rank
gimmik-max-nnz — cutoff for GiMMiK in terms of the number of non-zero entires in a constant matrix:
int
local-size-1d — local work size for one dimensional pointwise kernels:
int
local-size-2d — local work size for two dimensional pointwise kernels:
int
Example:
[backend-opencl]
platform-id = 0
device-type = gpu
device-id = local-rank
gimmik-max-nnz = 512
local-size-1d = 16
local-size-2d = 128
Parameterises the OpenMP backend with
cc — C compiler:
string
cflags — additional C compiler flags:
string
alignb — alignment requirement in bytes; must be a power of two and at least 32:
int
cblas — path to shared C BLAS library:
string
cblas-type — type of BLAS library:
serial | parallel
gimmik-max-nnz — cutoff for GiMMiK in terms of the number of non-zero entires in a constant matrix:
int
libxsmm-block-sz — blocking factor to use for libxsmm; must be a multiple of 16:
int
libxsmm-max-sz — cutoff for libxsmm in terms of the number of entires in a constant matrix:
int
Example:
[backend-openmp]
cc = gcc
cblas= example/path/libBLAS.dylib
cblas-type = parallel
Sets constants used in the simulation
gamma — ratio of specific heats for euler | navier-stokes:
float
mu — dynamic viscosity for navier-stokes:
float
nu — kinematic viscosity for ac-navier-stokes:
float
Pr — Prandtl number for navier-stokes:
float
cpTref — product of specific heat at constant pressure and reference temperature for navier-stokes with Sutherland’s Law:
float
cpTs — product of specific heat at constant pressure and Sutherland temperature for navier-stokes with Sutherland’s Law:
float
ac-zeta — artificial compressibility factor for ac-euler | ac-navier-stokes
float
Example:
[constants]
gamma = 1.4
mu = 0.001
Pr = 0.72
Parameterises the solver with
system — governing system:
euler | navier-stokes | ac-euler | ac-navier-stokes
where
navier-stokes requires
viscosity-correction — viscosity correction:
none | sutherland
shock-capturing — shock capturing scheme:
none | artificial-viscosity
order — order of polynomial solution basis:
int
anti-alias — type of anti-aliasing:
flux | surf-flux | flux, surf-flux
Example:
[solver]
system = navier-stokes
order = 3
anti-alias = flux
viscosity-correction = none
shock-capturing = artificial-viscosity
Parameterises the time-integration scheme used by the solver with
formulation — formulation:
std | dual
where
std requires
scheme — time-integration scheme
euler | rk34 | rk4 | rk45 | tvd-rk3
tstart — initial time
float
tend — final time
float
dt — time-step
float
controller — time-step controller
none | pi
where
pi only works with rk34 and rk45 and requires
atol — absolute error tolerance
float
rtol — relative error tolerance
float
errest-norm — norm to use for estimating the error
uniform | l2
safety-fact — safety factor for step size adjustment (suitable range 0.80-0.95)
float
min-fact — minimum factor by which the time-step can change between iterations (suitable range 0.1-0.5)
float
max-fact — maximum factor by which the time-step can change between iterations (suitable range 2.0-6.0)
float
dt-max — maximum permissible time-step
float
dual requires
scheme — time-integration scheme
backward-euler | bdf2 | bdf3
pseudo-scheme — pseudo time-integration scheme
euler | rk34 | rk4 | rk45 | tvd-rk3 | vermeire
tstart — initial time
float
tend — final time
float
dt — time-step
float
pseudo-dt — pseudo time-step
float
controller — pseudo time-step controller
none
pseudo-niters-max — minimum number of iterations
int
pseudo-niters-min — maximum number of iterations
int
pseudo-resid-tol — pseudo residual tolerance
float
pseudo-resid-norm — pseudo residual norm
uniform | l2
pseudo-controller — pseudo time-step controller
none | local-pi
where
local-pi only works with rk34 and rk45 and requires
atol — absolute error tolerance
float
safety-fact — safety factor for pseudo time-step size adjustment (suitable range 0.80-0.95)
float
min-fact — minimum factor by which the local pseudo time-step can change between iterations (suitable range 0.98-0.998)
float
max-fact — maximum factor by which the local pseudo time-step can change between iterations (suitable range 1.001-1.01)
float
pseudo-dt-max-mult — maximum permissible local pseudo time-step given as a multiplier of pseudo-dt (suitable range 2.0-5.0)
float
Example:
[solver-time-integrator]
formulation = std
scheme = rk45
controller = pi
tstart = 0.0
tend = 10.0
dt = 0.001
atol = 0.00001
rtol = 0.00001
errest-norm = l2
safety-fact = 0.9
min-fact = 0.3
max-fact = 2.5
Parameterises multi-p for dual time-stepping with
pseudo-dt-fact — factor by which the pseudo time-step size changes between multi-p levels:
float
cycle — nature of a single multi-p cycle:
[(order,nsteps), (order,nsteps), ... (order,nsteps)]
where order in the first and last bracketed pair must be the overall polynomial order used for the simulation, and order can only change by one between subsequent bracketed pairs
Example:
[solver-dual-time-integrator-multip]
pseudo-dt-fact = 2.3
cycle = [(3, 1), (2, 1), (1, 1), (0, 2), (1, 1), (2, 1), (3, 3)]
Parameterises the interfaces with
riemann-solver — type of Riemann solver:
rusanov | hll | hllc | roe | roem
where
hll | hllc | roe | roem do not work with ac-euler | ac-navier-stokes
ldg-beta — beta parameter used for LDG:
float
ldg-tau — tau parameter used for LDG:
float
Example:
[solver-interfaces]
riemann-solver = rusanov
ldg-beta = 0.5
ldg-tau = 0.1
Parameterises the line interfaces, or if -mg-porder is suffixed the line interfaces at multi-p level order, with
flux-pts — location of the flux points on a line interface:
gauss-legendre | gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing on a line interface:
int
quad-pts — name of quadrature rule for anti-aliasing on a line interface:
gauss-legendre | gauss-legendre-lobatto
Example:
[solver-interfaces-line]
flux-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre
Parameterises the triangular interfaces, or if -mg-porder is suffixed the triangular interfaces at multi-p level order, with
flux-pts — location of the flux points on a triangular interface:
williams-shunn
quad-deg — degree of quadrature rule for anti-aliasing on a triangular interface:
int
quad-pts — name of quadrature rule for anti-aliasing on a triangular interface:
williams-shunn | witherden-vincent
Example:
[solver-interfaces-tri]
flux-pts = williams-shunn
quad-deg = 10
quad-pts = williams-shunn
Parameterises the quadrilateral interfaces, or if -mg-porder is suffixed the quadrilateral interfaces at multi-p level order, with
flux-pts — location of the flux points on a quadrilateral interface:
gauss-legendre | gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing on a quadrilateral interface:
int
quad-pts — name of quadrature rule for anti-aliasing on a quadrilateral interface:
gauss-legendre | gauss-legendre-lobatto | witherden-vincent
Example:
[solver-interfaces-quad]
flux-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre
Parameterises the triangular elements, or if -mg-porder is suffixed the triangular elements at multi-p level order, with
soln-pts — location of the solution points in a triangular element:
williams-shunn
quad-deg — degree of quadrature rule for anti-aliasing in a triangular element:
int
quad-pts — name of quadrature rule for anti-aliasing in a triangular element:
williams-shunn | witherden-vincent
Example:
[solver-elements-tri]
soln-pts = williams-shunn
quad-deg = 10
quad-pts = williams-shunn
Parameterises the quadrilateral elements, or if -mg-porder is suffixed the quadrilateral elements at multi-p level order, with
soln-pts — location of the solution points in a quadrilateral element:
gauss-legendre | gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing in a quadrilateral element:
int
quad-pts — name of quadrature rule for anti-aliasing in a quadrilateral element:
gauss-legendre | gauss-legendre-lobatto | witherden-vincent
Example:
[solver-elements-quad]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre
Parameterises the hexahedral elements, or if -mg-porder is suffixed the hexahedral elements at multi-p level order, with
soln-pts — location of the solution points in a hexahedral element:
gauss-legendre | gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing in a hexahedral element:
int
quad-pts — name of quadrature rule for anti-aliasing in a hexahedral element:
gauss-legendre | gauss-legendre-lobatto | witherden-vincent
Example:
[solver-elements-hex]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre
Parameterises the tetrahedral elements, or if -mg-porder is suffixed the tetrahedral elements at multi-p level order, with
soln-pts — location of the solution points in a tetrahedral element:
shunn-ham
quad-deg — degree of quadrature rule for anti-aliasing in a tetrahedral element:
int
quad-pts — name of quadrature rule for anti-aliasing in a tetrahedral element:
shunn-ham | witherden-vincent
Example:
[solver-elements-tet]
soln-pts = shunn-ham
quad-deg = 10
quad-pts = shunn-ham
Parameterises the prismatic elements, or if -mg-porder is suffixed the prismatic elements at multi-p level order, with
soln-pts — location of the solution points in a prismatic element:
williams-shunn~gauss-legendre | williams-shunn~gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing in a prismatic element:
int
quad-pts — name of quadrature rule for anti-aliasing in a prismatic element:
williams-shunn~gauss-legendre | williams-shunn~gauss-legendre-lobatto | witherden-vincent
Example:
[solver-elements-pri]
soln-pts = williams-shunn~gauss-legendre
quad-deg = 10
quad-pts = williams-shunn~gauss-legendre
Parameterises the pyramidal elements, or if -mg-porder is suffixed the pyramidal elements at multi-p level order, with
soln-pts — location of the solution points in a pyramidal element:
gauss-legendre | gauss-legendre-lobatto
quad-deg — degree of quadrature rule for anti-aliasing in a pyramidal element:
int
quad-pts — name of quadrature rule for anti-aliasing in a pyramidal element:
witherden-vincent
Example:
[solver-elements-pyr]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = witherden-vincent
Parameterises solution, space (x, y, [z]), and time (t) dependent source terms with
rho — density source term for euler | navier-stokes:
string
rhou — x-momentum source term for euler | navier-stokes :
string
rhov — y-momentum source term for euler | navier-stokes :
string
rhow — z-momentum source term for euler | navier-stokes :
string
E — energy source term for euler | navier-stokes :
string
p — pressure source term for ac-euler | ac-navier-stokes:
string
u — x-velocity source term for ac-euler | ac-navier-stokes:
string
v — y-velocity source term for ac-euler | ac-navier-stokes:
string
w — w-velocity source term for ac-euler | ac-navier-stokes:
string
Example:
[solver-source-terms]
rho = t
rhou = x*y*sin(y)
rhov = z*rho
rhow = 1.0
E = 1.0/(1.0+x)
Parameterises artificial viscosity for shock capturing with
max-artvisc — maximum artificial viscosity:
float
s0 — sensor cut-off:
float
kappa — sensor range:
float
Example:
[solver-artificial-viscosity]
max-artvisc = 0.01
s0 = 0.01
kappa = 5.0
Parameterises an exponential solution filter with
nsteps — apply filter every nsteps:
int
alpha — strength of filter:
float
order — order of filter:
int
cutoff — cutoff frequency below which no filtering is applied:
int
Example:
[soln-filter]
nsteps = 10
alpha = 36.0
order = 16
cutoff = 1
Periodically write the solution to disk in the pyfrs format. Parameterised with
dt-out — write to disk every dt-out time units:
float
basedir — relative path to directory where outputs will be written:
string
basename — pattern of output names:
string
post-action — command to execute after writing the file:
string
post-action-mode — how the post-action command should be executed:
blocking | non-blocking
region — region to be written, specified as either the entire domain using *, a cuboidal sub-region via diametrically opposite vertices, or a sub-region of elements that have faces on a specific domain boundary via the name of the domain boundary
* | [(x, y, [z]), (x, y, [z])] | string
Example:
[soln-plugin-writer]
dt-out = 0.01
basedir = .
basename = files-{t:.2f}
post-action = echo "Wrote file {soln} at time {t} for mesh {mesh}."
post-action-mode = blocking
region = [(-5, -5, -5), (5, 5, 5)]
Periodically integrates the pressure and viscous stress on the boundary labelled name and writes out the resulting force vectors to a CSV file. Parameterised with
nsteps — integrate every nsteps:
int
file — output file path; should the file already exist it will be appended to:
string
header — if to output a header row or not:
boolean
Example:
[soln-plugin-fluidforce-wing]
nsteps = 10
file = wing-forces.csv
header = true
Periodically checks the solution for NaN values. Parameterised with
nsteps — check every nsteps:
int
Example:
[soln-plugin-nancheck]
nsteps = 10
Periodically calculates the residual and writes it out to a CSV file. Parameterised with
nsteps — calculate every nsteps:
int
file — output file path; should the file already exist it will be appended to:
string
header — if to output a header row or not:
boolean
Example:
[soln-plugin-residual]
nsteps = 10
file = residual.csv
header = true
Write time-step statistics out to a CSV file. Parameterised with
flushsteps — flush to disk every flushsteps:
int
file — output file path; should the file already exist it will be appended to:
string
header — if to output a header row or not:
boolean
Example:
[soln-plugin-dtstats]
flushsteps = 100
file = dtstats.csv
header = true
Write pseudo-step convergence history out to a CSV file. Parameterised with
flushsteps — flush to disk every flushsteps:
int
file — output file path; should the file already exist it will be appended to:
string
header — if to output a header row or not:
boolean
Example:
[soln-plugin-pseudostats]
flushsteps = 100
file = pseudostats.csv
header = true
Periodically samples specific points in the volume and writes them out to a CSV file. The plugin actually samples the solution point closest to each sample point, hence a slight discrepancy in the output sampling locations is to be expected. A nearest-neighbour search is used to locate the closest solution point to the sample point. The location process automatically takes advantage of scipy.spatial.cKDTree where available. Parameterised with
nsteps — sample every nsteps:
int
samp-pts — list of points to sample:
[(x, y), (x, y), ...] | [(x, y, z), (x, y, z), ...]
format — output variable format:
primitive | conservative
file — output file path; should the file already exist it will be appended to:
string
header — if to output a header row or not:
boolean
Example:
[soln-plugin-sampler]
nsteps = 10
samp-pts = [(1.0, 0.7, 0.0), (1.0, 0.8, 0.0)]
format = primative
file = point-data.csv
header = true
Time average quantities. Parameterised with
nsteps — accumulate the average every nsteps time steps:
int
dt-out — write to disk every dt-out time units:
float
tstart — time at which to start accumulating average data:
float
mode — output file accumulation mode:
continuous | windowed
basedir — relative path to directory where outputs will be written:
string
basename — pattern of output names:
string
precision — output file number precision:
single | double
region — region to be averaged, specified as either the entire domain using *, a cuboidal sub-region via diametrically opposite vertices, or a sub-region of elements that have faces on a specific domain boundary via the name of the domain boundary
* | [(x, y, [z]), (x, y, [z])] | string
avg-name — expression to time average, written as a function of the primitive variables and gradients thereof; multiple expressions, each with their own name, may be specified:
string
fun-avg-name — expression to compute at file output time, written as a function of any ordinary average terms; multiple expressions, each with their own name, may be specified:
string
Example:
[soln-plugin-tavg]
nsteps = 10
dt-out = 2.0
mode = windowed
basedir = .
basename = files-{t:06.2f}
avg-p = p
avg-p2 = p*p
fun-avg-varp = p2 - p*p
avg-vel = sqrt(u*u + v*v)
Parameterises constant, or if available space (x, y, [z]) and time (t) dependent, boundary condition labelled name in the .pyfrm file with
type — type of boundary condition:
ac-char-riem-inv | ac-in-fv | ac-out-fp | char-riem-inv | no-slp-adia-wall | no-slp-isot-wall | no-slp-wall | slp-adia-wall | slp-wall | sub-in-frv | sub-in-ftpttang | sub-out-fp | sup-in-fa | sup-out-fn
where
ac-char-riem-inv only works with ac-euler | ac-navier-stokes and requires
ac-zeta — artificial compressibility factor for boundary (increasing ac-zeta makes the boundary less reflective allowing larger deviation from the target state)
float
niters — number of Newton iterations
int
p — pressure
float | string
u — x-velocity
float | string
v — y-velocity
float | string
w — z-velocity
float | string
ac-in-fv only works with ac-euler | ac-navier-stokes and requires
u — x-velocity
float | string
v — y-velocity
float | string
w — z-velocity
float | string
ac-out-p only works with ac-euler | ac-navier-stokes and requires
p — pressure
float | string
char-riem-inv only works with euler | navier-stokes and requires
rho — density
float | string
u — x-velocity
float | string
v — y-velocity
float | string
w — z-velocity
float | string
p — static pressure
float | string
no-slp-adia-wall only works with navier-stokes
no-slp-isot-wall only works with navier-stokes and requires
u — x-velocity of wall
float
v — y-velocity of wall
float
w — z-velocity of wall
float
cpTw — product of specific heat capacity at constant pressure and temperature of wall
float
no-slp-wall only works with ac-navier-stokes and requires
u — x-velocity of wall
float
v — y-velocity of wall
float
w — z-velocity of wall
float
slp-adia-wall only works with euler | navier-stokes
slp-wall only works with ac-euler | ac-navier-stokes
sub-in-frv only works with navier-stokes and requires
rho — density
float | string
u — x-velocity
float | string
v — y-velocity
float | string
w — z-velocity
float | string
sub-in-ftpttang only works with navier-stokes and requires
pt — total pressure
float
cpTt — product of specific heat capacity at constant pressure and total temperature
float
theta — azimuth angle (in degrees) of inflow measured in the x-y plane relative to the positive x-axis
float
phi — inclination angle (in degrees) of inflow measured relative to the positive z-axis
float
sub-out-fp only works with navier-stokes and requires
p — static pressure
float | string
sup-in-fa only works with euler | navier-stokes and requires
rho — density
float | string
u — x-velocity
float | string
v — y-velocity
float | string
w — z-velocity
float | string
p — static pressure
float | string
sup-out-fn only works with euler | navier-stokes
Example:
[soln-bcs-bcwallupper]
type = no-slp-isot-wall
cpTw = 10.0
u = 1.0
Parameterises space (x, y, [z]) dependent initial conditions with
rho — initial density distribution for euler | navier-stokes:
string
u — initial x-velocity distribution for euler | navier-stokes | ac-euler | ac-navier-stokes:
string
v — initial y-velocity distribution for euler | navier-stokes | ac-euler | ac-navier-stokes:
string
w — initial z-velocity distribution for euler | navier-stokes | ac-euler | ac-navier-stokes:
string
p — initial static pressure distribution for euler | navier-stokes | ac-euler | ac-navier-stokes:
string
Example:
[soln-ics]
rho = 1.0
u = x*y*sin(y)
v = z
w = 1.0
p = 1.0/(1.0+x)
Proceed with the following steps to run a serial 2D Couette flow simulation on a mixed unstructured mesh:
Create a working directory called couette_flow_2d/
Copy the configuration file PyFR/examples/couette_flow_2d/couette_flow_2d.ini into couette_flow_2d/
Copy the Gmsh mesh file PyFR/examples/couette_flow_2d/couette_flow_2d.msh into couette_flow_2d/
Run pyfr to covert the Gmsh mesh file into a PyFR mesh file called couette_flow_2d.pyfrm:
pyfr import couette_flow_2d.msh couette_flow_2d.pyfrm
Run pyfr to solve the Navier-Stokes equations on the mesh, generating a series of PyFR solution files called couette_flow_2d-*.pyfrs:
pyfr run -b cuda -p couette_flow_2d.pyfrm couette_flow_2d.ini
Run pyfr on the solution file couette_flow_2d-040.pyfrs converting it into an unstructured VTK file called couette_flow_2d-040.vtu. Note that in order to visualise the high-order data, each high-order element is sub-divided into smaller linear elements. The level of sub-division is controlled by the integer at the end of the command:
pyfr export couette_flow_2d.pyfrm couette_flow_2d-040.pyfrs couette_flow_2d-040.vtu -d 4
Visualise the unstructured VTK file in Paraview
Proceed with the following steps to run a parallel 2D Euler vortex simulation on a structured mesh:
Create a working directory called euler_vortex_2d/
Copy the configuration file PyFR/examples/euler_vortex_2d/euler_vortex_2d.ini into euler_vortex_2d/
Copy the Gmsh file PyFR/examples/euler_vortex_2d/euler_vortex_2d.msh into euler_vortex_2d/
Run pyfr to convert the Gmsh mesh file into a PyFR mesh file called euler_vortex_2d.pyfrm:
pyfr import euler_vortex_2d.msh euler_vortex_2d.pyfrm
Run pyfr to partition the PyFR mesh file into two pieces:
pyfr partition 2 euler_vortex_2d.pyfrm .
Run pyfr to solve the Euler equations on the mesh, generating a series of PyFR solution files called euler_vortex_2d*.pyfrs:
mpiexec -n 2 pyfr run -b cuda -p euler_vortex_2d.pyfrm euler_vortex_2d.ini
Run pyfr on the solution file euler_vortex_2d-100.0.pyfrs converting it into an unstructured VTK file called euler_vortex_2d-100.0.vtu. Note that in order to visualise the high-order data, each high-order element is sub-divided into smaller linear elements. The level of sub-division is controlled by the integer at the end of the command:
pyfr export euler_vortex_2d.pyfrm euler_vortex_2d-100.0.pyfrs euler_vortex_2d-100.0.vtu -d 4
Visualise the unstructured VTK file in Paraview
Proceed with the following steps to run a serial 2D incompressible cylinder flow simulation on a mixed unstructured mesh:
Create a working directory called inc_cylinder_2d/
Copy the configuration file PyFR/examples/inc_cylinder_2d/inc_cylinder_2d.ini into inc_cylinder_2d/
Copy the compressed Gmsh mesh file PyFR/examples/inc_cylinder_2d/inc_cylinder_2d.msh.gz into inc_cylinder_2d/
Unzip the file and run pyfr to covert the Gmsh mesh file into a PyFR mesh file called inc_cylinder_2d.pyfrm:
zcat inc_cylinder_2d.msh.gz | pyfr import -tgmsh - inc_cylinder_2d.pyfrm
Run pyfr to solve the incompressible Navier-Stokes equations on the mesh, generating a series of PyFR solution files called inc_cylinder_2d-*.pyfrs:
pyfr run -b cuda -p inc_cylinder_2d.pyfrm inc_cylinder_2d.ini
Run pyfr on the solution file inc_cylinder_2d-75.00.pyfrs converting it into an unstructured VTK file called inc_cylinder_2d-75.00.vtu. Note that in order to visualise the high-order data, each high-order element is sub-divided into smaller linear elements. The level of sub-division is controlled by the integer at the end of the command:
pyfr export inc_cylinder_2d.pyfrm inc_cylinder_2d-75.00.pyfrs inc_cylinder_2d-75.00.vtu -d 4
Visualise the unstructured VTK file in Paraview