Difference between revisions of "Sig WindE/Validation Cases"

Latest revision as of 09:38, 2 November 2013

The aim of this page is to hold validation cases for ABL (Atmospheric Boundary Layer) flows that are of interest to anybody starting to model ABL flows in OF. The simulations I perform, and refer to mainly in this text are RANS with k-epsilon turbulence closure models, using simpleFoam. Most of the first simulations that are placed in this page were done according to the excellent Martinez 2011 Ms.C. thesis [1] All of the boundary conditions and mesh creation suggestions were followed by from his work and the references within it.

(0) Maintaining a logarithmic inlet profile over a 2D domain

The wall function boundary condition implemented in OF is the Nikuradse's grain of sand roughness. It produces an anomaly in the k profile in the first cell above the ground which propagates upwards and down-flow. This is discussed extensively in the literature, and alternative wall functions have been suggested such as Richards and Hoexy 1993 (RH). There has lately been a very good documentation of testing this anomaly by Martinez in his Ms.C. thesis [2].

File:BL2D.tar checked on OF 2.1

Attached case ("BL2Dtemplate") has python scripts using the pyFoam library to run different z0 values (all files included in this page require the pyFoam library). run with ./Allrun .

(1) Flow over Isolated 2D Hill (checked on OF 2.1)

The 2D simple hill from Martinez 2011 whose shape is given by

$y\left(x\right) = -\frac{H}{6.04844}\left[J_0(\Lambda)I_0(\Lambda\frac{x}{a})-I_0(\Lambda)J_0(\Lambda\frac{x}{a})\right] \; \mathrm{for} \; |x-x_0| \leq a$

$y\left(x\right) = 0 \; \mathrm{for} \; |x-x_0| > a$

The mesh is created with blockMesh, and a template file, and also a python script that writes the hill shape as a polyline object in the blockMeshDict. The mesh grading was done by breaking the control box into 3 segments - upwind, hill area and downwind (unlike the Martinez mesh grading tanh scheme). The case with the python script, ready to run (./Allrun) is attached here

File:2DBump.tar (checked on OF 2.1)

It compares the results with figures 21a,21b and 21c from Martinez 2011.

Another 2D case I tried simulating though unsuccessfully is mentioned here Sig_WindE_-_Validation_Cases_RUSHIL.

(2) Flow over Isolated 3D Hill (checked on OF 2.1.1)

The same geometry as used by Martinez 2011, with the same recommendations. The mesh is created by snappyHexMesh. The STL file is created in salome with this python script File:3DBump salome script.tar. The case with an Allrun file is attached here

File:3DBump.tar (checked on OF 2.1.1)

(3) Flow over real terrain - Askervein hill (checked on OF 2.1.1)

The Askervein hill [3] and [4] is a low isolated eliptic hill located in the island of south Uist in Scotland. In 1982-3 a measurement campaign was conducted there, and it has since been used as a bench mark case for linear and non-linear (CFD) models for flow and pollutant dispersion simulations. A simulation case is presented here with use of the meshing recommandations from Martinez 2011. The terrain is massaged (smoothed) using SAGA GIS open source software, as will be described below.

Creating the STL surface

In this case a real terrain is used. This means that the source of data is some DEM (Digital Elavation Model), and an additional step of smoothing the boundaries of the terrain has to be taken, in order to incorporate the boundary conditions used so far in the previous simulation, since they are appropriate only to flat terrain.

(a) For this simulation SRTM data was used, specifically this file [5]. It is probably not refined enough (90 meter resolution). In the OF7 workshop Carlos Peralta pointed out to another source of data here [6] (b) The DEM was uploaded into SAGA [7] through the Modules–>File–>Grid–>Import–>Import USGS SRTM Grid (picture of the result [8])

(c) Dfining the coordinate system - Moduls–>Projection–>Set Coordinate Reference System , under EPSG code insert 4326 (WGS 84)

(d) Projection - the lat/long information that the SRTM came in has to be projected unto a flat surface. I used the UTM projection (Universal Transverse Mercator). Modules–>Projection–>coordinate transformation grid . Under ESPG code chooce the one for the UTM zone that fits the Askervein hill. Finding this can be done through looking in UTM zone maps - [9]. Ours is V29. relating ESPG codes can be done according to this table for instance [10]. The ESPG code here is 32629

Smooting - the area around the simulation zone is extended and smoothed, that is - transformed into a plane of constant height. This is done here with a 1/r^2 transformation, that is continous at the height but creates a discontinouity in the slope.

