Under construction: Projected finish date: 5/18/2008 Sunday

CFD Simulation Scientific Paper (By: Jorge Rodriguez, Yong Sheng Khoo)

Title: Better Understanding of Flocculation Process through CFD Simulation

Abstract

Flocculation is an important process used by AguaClara to treat water. The process involves particle collisions and agglomeration to form flocs. Computational Fluid Dynamics was used to better understand the fluid dynamics in the reactor.  The standard K-ε model was used for every simulation model.  The pressure coefficient drop over one baffle turn is 3.75, which agrees with literature estimates.  After a clearance height of one baffle width or greater, the pressure coefficient drop and the maximum velocity become approximately constant.  Most of the energy dissipation occurs in the region after the turn over a distance of two baffle widths.  An area of flow recirculation occurs near the center wall immediately after the turn.  The pressure drop is not sensitive to the Reynolds number for a large range of inlet velocities.  Better understanding of the flocculation dynamics will enable optimized particle agglomeration and break-up. 

1. Introduction

Flocculation is the process by which particles collide and agglomerate. Past research has shown that shear gradients play an important role during flocculation. This process was simulated using Computational Fluid Dynamics (CFD).  The main task of this research is to find the optimum strain rate in the reactor to influence particle collision. Gambit and FLUENT were utilized to model one baffle turn. Gambit was used to create the geometry of the flocculator, and to generate the mesh.  FLUENT was used set up the boundary conditions and to obtain the results.

2. Methodology

The real life flocculation tank used by AguaClara involves 180 deg turns over few dozens baffles.  To save computational effort, a simple 180 deg turn over two baffles was modeled. The first step was to set up the geometry of the turn.  For future comparison with the experimental data, the design parameters for the pilot plant were used. The modeling approach was to create the geometry, mesh it, set boundary conditions, and solve it using FLUENT.

2.1 Creating Geometry 

Figure 1. Geometry of Flocculator

The design parameters used are:

Height: 1 m
Clearance: 0.15 m
Baffle width: 0.1 m
Velocity inlet: 0.1 m/s

With the geometry, the mesh for the model can then be set up.

2.2 Setting up Mesh 

Figure 2. Meshing Parameters (Click on the figure to see the original size)

 
Figure 2 shows the meshing parameters that were used. The boundary layers was first established at all the wall surfaces. The boundary layer was set such that the solution would provide a result of y+ less than 5. After that, the mesh edges were set up such that they will provide higher mesh resolution near the turn.  With the initial meshing conditions set up, all the faces were then meshed.

 
   
Figure 3. Mesh of the Model (Click on figure for original size)

Figure 3 shows the mesh of the flocculator model. As can be seen, the mesh is fine near the turn and at the walls.  The final step at this point was to set up the boundary conditions of the system.

2.3 Setting up Boundary Conditions 

Figure 4. Boundary Conditions

Figure 4 shows the boundary condition that was used for modeling. For a flocculator, there is an in flow and out flow of the fluids. Since inlet velocity inlet was known from the experimental data, the inlet was set to the Velocity Inlet type boundary condition. The outlet was set to Pressure Outlet boundary condition type, the atmospheric pressure.

The mesh was then saved and exported to FLUENT for further obtaining solution and further analysis.  

2.4 Solve using FLUENT 

At this stage, the Standard k-ε turbulence model was set up.  Water was defined as the working material from the FLUENT database. The discretization method for the momentum, turbulent kinetic energy, and turbulent dissipation rate were set to the 'Second Order Upwind' scheme to obtain a 'Second Order Accurate' solution.  The boundary conditions were set according to the values shown in the table below.  The solution was obtained by iterating until the residuals converged to 10e-6.  Results were then analyzed and plotted.        

TABLE IN BCS TABLE 1

2.5 Mesh Sensitivity Analysis

