Week 5 - Rayleigh Taylor Instability

AIM

     To conduct the Rayleigh Taylor CFD simulation.

 

OBJECTIVE

• Explain different CFD models used for solving Rayleigh Taylor Instability problems
• Perform a Rayleigh Taylor instability analysis using the Volume of Fluid (VOF) multiphase model for air and water for two different mesh element sizes.
• Repeat the process for water and user-defined material with density = 400 kg/m3, viscosity = 0.001 kg/m-s, and compare the behaviour
• Define the Atwood number and derive it for each case. Explain the effect of Atwood number on instability.
• So there will be total of 3 cases

1. Case 1 : Mesh size- 0.5mm, Water and air
2. Case 2 : Mesh size- 0.1mm, Water and air
3. Case 3 : Mesh size- 0.1mm, Water and user defined fluid

 

INTRODUCTION

            In this project, we are going to perform the two dimensional simulation for two immiscible fluids that have different densities. Initially, the fluid with high density is placed at the top region while the fluid having lesser density at the bottom. Due to gravity, the heavy dense fluid flows downwards making the less dense fluid to go up. This causes the Rayleigh Taylor instability between the interface of those two fluids. This instability is captured from this simulation which set to run at transient state.

 

THEORY

 

Rayleigh–Taylor instability

            The Rayleigh–Taylor instability, or RT instability is an instability of an interface between two fluids of different densities which occurs when the lighter fluid is pushing the heavier fluid. Examples include the behaviour of water suspended above oil in the gravity of Earth, mushroom clouds like those from volcanic eruptions and atmospheric nuclear explosions, supernova explosions in which expanding core gas is accelerated into denser shell gas, instabilities in plasma fusion reactors and inertial confinement fusion.

              As the RT instability develops, the initial perturbations progress from a linear growth phase into a non-linear growth phase, eventually developing "plumes" flowing upwards (in the gravitational buoyancy sense) and "spikes" falling downwards. In the linear phase, the fluid movement can be closely approximated by linear equations, and the amplitude of perturbations is growing exponentially with time. In the non-linear phase, perturbation amplitude is too large for a linear approximation, and non-linear equations are required to describe fluid motions. In general, the density disparity between the fluids determines the structure of the subsequent non-linear RT instability flows (assuming other variables such as surface tension and viscosity are negligible here). The difference in the fluid densities divided by their sum is defined as the Atwood number, A. For A close to 0, RT instability flows take the form of symmetric "fingers" of fluid; for A close to 1, the much lighter fluid "below" the heavier fluid takes the form of larger bubble-like plumes.

            This process is evident not only in many terrestrial examples, from salt domes to weather inversions, but also in astrophysics and electrohydrodynamic. For example, RT instability structure is evident in the Crab Nebula, in which the expanding pulsar wind nebula powered by the Crab pulsar is sweeping up ejected material from the supernova explosion 1000 years ago. The RT instability has also recently been discovered in the Sun's outer atmosphere, or solar corona, when a relatively dense solar prominence overlies a less dense plasma bubble. This latter case resembles magnetically modulated RT instabilities.

 

Stages of development and eventual evolution into turbulent mixing

                       The evolution of the RTI follows four main stages. In the first stage, the perturbation amplitudes are small when compared to their wavelengths, the equations of motion can be linearized, resulting in exponential instability growth. In the early portion of this stage, a sinusoidal initial perturbation retains its sinusoidal shape. However, after the end of this first stage, when non-linear effects begin to appear, one observes the beginnings of the formation of the ubiquitous mushroom-shaped spikes (fluid structures of heavy fluid growing into light fluid) and bubbles (fluid structures of light fluid growing into heavy fluid). The growth of the mushroom structures continues in the second stage and can be modelled using buoyancy drag models, resulting in a growth rate that is approximately constant in time. At this point, nonlinear terms in the equations of motion can no longer be ignored. The spikes and bubbles then begin to interact with one another in the third stage. Bubble merging takes place, where the nonlinear interaction of mode coupling acts to combine smaller spikes and bubbles to produce larger ones. Also, bubble competition takes places, where spikes and bubbles of smaller wavelength that have become saturated are enveloped by larger ones that have not yet saturated. This eventually develops into a region of turbulent mixing, which is the fourth and final stage in the evolution. It is generally assumed that the mixing region that finally develops is self-similar and turbulent, provided that the Reynolds number is sufficiently large.

 

