• No results found

Modelling of wind flow over complex terrain using OpenFoam

N/A
N/A
Protected

Academic year: 2023

Share "Modelling of wind flow over complex terrain using OpenFoam"

Copied!
115
0
0

Loading.... (view fulltext now)

Full text

New wall functions including roughness have been added to OpenFoam's current version (OF-1.5.) to be able to establish non-uniform variable roughness along the domain. It was concluded that more research is needed to get the best out of the meshing tool, as the large amount of variables to be set makes it quite difficult to control the process correctly.

INTRODUCTION

  • Energy background
  • Vattenfall
  • CFD for wind farm siting
  • Aims and tasks
  • Method
  • Lay out

The following image (Figure 1.3) shows Vattenfall's wind energy assets spread across the Nordic region, the United Kingdom and Central Europe. In the final part of the project, SnappyHexMesh (OF's post-processing tool) will be tried out to analyze the ease of using it to piece together the areas of interest.

Figure 1.3. Vattenfall’s wind assets in June 2009. 3  
Figure 1.3. Vattenfall’s wind assets in June 2009. 3  

THEORY

Atmospheric boundary layer (ABL)

  • ABL structure and depth
  • Stratification and stability
  • Turbulence
  • Governing equations
    • Coriolis force
  • Time‐Average Transport Equations
  • Roughness
  • Velocity profile approximations in the dynamic sublayer
    • Logarithmic velocity profile
    • Power‐law expression
    • Alexandrou’s formulation

In the future, stability effects should be taken into account in order to obtain more accurate results. Its components can be derived from the following diagram, where φ is the latitude of the study area.

Figure 2.3 Influence of the stability on the flow over a hill.[9]
Figure 2.3 Influence of the stability on the flow over a hill.[9]

Numerical Modeling

  • Turbulence modeling
    • Boussinesq approximation
    • K‐ε standard model
    • K‐ε RNG model
  • Discretization schemes
  • Near wall treatment: Wall functions
    • Smooth surfaces
    • Rough surfaces

For rough surfaces, some changes need to be made to account for the roughness. The current version of OpenFoam (OF-1.5) does not include roughness, so it was decided to include it in an attempt to base it on the formulation used in Fluent, as it has already been validated and proven to work well.

Table  2.2.  Standard k-ε model constants.
Table 2.2. Standard k-ε model constants.

CFD: OpenFoam

  • OpenFoam Introduction
  • Programming language of OpenFoam
    • File structure & compiling
    • Equation representation
  • OpenFoam cases
    • System folder
    • Constant folder
    • Time directories
  • Standard solvers and libraries
  • Pre‐processing
    • Gambit
    • SnappyHexMesh
  • Post‐processing with Fluent

To do this, we will need to include a folder called "Make" in the class directory, which contains two files, "options" and "files". To do this, wclean is entered in the command line and once all the old dependencies are deleted, it can be recompiled. One of the most important advantages of OpenFoam is that the equation syntax in solver applications strongly resembles differential equations.

FOAM_TUTORIALS and are divided into folders depending on the solver used in the example. Before typing the command, you must create the controlDict dictionary in the system subfolder or the command will not work. In the system folder is the SnappyHexMesh dictionary, where all variables related to mesh characteristics and quality are set.

The meshing process consists of three steps that are saved in the root folder as 1, 2 and 3 files respectively. In order to run the application like the rest of the applications in OpenFoam, all you need to do is type its name (in this case snappyHexMesh) on the command line in the root of the case and the meshing will start. An example of the foamDataToFluentDict file can be found on the attached CD in the directory “cases/Coriolis/system”).

Figure 3.2.  OpenFoam class structure.
Figure 3.2. OpenFoam class structure.

PROCESS AND RESULTS

Comparison OpenFoam vs Fluent

  • Flat plane
    • Velocity comparison OpenFoam vs Fluent 1 st order
    • Velocity comparison OpenFoam vs Fluent using 2 nd order
    • Velocity profiles
  • Bump case
    • Comparison Experimental vs OpenFoam and Fluent
  • Conclusion

The velocity profiles will be compared using 1st order and 2nd order discretization for the momentum, K and epsilon equations to be sure of the reliability of the results. In the following diagram, the velocity profiles for 1st and 2nd order are represented by the previously calculated vertical lines. The results obtained to fit the velocity profile 400 meters far from the inlet are summarized in Table 4.3.

It can be appreciated that both formulations match the results obtained with OpenFoam very well. The hump network and its surroundings can be seen in the figure below. The results obtained with OpenFoam and Fluent will be compared close to the hump, analyzing the separation along the centerline with first and second order discretization orders.

In the following images, the separation along the centerline with 1st discretization order can be seen. A final simulation is done to try to improve the results obtained with OpenFoam in the recirculation region. The results using the new turbulence model calculated by OF are presented in the following image.

Figure 4.3. Bump case.
Figure 4.3. Bump case.

Roughness implementation

  • Adding roughness as a variable in OF
    • Wall function validation
  • Pre‐processing roughness data application
    • Roughness map files
    • Roughness application

To try it out and see that it works well, run the bumpbox with the new solver (windFoam) and the rest of the edits, including the roughness set to zero. According to Fluent, to get accurate results, the roughness can never be such that the wall-adjacent cell is smaller than this. In addition, the velocity profiles for the smooth casing and roughness equal to 0.0024 meters are also plotted to see how the roughness affects the velocity profile.

