the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
SANDWake3D: A 3D parabolic RANS solver for atmospheric boundary layers and turbine wakes
Abstract. Despite many recent advances, modeling wind turbine wakes using semi-empirical and analytical models still face challenges when dealing with more complicated situations involving wind shear, veer, atmospheric stratification, and wake superposition. To address these limitations, this study introduces a three-dimensional, parabolic Reynolds Averaged Navier Stokes (RANS) k-epsilon formulation which includes an atmospheric boundary layer model and an actuator disk model for turbine wakes. The full three-dimensional solution for the velocity, temperature, and turbulence variables are efficiently solved through an alternating direction implicit scheme that requires orders of magnitude less computational resources than traditional high fidelity approaches. The results of the parabolic RANS model are compared to the equivalent large-eddy simulations (LES) and semi-empirical wake models at different wind speeds under stable atmospheric conditions with veer and shear. For the single turbine wake the RANS model was able to capture the wake deficit behavior, including the wake stretching and skewing that was observed in the LES. The distribution of the wake turbulence in the RANS model also agreed with results from the higher fidelity simulations. In simulations of a two-turbine, directly waked configuration, the new RANS model was able to handle the wake superposition behavior without difficulty, and also correctly modeled the corresponding increase in wake turbulence when compared to LES. Lastly, a demonstration of the RANS model on a 9-turbine, 3 row wind farm is shown and compared to LES.
- Preprint
(3030 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
-
CC1: 'Comment on wes-2025-249', Blondel Frédéric, 27 Nov 2025
-
AC1: 'Reply on CC1', Lawrence Cheung, 27 Nov 2025
Dear Frédéric,
That is correct, for the parabolic approach we studied, the induction zone upstream of the turbine is not included. We briefly considered including an induction model similar to the one in https://onlinelibrary.wiley.com/doi/abs/10.1002/we.2956, but that is something to be done in a separate study. In our formulation, we are modeling the effects of the turbine rotor through an actuator disk model, where the disk forces are applied at exactly where the turbine rotor is placed (i.e., at x/D=0). That avoids having to prescribe an explicit momentum deficit in the velocity profile themselves, and should be general enough to handle partially waked or fully waked scenarios. The details of this are described in section 2.4 of the manuscript.
Regarding the Monin-Obukhov length L, there are two limits that should be considered when using the similiarty theory for the ABL profiles. The first case is when L approaches ∞ (either positive or negative infinity), which corresponds to when the ABL is neutrally stratified. In that situation, the relations given in section 2.2 are well defined, as the non-dimensional profiles collapse to unity. The second limit to consider is when L approaches zero, which corresponds to when the ABL is unstably or stably stratified to a inifinite degree, possibly to due to the surface heat flux approaching very large, unreasonable values. In that case, I suspect the Monin-Obukhov similarity theory would break down and the profiles would not really be classified as atmospheric boundary layers anymore, so the equations in section 2.2 shouldn't be used.
Hopefully this helps clarify things, but let us know if you have any other questions.
Cheers,
Lawrence
Citation: https://doi.org/10.5194/wes-2025-249-AC1
-
AC1: 'Reply on CC1', Lawrence Cheung, 27 Nov 2025
-
RC1: 'Comment on wes-2025-249', Anonymous Referee #1, 09 Dec 2025
The comment was uploaded in the form of a supplement: https://wes.copernicus.org/preprints/wes-2025-249/wes-2025-249-RC1-supplement.pdf
-
RC2: 'Comment on wes-2025-249', Anonymous Referee #2, 29 Dec 2025
The submitted paper develops a parabolised Reynolds Averaged Navier Stokes (RANS) model and code for predicting wind turbine wakes. The approach follows similar philosophy to WakeBlaster and the curled wake model, both published several years ago. A main difference is how the pressure Poisson is formulated. The model is calibrated to LES data and then compared to the same LES data for a single turbine wake. The model is then compared to similar LES data for two turbines and for a wind farm. Overall, some qualitative results appear promising, but more comparison ABL cases as well as more detailed and clearer results and comparisons are required.
Main comments
- The pressure Poisson equation neglects streamwise derivatives. This appears to be a major limitation, as it will not correctly resolve the induced velocity in the streamwise direction which is critical to predict the power production and the wake velocity at the outlet of the inviscid streamtube. The authors should calculate the resolved induction in their model and compare it to LES. Also the maximum wake velocity as a function of x should be plotted from x/D=-5 to 15 for each single turbine LES case and model.
- ABL physics inconsistency: A substantial portion of this paper focuses on veer effects, but the presented equations neglect Coriolis forces which are largely responsible for veer. The authors also use a turbulence model that is for surface layer, which neglects Coriolis forces and veer.
- The wake models for comparison are too simple. The authors should compare their results to the WakeBlaster or the curled wake model. Also, given their focus on veer the authors should compared to the Gaussian model with the Abkar approach for veer (Abkar et al (2018) https://doi.org/10.3390/en11071838). Also, quantitative error metrics should be presented.
- It appears that the RANS model is calibrated to the LES data and then compared back to the same LES data as validation. The RANS model must be compared to unseen test data. Some very brief comparisons are shown to a wind farm case but these comparisons are brief and not sufficiently quantitative. Additional unseen single turbine LES cases with different ABL conditions should be simulated and tested against the model.
- Section 2.2: The authors describe a boundary condition methodology which is something like a hybrid between wake modeling and RANS modeling approaches and must be more clearly explained and justified. Specifically, with a prescribed turbulence model and boundary conditions for k and \epsilon, this will imply a particular solution to the RANS equations with some velocity. However, the authors also claim they can achieve a “desired” veer profile. This is inconsistent. Any prescribed profile that is not the solution to the RANS & turbulence model combination at the inlet will evolve within the domain away from the prescribed profile, as is well established in RANS (Parente et al (2011) https://doi.org/10.1016/j.jweia.2010.12.017)
Minor comments
- Line 5: “The full three-dimensional solution for the velocity, temperature, and turbulence variables are efficiently solved through an alternating direction implicit scheme that requires orders of magnitude less computational resources than traditional high fidelity approaches.”
The authors describe that their method is orders of magnitude cheaper than “high fidelity approaches” but it is not clear what approaches are referred to. Non-parabolic RANS or LES? - Line 30: “More recent work (Narasimhan et al., 2022, 2025) has extended analytical wake models to include atmospheric shear and veer, but consistently accounting for these effects in interacting wakes or in the wake-added turbulence behavior remains an open question.”
Reference to original (Abkar et al (2018) https://doi.org/10.3390/en11071838) paper on wake modeling in shear and veer is absent. - Line 49: “This restricts the ability of the model to handle effects such as wake-ABL interactions or wake skewing and veering.”
The authors state that initializing the wakes in the velocity field, rather than a body forcing, restricts the ability to handle “wake-ABL interactions or wake skewing and veering.” This statement is not sufficient clear about what specific mechanisms cannot be handled by the velocity field approach. This must be rephrased to be clearer. Wake skewing and veering can be capturing via the velocity field (e.g. Abkar & Porte-Agel (2018)) and “wake-ABL interactions” is far too vague. - Line 65: “Previous studies have shown that the k −ϵ model can accurately simulate stratified ABL conditions (Alinot and Masson, 2005), [...]”
I do not agree that the k-\epsilon model is accurate in stratified ABL conditions. This is too sweeping. It can be accurate, but has shown clear limitations as well (van der Laan et al (2017) https://doi.org/10.1002/we.2017) - Equations (1b)-(1d): The Coriolis forcing and geostrophic pressure gradient force are absent from these equations. This is concerning because the introduction motivated this paper by seeking to capture stratified ABL conditions and “veer.” This RANS model will not produce veer from Coriolis effects.
- Line 98: “values of An are given in table 1 of Alinot and Masson (2005).”
The 5 tunable parameters for An must be shown in this paper, not in a reference to a paper from 2005. - Equation 6: the authors solve a “parabolic” RANS equation that solves an elliptic pressure Poisson equation where propagation of forcing in x is neglected. This approximation of the Poisson equation to neglect x-derivatives may lead to major errors in the degree to which the induction is resolved by this model. This is a concern of this approach and must be tested in detail.
- Equation (19): this actuator disk model is not well posed for any simulation besides a single turbine, because the freestream velocity (U_\infty) is not known in wind turbine arrays. Please clarify how this is done. Two turbine and wind farm cases are presented and U_\infty should not have been used to calculate C_T in the RANS model.
- Section 3.1: AMR-Wind with ALM turbines are used for comparison. The actuator disk model should be used first, to eliminate the drivers of error from differences between uniformly loaded ADM and ALM, instead focusing on parabolised RANS vs. LES
- Figure 6: The results seem to indicate that the initial wake magnitude is too weak in the parabolised RANS, and then the wake recovers too slowly. Again, the RANS was calibrated to the LES for x/D between 4 and 6. So the model is likely just overfit to agree well in these regimes, but will have the wrong wake recovery. This will cause large error outside of 4-6 D downwind.
- Figure 8: The same results are seen, the wake is initially too weak, and then recovery is too slow. This would be clear if the velocity was plotted for all x values, similar to Figure 5 from Bastankhah (2014) https://doi.org/10.1016/j.renene.2014.01.002. The authors should make this plot for their comparisons with LES and other wake models.
Citation: https://doi.org/10.5194/wes-2025-249-RC2 - The pressure Poisson equation neglects streamwise derivatives. This appears to be a major limitation, as it will not correctly resolve the induced velocity in the streamwise direction which is critical to predict the power production and the wake velocity at the outlet of the inviscid streamtube. The authors should calculate the resolved induction in their model and compare it to LES. Also the maximum wake velocity as a function of x should be plotted from x/D=-5 to 15 for each single turbine LES case and model.
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 254 | 66 | 23 | 343 | 13 | 12 |
- HTML: 254
- PDF: 66
- XML: 23
- Total: 343
- BibTeX: 13
- EndNote: 12
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
Dear Authors,
I have a question regarding your use of a forward-in-space numerical algorithm in your preprint. Since such algorithms only allow downstream information propagation, I assume there is no induction zone upstream of the actuator disk. In the absence of an induction zone, it is standard to inject the momentum deficit at the end of the near wake (as in Bradstock and Schlez, WES, 2020), or at the rotor plane if pressure is being solved for. Could you clarify whether this is the approach you have taken, or if I have misunderstood your method?
Additionally, I noticed that Equations (9) and (10)—and possibly others—are ill-posed when L equals zero, as this would result in division by zero.
Thank you for your time and clarification.