Practical CFD models that have been based on mathematical analysis of RT waves

1. Kelvin-Helmholtz instability
2. Richtmyer-Meshkov instability
3. Plateau-Rayleigh instability

 

Kelvin-Helmholtz instability

            This instability typically occurs when there is velocity shear in a single continuous fluid, or additionally where there is a velocity difference across the interface between two fluids. A common example is seen when wind blows across a water surface; the instability constant is able to manifest itself through waves on water surface. The Kelvin-Helmholtz instability is not only restricted to a water surface as clouds, but is evident through other natural phenomena as ocean and the sun's corona. The theory predicts the onset of instability and transition to turbulent flow in fluids of different densities moving at various speeds. This instability occurs where there is a velocity shear or velocity difference between the interface of the two fluids.

Richtmyer-Meshkov instability

                        This occurs when two fluids of different density are impulsively accelerated. Normally this is by the passage of a shock wave. The development of the instability begins with small amplitude perturbations which initially grow linearly with time. This is followed by a nonlinear regime with bubbles appearing in the case of a light fluid penetrating a heavy fluid, and with spikes appearing in the case of a heavy fluid penetrating a light fluid. A chaotic regime eventually is reached and the two fluids mix. This instability can be considered the impulsive-acceleration limit of the Rayleigh Taylor instability.

Plateau-Rayleigh instability

            This instability, often just called the Rayleigh instability, explains why and how a falling stream of fluid breaks up into smaller packets with the same volume but less surface area. It is related to the RTI and is part of a greater branch of fluid dynamics concerned with fluid thread breakup. This fluid instability is exploited in the design of a particular type of ink jet technology whereby a jet of liquid is perturbed into a steady stream of droplets. The driving force of the Plateau–Rayleigh instability is that liquids, by virtue of their surface tension, tend to minimize their surface area.

Atwood number

            The Atwood number (A) is a dimensionless number in fluid dynamics used in the study of hydrodynamic instabilities in density stratified flows. It is a dimensionless density ratio defined as

where

ρ₁= density of heavier fluid

ρ₂= density of lighter fluid

 

SOLVING & MODELLING APPROACH

1. Construct a 2D geometry for the 2 fluid phases.
2. Mesh independence study is carried out for different mesh length.
3. Laminar Transient state with volume of fluid mode is used
4. The simulation is calculated for very small time steps to avoid floating point error.
5. Animations are created to analyse the flow patterns

                                                   Total Cases.

 

PRE PROCESSING AND SOLVER SETTING

            In our challenge we will create a flow simulation of Rayleigh Taylor Instability. In general I will be explaining only one case and posting the screenshots of the other cases.

 

Case 1

                        Mesh size        - 0.5mm,

                        Fluid                - Water and air

 

            Step 1 : Open the geometry Space Claim. Two square with 20x20mm is drawn one over the other for the two fluids. Once this is done the topology is shared  so the meshing will be smooth. 

Once this process is done. We will close the space claim and open the mesh module.

 

            Step 2 : Open the mesh module and under the mesh details give CFD Fluent as preference and Element size of 0.5mm. But since we are going to go for mesh independence study uniform meshing size is maintained throughout the volume. Note the number of nodes and elements for the study.

                        After the meshing is done check for the quality criterion. Check whether the mesh quality is above 5%. Once this is done name the faces of the volume as Water and air.

 

 

Step 3 : After the meshing and face naming are done move on to Fluent Solver. In the fluent launcher select double precision, display mesh after reading and give the appropriate solver processors and GPUs.