However, it is not recommended to work with cases where the roughness height is much greater than the height of the first cell. The roughness data is provided in ".map" extension files and an example can be seen in Appendix A.4. When the last polygon and the last central point are reached, the Cs and Ks OF boundary condition files are rewritten with the new roughness values.

If the point meets the requirements for all the pair of consecutive vertices, it means that the point is inside the polygon, so the roughness values ​​of the polygon will be set for such a point. However, the roughness data is going to be delivered in 2D (plane x-y), so this advantage is useless. To test this, a new case is run, this will check if the roughness data has been correctly transferred to OpenFoam's boundary condition files.

Figure 4.19. Velocity distribution computed by OF in m/s.
Figure 4.19. Velocity distribution computed by OF in m/s.

Coriolis force

  • Validation of the implemented coriolis forces

Such restrictions were taken into account in the example above, so that the entrance is in the west, the outlet is in the east, the upper wall is in the north, and the bottom is in the south. The following figure shows the wind path, particles entering through the inlet wall and exiting through the bottom wall. Due to the Coriolis force, it can be seen that the wind flow does not go out through the outlet, but through the bottom wall.

If instead of being located at the northern hemisphere, the case would be at the southern hemisphere, coriolis force would affect it in the opposite direction and the wind flow would instead tend to go out through the upper wall. We will compare the velocity in x- and y-direction of the particle drawn in red in Figure 4.27, the one entering the domain at the top of the inlet wall. Therefore, taking into account that the coriolis force is the only force applied to the wind flow, the velocity for the particles is calculated.

During time intervals of 10 seconds (a period of time that is short enough compared to the time required for changes in acceleration and velocity to occur), the velocity and acceleration of the wind stream are considered to be constant, i.e. uniformly accelerated during such intervals you can use motion formulas. The following figure shows the path of the particles along the domain according to the theoretical formulas and with OpenFoam. We can see that the velocities and trajectory of the particle are exactly the same, so the Coriolis force formulation is considered to be correct and the solver is working correctly.

Figure 4.28. Coriolis validation case.
Figure 4.28. Coriolis validation case.

Gravity force

As might be expected, the pressure varies linearly with respect to the vertical and reaches its highest value at the inlet.

Figure 4.32. Static pressure depending on height (z).
Figure 4.32. Static pressure depending on height (z).

SnappyHexMesh

  • Bump case with SnappyHexMesh
  • Askervein

Then the value of the "locationInMesh" variable in the SnappyHexMeshDict will determine the mesh volume, which will be the one containing such a point. In our case, the locanInMesh point will be located in the top volume, that is, volume A. As can be seen in the figures above, the mesh looks very good and the seven surface layers were correctly created.

Many research works have been carried out on such a hill due to the extent of measures taken in the 80s. A fairly coarse background mesh (about 800,000 cells) is created that meets the requirement of being split into two volumes by the stl surface as commented in the previous case. Such a net is composed of 20 meter high, long and wide cells, and its location in relation to the hill can be seen in the following plan view.

The mesh is analyzed in five planes (a, b, c, d and e) which are defined in the following image. Such a region can be seen and compared in the following figures cut through the Plane e for 3 and 7 added layers respectively. 15 The OpenFoam user guide barely explains how to use it and no tutorials are available anymore.

Figure 4.33. Bump case displayed with Global mapper and background mesh wireframe.
Figure 4.33. Bump case displayed with Global mapper and background mesh wireframe.

CONCLUSION

FUTURE WORK

The energy contained in the flowing wind in a given period can be calculated by multiplying the above equation by the time period in question. Several efficiencies must be factored into the equation, reducing the amount of energy obtained. As commented in the theoretical part, from the value of y+ to (y_lam), there is a change in flow behavior from laminar flow (viscous sublayer) to fully turbulent (turbulent sublayer).

In the table below the y+ values ​​for which the turbid layer starts (depending on the E value) calculated from Eq. After typing roughnessToFoam on the command line, the main file (roughnessToFoam.C) will call other files to first read the mesh, read the map file, and finally write the roughness values ​​according to the map file to the files border conditions. The application name (roughnessToFoam) will be typed on the command line in the main folder of the case where we want to add the roughness.

It is recommended to set all roughness values ​​to zero, because the application will only write the values ​​for the areas defined in the map file, so the rest of the cells will not be changed and if before running the application a different value is used than zero was set, it will remain that way. We need the points, faces, etc. in the. Coordinates relative to the mesh coordinate system for the point (1,0) in the map file coordinate system (a1 and b1 in Figure B.1 respectively).

If necessary, this value can be modified by changing the "file" value in the ReadmapFile.H line 6. If necessary, this value can be modified by changing the "nam" value in the ReadmapFile.H 12 line.

Figure A.1. (a) Weibull wind distribution f(v) and (b) Wind rose.  16
Figure A.1. (a) Weibull wind distribution f(v) and (b) Wind rose. 16

Figure

Figure 1.3. Vattenfall’s wind assets in June 2009. 3  
Figure 2.1.  Atmospheric boundary layer structure.
Figure 2.3 Influence of the stability on the flow over a hill.[9]
Figure 4.1. Flat plane case.
+7

References

Related documents