Input setting for wind simulation

Hi there,

Hope you are doing well.

Recently, I used PyFR to simulate the aerodynamics of a wind farm. First, I have generated a straightforward structural grid (500m*300m*300m) with a cell size of (10m*10m*10m). In my input .ini file, I specified the boundary conditions for the top surface (z=300) and the bottom surface (z=0) as ‘slp-adia-wall’ and ‘slp-adia-wall’, respectively. The boundary conditions for the inflow surface (x=0), the back surface (y=300), and the front surface (y=0) are all set to ‘char-riem-inv’. The outflow surface (x=1000) has a boundary condition of ‘sub-out-fp’.

  1. When I set the air static pressure Pc = 101325 Pa, PyFR crashes very quickly (e.g., at t=0.19s), indicating the appearance of NaN values in the calculations. However, if Pc is changed to 1, PyFR runs smoothly. Could you help me understand the cause of this issue?

  2. I’m uncertain whether the crash is related to the boundary conditions. Could you please verify if the boundary conditions I set are correct? For example, should the inflow condition be ‘sub-in-frv’ or ‘char-riem-inv’? What about the other boundaries?

I have attached the input .ini file.

Looking forward to your reply

Many thanks,

Yubo

The input .ini file:

[backend]
precision = double

[constants]

;; air
gamma = 1.4
mu = 0.0000181
Pr = 0.72

cp = 10.0
Uu = 8.0
Uv = 0
Uw = 0
Pc = 101325
; Pc = 1
Tw = 1.0


[solver]
system = navier-stokes
order = 2
anti-alias = flux
viscosity-correction = none
shock-capturing = artificial-viscosity

[solver-time-integrator]
formulation = std
scheme = rk4
controller = none
tstart = 0.0
tend = 100.0
dt = 0.01

[solver-interfaces]
riemann-solver = rusanov
ldg-beta = 0.5
ldg-tau = 0.1

[solver-artificial-viscosity]
max-artvisc = 0.01
s0 = 0.01
kappa = 5.0

[solver-interfaces-line]
flux-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre

[solver-interfaces-tri]
soln-pts = williams-shunn
quad-deg = 10
quad-pts = williams-shunn

[solver-interfaces-quad]
flux-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre

[solver-elements-tri]
soln-pts = williams-shunn
quad-deg = 10
quad-pts = williams-shunn

[solver-elements-quad]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre

[solver-elements-hex]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = gauss-legendre

[solver-elements-tet]
soln-pts = shunn-ham
quad-deg = 10
quad-pts = shunn-ham

[solver-elements-pri]
soln-pts = williams-shunn~gauss-legendre
quad-deg = 10
quad-pts = williams-shunn~gauss-legendre

[solver-elements-pyr]
soln-pts = gauss-legendre
quad-deg = 10
quad-pts = witherden-vincent

[soln-plugin-nancheck]
nsteps = 10

[soln-plugin-writer]
dt-out = 0.1
basedir = .
basename = wind_farm_3d-{n:03d}
region = [(0, 0, 0), (500, 300, 300)]

[soln-bcs-bottom]
type = no-slp-isot-wall
cpTw = cp*Tw
u = Uu
v = Uv
w = Uw

[soln-bcs-top]
type = slp-adia-wall

[soln-bcs-back]
type = char-riem-inv
rho  = 1.225
u    = Uu
v    = Uv
w    = Uw
p    = Pc

[soln-bcs-front]
type = char-riem-inv
rho  = 1.225
u    = Uu
v    = Uv
w    = Uw
p    = Pc

[soln-bcs-inflow]
type = char-riem-inv
rho  = 1.225
u    = Uu
v    = Uv
w    = Uw
p    = Pc

[soln-bcs-outflow]
type = sub-out-fp
p = Pc


[soln-ics]
rho = 1.225
u = Uu
v = Uv
w = Uw
p = Pc

There are a couple of things that I can see that you probably want to try.

  1. Based on your inlet conditions you are running at a Mach number of ~0.02, ie you are pretty firmly in the incompressibel regime. You should problebly try using ac-navier-stokes.
  2. If you want to continue using navier-stokes you should probably use adaptive time stepping, at such a low Mach number the system will be very stiff.
  3. Using such a large value of pressure can lead to numerical precision issues, you might want to try non-dimensionalising things. So that the numbers are smaller.
  4. It looks like you are trying to use the artificial viscosity to stabilise the system, it is intended to be used only as a shock capturing approach and not as a smoother like you see in many FV codes.

Another suggestion for you, given the grid size it seems like you are really trying to simulate atmospheric turbulence, as you don’t seem to have actually meshed the wind turbines. If you are just trying to simulate the atmosphere/atmospheric turbulence, PyFR maybe isn’t the best tool for this. You should maybe try using the open source tool WRF: GitHub - wrf-model/WRF: The official repository for the Weather Research and Forecasting (WRF) model

I’ve used it quite a lot, and for your use case you should be able to do that on a laptop.

Thanks very much for your useful suggestions, Will. I will first try to use the ac-navier-stokes, and then I will feedback my results.

Best wishes,
Yubo