In the ANSYS CFD we need to give all the conditions, parameters and models. To start with go to physics menu and click general settings. In that select Pressure based type solver, Absolute velocity formulation and transient state flow. Next we need to give the model. For this case Laminar is selected.

Next is the flow material selection. By default air is selected. We have two fluid phases. Hence select water from the fluent database.

Click on the multiphase and select volume of fluids with two fluids. In the phases tab rename and assign the primary and secondary phases.

 

After the physics part is done move to the solution part. Initialise the simulation first. click Initialize to initialize the boundary conditions. Hybrid method is selected before initialization. Click on Autosave to obtain a animation of the flow in the post results.

 

Patch the created region with the lubricating fluid. Set the fluid and select volume fraction. Set the values to 1 and 0. 

      

Then run the calculation for a given number of iterations till the convergence is obtained. Depending upon the number of elements and model selected the time required for convergence varies. In the CFD module itself we can compute, measure, plot, animate, etc., if needed.

 

Step 4 : After the solutions and calculations move to the result module to get different graphs, plots, contours, animations etc.,. Sectional views can be created if required.

 

 

RESULTS

 

Case 1

                              Mesh size        - 0.5mm,

                              Fluid                - Water and air

                                                                        Residuals

                                                                    Volume Fraction

                                                                   

 

                                                Animation

 

Case 2

                              Mesh size        - 0.1mm,

                              Fluid                - Water and air

 

                                                                         Residuals

                                                                     Volume Fraction

 

 

                                                 Animation

 

Case 3

                              Mesh size        - 0.1mm,

                              Fluid                - Water and User defined fluid

 

                                                                       Residuals

 

                                                                    Volume Fraction

 

 

                                                    Animation

 

 

Atwood Number

It is a dimensionless number that is used to study hydrodynamic instabilities in density stratified flows.

Where,

ρ1 is density of heavier fluid 

ρ2 is density of lighter fluid

For Air-water RT instability,

ρ1 = 998.2 kg/m³ [Water density]

ρ2 = 1.225 kg/m³ [Air density]

 

         A = 0.998

 

For water- User Defined material RT instability,

ρ1 =  998.2 kg/m³ [Water density]

ρ2 =  400 kg/m³ [UD material density]

 

        A = 0.428

 

 

                            RTI evolves in three distinct stages: Linear stability , mushroom head(falling) or bubble(rising) formation , long term evolution due to bubble merging and mixing. For low Atwood number, bubble and mushroom head shape is relatively symmetrical, the disturbances or perturbations are linear for a larger time. In contrast, for large Atwood number, bubble and mushroom head shape is less symmetrical, also it occurs very rapidly at an accelerating rate which leads to a nonlinear growth rate.

                            Finer mesh captures the instability better than the coarse mesh , but it takes more time for the heavier fluid to settle down. So it is computationally expensive but produces results with better accuracy.

                            Since the simulation involves capturing instabilities which are changing in an unpredictable manner as a function of time , transient approach is preferred over steady state approach. 

 

Why Steady State simulation is not suitable ?

• We know that steady state simulation is a time independent approach. 
• Based on the simulation we know there is high turbulence involve and the change in the flow occurs between very small time variation. This is not possible to capture in steady state solution
• Also steady state solver ignore the time differential terms involved in the simulation and thus even if it may be able to capture the state of equilibrium, it skips to capture the non-linearities in the simulation precisely.

 

 

 

CONCLUSION

• From the above cases, we observe that finer mesh captures the instability better than the coarse mesh, but it takes more time for the heavier fluid to settle down. As computational time is long but produces results with better accuracy.
• Since the simulation involves capturing instabilities that are changing in an unpredictable manner as a function of time, the transient approach is preferred over the steady-state approach. 

• For Atwood number close to 0, RT instability flows take the form of symmetric fingers of fluid.
• when Atwood number close to 1, the much lighter fluid below the heavier fluid takes the form of larger bubble-like plumes.