The effect of the number of mesh elements on the result was also carried out. Coarse, medium and fine meshes were created and the pressure drop across the turn from each mesh was compared. This analysis will provide confidence on the accuracy of certain mesh. If the changes in mesh elements does not result in a lot of change in pressure coefficient drop, it is concluded that the mesh elements were refined enough that the truncation and discretization errors can be neglected. Table 2 shows the summary of 3 meshes created for mesh sensitivity analysis. Please refer back to figure 2 for corresponding meshing parameters.

Table 2. Mesh Meshing Parameters for Coarse, Medium and Fine Meshes

2.6 Effect of Reynolds Number

Since the inlet flow rate can vary substantially across AguaClara plants the effects of inlet Reynolds number on pressure coefficient drop were also examined.

2.7 Parameterization

At the later stage of project, after the confident on result of the model was built up, the effect of geometry on results was analyzed. Different clearance heights were used for analyzing pressure drops and maximum velocities. It would be tedious to individually recreate each geometry and mesh for different clearance height from scratch. For this reason, a parameterization technique was used. The original Gambit journal file was modified to include the variable clearance height. Using this method, changes in corresponding clearance height were plugged into the journal file to obtain the desired mesh/geometry. The journal file is included in the Appendix.

2.8 Comparing Turbulence Model

The pressure coefficient results were compared using Standard K-ε, K-ε Realizable and K-ω turbulence models.

3. Results and Discussions

The results considered were plots of the velocity vectors, pressure coefficient contours, contours of strain rate and contours of turbulence dissipation rate.

Figure 5. Velocity Vectors (Click on figure for original size)
 
Velocity vector plot shows the velocity of the fluids throughout the flocculator.  As can be seen, there is a region of high velocity at the outer turn and recirculation at the inner turn. At the bottom of the flocculator, there is region of stagnant fluid.
 

Figure 6. Contours of Stream Function (Click on figure for original size)
 
Contours of stream function tell us how particles of fluid travel in the flocculator. As can be seen from figure 6, there is an enclosed streamline at the inner turn. This means there is recirculating fluid 'trapped' in that region.

 
 

Figure 7. Contours of Pressure Coefficient

Figure 7 shows most of the pressure coefficient drop occurs around the bend. There is a pressure coefficient drop of about 3.7 across the bend. This is in excellent agreement with literature estimates.

Figure 8. Contours of Strain Rate

Contours of strain rate show high velocity gradients around the turn. There are also high strain rates in the boundary layer near the wall.  It is postulated that flocculation is directly proportional to the strain rate.

Figure 9. Contours of Turbulent Disssipation Rate

Contours of turbulent dissipation rate show a similar trend as the contours the strain rate right after the turn. The region of highest turbulence dissipation occurs after the turn.  As can be seen from figure 9, this region is roughly twice the length of baffle spacing. 

Figure 10. Wall Yplus

Figure 10 shows that the yplus values were consistently less than 5.  According to FLUENT documentation "the mesh should be made either coarse or fine enough to prevent the wall-adjacent cells from being placed in the buffer layer (yplus = 5~30)". Since the yplus from the model was consistently less than five (inside the viscous sublayer) the turbulence flow near the walls was resolved properly.

Figure 11. Mesh Sensitivity Analysis

Figure 11 shows the pressure coefficient drop over one turn for different mesh densities.  In general finer meshes provide more accurate results. However, as the mesh was refined the pressure drop remained constant as can be seen in figure 11.  Hence, it was concluded that results were not sensitive to mesh density and the coarse mesh was sufficient.

Figure 12. Reynolds Number Effect on Pressure Coefficient Drop

By changing the inlet velocity a range of Reynolds numbers were analyzed. The normal flow rate has a Reynolds number of 10,000. From figure 12, it is seen that the pressure coefficient drop changes only a little with big changes in Reynolds number. In other words, the pressure coefficient drop is not sensitive to the Reynolds number at the inlet. This is a good thing because in the design of the flocculator, the inlet flow rate can be neglected. This is to say that one flocculator design can be used for different flow rates.
 
Figure 13. Clearance Height Effect on Pressure Coefficient Drop and Maximum Velocity

