The Turbulent Pipe Flow Example - High-Re
The purpose of the Cornell tutorials with OpenFOAM® is to give a few easy examples which are well known for the CFD users. The original tutorials were written by Rajesh Bhaskaran. They give the user a practical guide on software usage (FLUENT / ANSYS) and the bacground is also explained well. Therefore it is recommended to read the particular Cornell tutorial first. In my tutorials I show the same steps using OpenFOAM-1.7 and other open software so that the reader can produce comparable results. I made these calculations on a Linux operating system (Ubuntu 10.04) using OpenFOAM-1.7, Octave, gnuplot and xmgrace.
The original tutorial can be found here.
The case will be calculated with a High-Re turbulence model. An other turbulent pipe flow tutorial calculates the same example with a Low-Re model. One may find information on the two calculation types here and here.
A similar mesh can be used as in the laminar case. The only difference is that we will use a higher resolution in radial direction (10 cells). First copy the full laminar case:
cp -rv <laminar_case_directory> <turbulent_case_directory>
Change the blockMeshDict file:
Run the meshing again:
The mesher will create a default patch with no elements in it. Delete this patch and change the number of patches to 5 in the constant/polyMesh/boundary file:
The boundary and initial conditions can be set within the field files in the directory of the first iteration/timestep (counting time/iteration starts with 0). This time we will need the turbulent field variables as well, so they should be copied new from the tutorials if they have been deleted in the laminar flow tutorial. Copy the k field first:
cp -rv $FOAM_TUTORIALS/incompressible/simpleFoam/pitzDaily/0/k <turbulent_case_directory>/0/
Repeat this for the epsilon and nut fields also.
When using a High-Re model a wall function should be used. The main idea is that the boundary layer is not resolved by the mesh (we do not have to use a really fine mesh near the wall), but the boundary layer is calculated by a semi-empirical formula in one cell: the boundary layer is in the first cell. The type of wall function to be used can be defined in the nut field file. First We will use the so called standard wall function which is called nutWallFunction in OpenFOAM®.
Additionally we have to tell the solver that we use wall functions for k and epsilon.
Initialize the velocity field to the inlet value (according to the shown turbulence fields).
This is similar case to the simpleFoam pitzDaily tutorials where wall functions are used also.
These values should be put to the inlet value and to the value at the wall functions also.
Summary of all the boundary conditions:
We shall change now the value of viscosity in the transportProperties file (located in the constant directory):
Set the followings in the file (keep the fluid type Newtonian):
We will use the k-epsilon turbulence model (called kEpsilon in OpenFOAM®, see available turbulence models here) therefore the RASProperties file should be modified:
fvSchemes: keep it unchanged.
Change only the endTime and writeInterval in the controlDict file (both to 250):
Launch the solver the same way as in the laminar flow example:
user@machine:~$ simpleFoam > log &
Monitor the log file:
user@machine:~$ tail -200f log
A similar gnuplot script can be used as in the laminar flow example, but it should be changed to monitor the k and epsilon residuals as well:
Residual history for this calculation:
Download the Residuals-turb gnuplot script.
Start with checking the axial velocity profile. Exporting the profile using the sample utility (use the sampleDict1 file from the laminar flow example files) one will get the following (as datasets were written to the format of Grace, it is enough to double-click on the axial_profile_U.agr in the file manager e.g. in nautilus when using Ubuntu):
We can not be sure if that axial velocity is converged until the end of the pipe. Therefore launch the case on a longer mesh. Lets change the mesh length to 16m. Follow these steps:
This is the new axial velocity plot:
We can see that the axial velocity finally converged to the value of 1.2 m/s. The comparable Fluent result can be checked here. The converged axial velocity obtained with ANSYS 12 is 1.19 m/s (according to the chart shown).
We should check the height of the near wall cells now. The non-dimensional measure for this is the y+ (definition @ cfd-online). Using a High-Re turbulence model y+ should be ~=30. The yPlus field is not written by default, we should generate it using the yPlusRAS utility:
user@machine:~$ yPlusRAS -latestTime
We will need the wallShearStress for checking the skin friction coefficient (definition @ cfd-online), therefore create it now:
user@machine:~$ wallShearStress -latestTime
One may find the definition of skin friction coefficient
The following new files are created:
Modify the sampleDict2 file from the laminar example files to have also the yPlus sampled::
Load yPlus to Grace:
This is finally the yPlus plot:
We can see that y+ is too low. Therefore we should recalculate the case with a lower radial mesh resolution. Follow these steps:
Loading the yPlus plot again:
y+ is closer to 30 with the rough mesh!
Create the skin friction coefficient now from the wallShearStress plot:
This is what you should get:
The value at the developed flow range (end of pipe) is ~0.00734 for the both meshes. The corresponding Fluent result is close to 0.01 (according to the graph shown).
The radial profile is already created together with the axial profile when using the sampleDict1 file with the sample utility. As sample created the output the default formating is: velocity on the y axis and distance on the x axis. We can change this to have an easier comparison with the Fluent results:
After some formatting (and loading results obtained on both meshes):
Comparing it with the Fluent radial profile we can see that the max velocity is better captured on the rough mesh (200x6). Please note that the Fluent calcultion was done with a Low-Re model, therefore the distribution is finer due to the higher number of points displayed.
Version: 1.7, 2011-11-06