(e) The area is cut and interpolated in one step: Modules–>Grid–>Gridding–>Spline Interpolation–>Multilevel B spline Interpolation (from grid) , with 5 km and 30 km grid extents, with a cell size of 50 meter. Next, a constant grid is created - the one the smoothing is done unto. A grid of height of 5 meter is chosen (according to the Askervein map and Martinez's experience) Modules–>Grid–>Construction–>create constant grid, with the 30 km base created in the last step.

(f) Creating two circular areas. Create a text file with the center coordinate - 598060 6339670 5. Next import it into SAGA by Modules–>File–>Shapes–>Import–>import point cloud from text file and choose the seperator accordingly. Modules–>Shapes–>Tools–>Shape buffer, choose 2000 meter and the central point. Repeat again with 4000 meter.

(g) Interpolating. The interpolation is done with Modules–>Grid–>Calculus–>Grid generation–>function, inserting the grid limits, and the following equation: (sqrt((x-598060)^2+(y-6339670)^2)-2000)/(4000-2000). Next, using grid–>calculus–>grid calculator, using the bispline (initial grid), constant height grid and finally the interpolation function grid, with the following formula: g1*(1-g3)+g2*g3 (works in Saga 2.08)

(h) connecting it together. Shapes–>grid–>spatial extent–>clip grid with polygon. choosing the original grid and a 2000 meter radius (thus creating the inner circle)). Repeating the same with the grid that went the interpolation (it's name should end with "... calculation") and 4000 meter radius (thus creating the outer circle). Finally - merging: Grid–>Construction–>merging, choosing the inner circle, outer circle and the constant height grid. Make sure in the "overlapping cells" line you choose "first value in order of grid list".

(i) Finally - exporting to STL: Two steps - 1. Modules–>TIN–>conversions–>Grid to TIN, 2. Modules–>file–>shapes–>export–>export TIN to stereo litography file (STL)

The case: The following case includes a "Allrun" file, calling a python script twice, first for solving a mesh with about 500,000 cells, and then with 5,000,000 . The results are compared to the observations. They are not perfect, but comparing to the results Martinez 2011 got they are quite similar. With one discrepency though - downhill the recirculation is not predicted well by this mesh and discritization scheme (which is QUICKV for the div(phi,U)).

File:Askervein template.tar.gz (checked on OF 2.1.1)

Page under construction

@@ Line 1: / Line 1: @@
-The aim of this page is to hold validation cases for ABL (Atmospheric Boundary Layer) flows that are of interest to anybody starting to model in OF.
+The aim of this page is to hold validation cases for ABL (Atmospheric Boundary Layer) flows that are of interest to anybody starting to model ABL flows in OF.
-The page will hold links to databases and articles describing experiments, and also instructions on how to make meshes and set boundary conditions that give fair results, hopefully the end result will be complete cases.
+The simulations I perform, and refer to mainly in this text are RANS with k-epsilon turbulence closure models, using simpleFoam. Most of the first simulations that are placed in this page were done according to the excellent Martinez 2011 Ms.C. thesis [http://windenergyresearch.org/2011/09/wind-resource-in-complex-terrain-with-openfoam/] All of the boundary conditions and mesh creation suggestions were followed by from his work and the references within it.
-The simulations I perform, and refer to mainly in this text are RANS with k-epsilon turbulence closure models, using simpleFoam.
 ----
@@ Line 7: / Line 6: @@
 '''(0) Maintaining a logarithmic inlet profile over a 2D domain
 The wall function boundary condition implemented in OF is the Nikuradse's grain of sand roughness. It produces an anomaly in the k profile in the first cell above the ground which propagates upwards and down-flow. This is discussed extensively in the literature, and alternative wall functions have been suggested such as Richards and Hoexy 1993 (RH). There has lately been a very good documentation of testing this anomaly by Martinez in his Ms.C. thesis [http://windenergyresearch.org/2011/09/wind-resource-in-complex-terrain-with-openfoam/].
-Attached [[File:HomogenousTerrain.tar.gz]] is a case ("_BL2Dtemplate") and two python scripts using the pyFoam library to run different z0 values. To reproduce these graphs (p. 24 [http://windenergyresearch.org/2011/09/wind-resource-in-complex-terrain-with-openfoam/]) use "runZ0.py _BL2Dtemplate BL2Dvalidation_ 1 0.001 0.005 0.01 0.03 0.05 0.1"
+ [[File:BL2D.tar]] checked on OF 2.1
+Attached case ("BL2Dtemplate") has python scripts using the pyFoam library to run different z0 values (all files included in this page require the pyFoam library). run with ./Allrun .
-'''(1) Flow over Isolated 2D Hill
-Experiments by Khurshudyan et al.'''[http://cfd.mace.manchester.ac.uk/cgi-bin/cfddb/prpage.cgi?69&EXP&&database/cases/case69&cas69_head.html&cas69_desc.html&cas69_meth.html&cas69_data.html&cas69_refs.html&cas69_rsol.html&1&1&1&0&1&unknown]
+'''(1) Flow over Isolated 2D Hill (checked on OF 2.1)
-This is a wind tunnel experiment with a 2D hill of several aspect ratios, named also RUSHIL experiment. The latest comparison to it was published by Kasmi and Mason 2010 [http://onlinelibrary.wiley.com/doi/10.1002/we.390/full]
+The 2D simple hill from Martinez 2011 whose shape is given by
-'''Creating the STL surface'''
+<math>y\left(x\right) = -\frac{H}{6.04844}\left[J_0(\Lambda)I_0(\Lambda\frac{x}{a})-I_0(\Lambda)J_0(\Lambda\frac{x}{a})\right] \; \mathrm{for} \; |x-x_0| \leq a</math>
-The shape of the 2D hill is an analytical function described in [http://onlinelibrary.wiley.com/doi/10.1002/we.390/full] as:
+<math>y\left(x\right) = 0 \; \mathrm{for} \; |x-x_0| > a </math>
-<math>x = \frac{1}{2} \xi \left[ 1+\frac{a^2}{\xi^2+m^2(a^2- \xi^2)} \right]   </math> for <math> |x| \leq a</math>
+The mesh is created with blockMesh, and a template file, and also a python script that writes the hill shape as a polyline object in the blockMeshDict. The mesh grading was done by breaking the control box into 3 segments - upwind, hill area and downwind (unlike the Martinez mesh grading tanh scheme).
+The case with the python script, ready to run (./Allrun) is attached here
-<math>z = \frac{1}{2} m \sqrt{a^2-\xi^2} \cdot \left[ 1- \frac{a^2}{\xi^2+m^2(a^2- \xi^2)} \right] </math>
+ [[File:2DBump.tar]] (checked on OF 2.1)
-where <math> m = \frac{h}{a}+\sqrt{\frac{h}{a}+1}</math> and h is the height of the hill (<math>h = 0.117 [m]</math>) and a is the length of the hill. <math>\xi</math> is a parameter that changes from 0 to a. The aspect ratio of the hill is 3, 5 and 8.
+It compares the results with figures 21a,21b and 21c from Martinez 2011.
-The experimental setting was:
+----
+Another 2D case I tried simulating though unsuccessfully is mentioned here [[Sig_WindE_-_Validation_Cases_RUSHIL]].
+----
-logarithmic inlet profile with <math>z_0=0.157 \cdot 10^-3 [m]</math>, <math>u_{*} = 0.178 \left[ \frac{m}{s}\right]</math> which gives for instance <math>U_\infty = 3.9 \left[ \frac{m}{s} \right]  </math > at <math> z = 1 [m] </math>
+'''(2) Flow over Isolated 3D Hill (checked on OF 2.1.1)
-'''1''' The Profile was created with the desired discretization in a spreadsheet program. Column x (1st column) was the width of the hill (arbitrary width), 2nd column is x and the third z. Finally the 3 columns where exported as a csv file [[File:RUSHIL_8.csv]] (This should be according to the .xyz format).
+The same geometry as used by Martinez 2011, with the same recommendations. The mesh is created by snappyHexMesh. The STL file is created in salome with this python script [[File:3DBump_salome_script.tar‎]]. The case with an Allrun file is attached here
-Next, paraview is used to transform the csv into a STL surface, as explained  [http://www.cfd-online.com/Forums/openfoam/79032-movedynamicmesh.html in this thread], an reiterated here:
+  [[File:3DBump.tar]] (checked on OF 2.1.1)
-'''2''' The profile is uploaded in paraview.
+'''(3) Flow over real terrain - Askervein hill (checked on OF 2.1.1)
-'''3''' Open the csv in paraview using the csv reader, choose 1 column for each coordinate.
+The Askervein hill [http://www.yorku.ca/pat/research/Askervein/ASK82.pdf] and [http://www.yorku.ca/pat/research/Askervein/ASK83.pdf] is a low isolated eliptic hill located in the island of south Uist in Scotland. In 1982-3 a measurement campaign was conducted there, and it has since been used as a bench mark case for linear and non-linear (CFD) models for flow and pollutant dispersion simulations. A simulation case is presented here with use of the meshing recommandations from Martinez 2011. The terrain is massaged (smoothed) using SAGA GIS open source software, as will be described below.
-'''4''' Use the "TableToPoints' filter to obtain an array of points. The columns choice here is important so that the result will be a right hand side coordinate system. For the file above the order is '''y''' - '''x''' - '''z'''
+''Creating the STL surface''
-'''5''' Use the delaunay tool to "map" a suface from the point (The Delaunay 2D filter)
+In this case a real terrain is used. This means that the source of data is some DEM (Digital Elavation Model), and an additional step of smoothing the boundaries of the terrain has to be taken, in order to incorporate the boundary conditions used so far in the previous simulation, since they are appropriate only to flat terrain.
-'''6''' Save the data, you 'll be able to save it as an stl
+(a) For this simulation SRTM data was used, specifically this file [http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Eurasia/N57E008.hgt.zip]. It is probably not refined enough (90 meter resolution). In the OF7 workshop [http://www.cfd-online.com/Forums/openfoam/85128-wind-simulation-askervein-hill-2.html Carlos Peralta]
+pointed out to another source of data here [http://code.google.com/p/fds-smv/source/browse/trunk/FDS/?r=6478#FDS%2Ftrunk%2FValidation%2FAskervein_Hill%2FFDS_Input_Files]
+(b) The DEM was uploaded into SAGA [http://www.saga-gis.org/en/index.html] through the Modules–>File–>Grid–>Import–>Import USGS SRTM Grid (picture of the result [http://www.flickr.com/photos/hanan_levy/6817831074/])
-'''Creating the mesh'''
+(c) Dfining the coordinate system - Moduls–>Projection–>Set Coordinate Reference System , under EPSG code insert 4326 (WGS 84)
+(d) Projection - the lat/long information that the SRTM came in has to be projected unto a flat surface. I used the UTM projection (Universal Transverse Mercator). Modules–>Projection–>coordinate transformation grid . Under ESPG code chooce the one for the UTM zone that fits the Askervein hill. Finding this can be done through looking in UTM zone maps - [http://www.dmap.co.uk/utmworld.htm]. Ours is V29. relating ESPG codes can be done according to this table for instance [http://reference.mapinfo.com/common/docs/mapxtend-dev-web-none-eng/miaware/doc/guide/xmlapi/coordsys/systems.htm]. The ESPG code here is 32629
+Smooting - the area around the simulation zone is extended and smoothed, that is - transformed into a plane of constant height. This is done here with a 1/r^2 transformation, that is continous at the height but creates a discontinouity in the slope.
+(e) The area is cut and interpolated in one step: Modules–>Grid–>Gridding–>Spline Interpolation–>Multilevel B spline Interpolation (from grid) , with 5 km and 30 km grid extents, with a cell size of 50 meter. Next, a constant grid is created - the one the smoothing is done unto. A grid of height of 5 meter is chosen (according to the Askervein map and Martinez's experience) Modules–>Grid–>Construction–>create constant grid, with the 30 km base created in the last step.
+(f) Creating two circular areas. Create a text file with the center coordinate - 598060  6339670  5. Next import it into SAGA by Modules–>File–>Shapes–>Import–>import point cloud from text file and choose the seperator accordingly. Modules–>Shapes–>Tools–>Shape buffer, choose 2000 meter and the central point. Repeat again with 4000 meter.
+(g) Interpolating. The interpolation is done with Modules–>Grid–>Calculus–>Grid generation–>function, inserting the grid limits, and the following equation: (sqrt((x-598060)^2+(y-6339670)^2)-2000)/(4000-2000). Next, using grid–>calculus–>grid calculator, using the bispline (initial grid), constant height grid and finally the interpolation function grid, with the following formula: g1*(1-g3)+g2*g3 (works in Saga 2.08)
-----
+(h) connecting it together. Shapes–>grid–>spatial extent–>clip grid with polygon. choosing the original grid and a 2000 meter radius (thus creating the inner circle)). Repeating the same with the grid that went the interpolation (it's name should end with "... calculation") and 4000 meter radius (thus creating the outer circle). Finally - merging: Grid–>Construction–>merging, choosing the inner circle, outer circle and the constant height grid. Make sure in the "overlapping cells" line you choose "first value in order of grid list".
+(i) Finally - exporting to STL: Two steps - 1. Modules–>TIN–>conversions–>Grid to TIN, 2. Modules–>file–>shapes–>export–>export TIN to stereo litography file (STL)
+''The case'': The following case includes a "Allrun" file, calling a python script twice, first for solving a mesh with about 500,000 cells, and then with 5,000,000 . The results are compared to the observations. They are not perfect, but comparing to the results Martinez 2011 got they are quite similar. With one discrepency though - downhill the recirculation is not predicted well by this mesh and discritization scheme (which is QUICKV for the div(phi,U)).
+ [[File:Askervein_template.tar.gz]] (checked on OF 2.1.1)
+'''Page under construction
+[[Category:Incomplete pages]]