By adjusting clearance height, the effect on the pressure coefficient drop was also analyzed. It can be seen that the pressure coefficient drop is independent of the change in clearance height as long as the clearance height is greater than a certain critical value. Figure 13 shows that after a critical value of 1, the pressure coefficient drop is constant. This phenomena can be explained using figure 5. Figure 5 shows the clearance height of 0.15 m. However, from 0.1 m onward, the flow is mostly stagnant in the flocculator. This mean that 0.1 m is needed for the flow to navigate through the turn and after this point onward, there is not much activity happening. Clearance height of less than 0.1 m gave higher pressure coefficient drop as it created a constriction of flow and the frictional loss was increased. With this result, it is recommended for the design team that the clearance height must be at least the same of bigger than the baffle width to produce predictable pressure coefficient drop.

Figure 14. Comparison of Turbulent Dissipation Rate for Clearance height of 0.1 m and 0.15 m

To further validate that the result is not sensitive to the change in clearance height, the contours of turbulence dissipation rate of clearance height 0.1 m and 0.15 m was compared. The result showed that the region of active turbulent dissipation was the same, about two times the length of baffle spacing. With this result, it is concluded that the design team has the freedom of choosing clearance height according to their design constraint and not theoretical constraint as long as the clearance height is greater than the baffle spacing.

Figure 14. Effect of Turbulence Model on Pressure Coefficient Drop

Pressure coefficient drop was least in the K-ε model and the most in the K-ω model. 
 

Figure 15. Contours of Velocity Magnitude for Different Turbulence Model

K-ε has smallest high and low velocity region. This explain why it has the lowest pressure coefficient drop. K-ω has the biggest region of high and low velocity. K-ω has highest pressure coefficient drop.
 
 

4. Conclusions

1. An area of recirculation occurs near the center wall immediately after the turn.

2. Pressure coefficient drop over one baffle turn is 3.75.

3. After a clearance height of one baffle width or greater, the pressure coefficient drop and the maximum velocity becomes constant. 

4. Most of the energy dissipation occurs in the region after the turn over a distance of about two baffle width.

5. The pressure drop is insensitive to the Reynolds number for a large range of inlet velocities.  

6. About 20,000 mesh elements is sufficient to obtain accurate modeling results.

7. Different turbulence model resulted in fairly different results.

5. Acknowledgements

During the course of this project, we received invaluable advice and direction from Dr.  Monroe Weber-Shirk, whose direction of the AguaClara team project spurred this research and technical investigation.  All questions we had were directed towards Dr. Rajesh Bhaskaran, his expertise and academic guidance kept us on schedule and on task.  We would also like to thank Prof. Brian Kirby for elucidating some technical aspects we encountered.  For their help, we are grateful. 

Appendix A

/ Journal File for GAMBIT 2.4.6, Database 2.4.4, ntx86 SP2007051421
/ Identifier "Clearance 1.5W"
/ File opened for write Fri Apr 18 10:07:42 2008.
undo begingroup

$clearance = 0.1

coordinate modify "c_sys.1" xyplane xaxis add 0 AND 0.1 AND 0.2 reset snap \
  lines
coordinate modify "c_sys.1" xyplane yaxis reset snap lines
window modify coordinates "c_sys.1" xyplane grid
undo endgroup
undo begingroup
coordinate modify "c_sys.1" xyplane xaxis add 0 AND 0.1 AND 0.2 reset snap \
  lines
coordinate modify "c_sys.1" xyplane yaxis add 0 AND 0.05 AND 0.1 AND 0.15 AND \
  0.2 AND 0.25 AND 0.3 AND 0.35 AND 0.4 AND 0.45 AND 0.5 AND 0.55 AND 0.6 AND \
  0.65 AND 0.7 AND 0.75 AND 0.8 AND 0.85 AND 0.9 AND 0.95 AND 1 reset snap \
  lines
window modify coordinates "c_sys.1" xyplane grid
undo endgroup
vertex create coordinates 0 0 0
vertex create coordinates 0.1 0 0
vertex create coordinates 0.2 0 0
vertex create coordinates 0.2 $clearance 0
vertex create coordinates 0.1 $clearance 0
vertex create coordinates 0 $clearance 0
vertex create coordinates 0.2 1 0
vertex create coordinates 0.1 1 0
vertex create coordinates 0 1 0
edge create straight "vertex.4" "vertex.7"
edge create straight "vertex.5" "vertex.8"
edge create straight "vertex.6" "vertex.9"
edge create straight "vertex.8" "vertex.9"
edge create straight "vertex.7" "vertex.8"
edge create straight "vertex.3" "vertex.4"
edge create straight "vertex.2" "vertex.5"
edge create straight "vertex.1" "vertex.6"
edge create straight "vertex.5" "vertex.6"
edge create straight "vertex.4" "vertex.5"
edge create straight "vertex.2" "vertex.3"
edge create straight "vertex.1" "vertex.2"
undo begingroup
coordinate modify "c_sys.1" xyplane xaxis add 0 AND 0.1 AND 0.2 reset snap \
  lines
coordinate modify "c_sys.1" xyplane yaxis add 0 AND 0.05 AND 0.1 AND 0.15 AND \
  0.2 AND 0.25 AND 0.3 AND 0.35 AND 0.4 AND 0.45 AND 0.5 AND 0.55 AND 0.6 AND \
  0.65 AND 0.7 AND 0.75 AND 0.8 AND 0.85 AND 0.9 AND 0.95 AND 1 reset snap \
  lines
window modify coordinates "c_sys.1" xyplane nogrid
undo endgroup
face create wireframe "edge.4" "edge.2" "edge.9" "edge.3" real
face create wireframe "edge.5" "edge.1" "edge.10" "edge.2" real
face create wireframe "edge.6" "edge.11" "edge.7" "edge.10" real
face create wireframe "edge.12" "edge.7" "edge.9" "edge.8" real
undo begingroup
/ERROR occurred in the next command!
blayer create first 0 growth 1.2 total 0 rows 4 transition 1 trows 0 uniform
undo endgroup
undo begingroup
blayer create first 0.0003 growth 1.25 total 0.0077407 rows 9 transition 1 \
  trows 0 uniform
blayer attach "b_layer.1" face "face.2" "face.3" "face.3" "face.4" "face.1" \
  "face.2" "face.1" "face.4" "face.4" "face.3" edge "edge.1" "edge.6" \
  "edge.7" "edge.7" "edge.2" "edge.2" "edge.3" "edge.8" "edge.12" "edge.11" \
  add
undo endgroup
undo begingroup
edge picklink "edge.3" "edge.2" "edge.1"
edge mesh "edge.1" "edge.2" "edge.3" successive ratio1 1.01 size 0.005
undo endgroup
undo begingroup
edge picklink "edge.6" "edge.7" "edge.8"
edge mesh "edge.8" "edge.7" "edge.6" successive ratio1 1 size 0.002
undo endgroup
undo begingroup
edge picklink "edge.4" "edge.5" "edge.10" "edge.9" "edge.11" "edge.12"
edge mesh "edge.12" "edge.11" "edge.9" "edge.10" "edge.5" "edge.4" successive \
  ratio1 1 size 0.002
undo endgroup
face mesh "face.2" map size 1
face mesh "face.1" map size 1
face mesh "face.3" map size 1
face mesh "face.4" map size 1
physics create "Inlet" btype "VELOCITY_INLET" edge "edge.4"
physics create "Outlet" btype "PRESSURE_OUTLET" edge "edge.5"
physics create "Left_Wall" btype "WALL" edge "edge.3" "edge.8"
physics create "Middle_Wall" btype "WALL" edge "edge.2"
physics create "Right_Wall" btype "WALL" edge "edge.1" "edge.6"
physics create "Bottom_Wall" btype "WALL" edge "edge.12" "edge.11"
export fluent5 "Clearance 1.5W.msh" nozval
save