Rayleigh-Taylor Instability Analysis using VOF Multi-Phase Method

Aim: To perform Rayleigh Taylor Instability simulation for 3 different mesh Sizes.

Theory:

Rayleigh Taylor Instability is the instability at the interface between two fluids having different densities in hydrostatic equilibrium under acceleration. The lighter fluid pushes the heavier fluid which would result in the formation of spikes. These spikes eventually develops in a non-linear fashion and changes into plumes like structures. Spike like structures are also developed by the denser fluid as it falls downward which also grows non-linearly.

In the linear phase, equation can be linearized and the amplitude of pertubations is growing exponentially with time. In the non-linear phase, pertubation amplitude is too large for the non-linear terms to be neglected.

Some of the example for Rayleigh Taylor instability:

• Behavior of water suspended above oil against the earth\'s gravity.
• Supernova explosions in an expanding core gas is accelerated into denser shell.
• Instabilites in Plasma fusion reactors and inertial confinement fusion.
• Mushroom cloud formation from an volcanic eruptions or nuclear explosion.

Some practical CFD models that have been based on the mathematical analysis of Rayleigh Taylor waves are:

• Kelvin-Helmholtz or Rayleigh Taylor Hybrid Model : This model is used to predict primary break up of the intact liquid core of the jet and the Rayleigh Taylor model determines the effects of acceleration of different phases. This model predicts the shape of the spray taking into consideration of temperature at different phases.
• Granular Rayleigh Taylor Instability: A granular instability driven by gravity is studied both numerically and experimentally. The instability occurs as grains fall in a closed Hele-Shaw cell where a layer of dense granular material is positioned above a layer air. The initially flat front defined by the grains subsequently develops into a pattern of falling granular fingers seperated by rising bubbles of air. A transient coarsening of the front is observed right from the start by a finger merging process. The coarsening is later stabilized by new fingers growing from the center of the rising bubbles.
• Potential Flow Models: It is a well known framework for the study of surface modes and has been the basis of many classical stability results. This model is applied to two dimensional single mode evolution as well as two bubble competation,for both the Rayleigh-Taylor and Richtmyer-Meshkov Instabilities, latter treated in an impulse approximation.This model is also applied to the evolution of three dimensional modes with an initial rectangular geometry. The model yeilds the aspect ratio dependence of the early non-linear stages, in agreement with third order pertubation theory.
• Lattice Boltzmann Model: This is developed for two dimensional combustion simulation which is an efficient numerical method. They are also suitable for incompressible flows.
• K L Turbulence Model: It is used to describe self similar growth of Rayleigh-Taylor and Richtmyer-Meshkov instabilies. The model describes the dominant eddies in the mixing zone with evolutionary equations for their characteristic dimension (L) and energy per unit mass. These equations are based on buoyancy-drag models for RT and RM flows.

We will be going through the common steps used in the CFD analysis.

• Geometry Preparation:

In this case, we will be considering a two dimensional domain whose dimensions are 20 x 40 mm where initially the lower domain is filled with air and the upper domain is filled with water seperated by an interface line. Share topology option was used to define the inetrface.

Meshing :

Case 1:

Number of elements : 3200

Element Size: 0.0005 m

The image of the mesh for the baseline case can be seen below:

Setup :

We will be using a pressure based solver and will be proceeding with transient simulation. We will be enabling the gravity and will be entering the value of -9.81 m/s as our water is in the upper domain. Then we will be selecting Multi-phase method as we have 2 different fluids. Once we select the multi-phase method we will be using a volume of fluid method. The VOF model is a surface tracking technique applied to a fixed Eulerian mesh. Then we select the material used for the fluid which would be the standard air and water from the ansys material libary.

Once these has been set up we move on to solve by initializing the solution using standard initialization. After standard initialization, patch option is used to allocate the fluid domain with their respective fluids (i.e. water is allocated to the upper domain and air is allocated to the lower domain). In the patch volume fraction of waterin the upper domain is set to 1 and in the lower domain it is set to 0. Once the patching of the fluids is done we can check this using the contour of volume fraction of water as shown below:

Results:

Residual Plots:

Rayleigh-Taylor Contour:

Volume fraction of water contour animation:

Time step: 0.005

Case 2:

Number of elements : 20000

Element Size: 0.0002 m

The image of the mesh for the baseline case can be seen below:

Setup :

We will be using a pressure based solver and will be proceeding with transient simulation. We will be enabling the gravity and will be entering the value of -9.81 m/s as our water is in the upper domain. Then we will be selecting Multi-phase method as we have 2 different fluids. Once we select the multi-phase method we will be using a volume of fluid method. The VOF model is a surface tracking technique applied to a fixed Eulerian mesh. Then we select the material used for the fluid which would be the standard air and water from the ansys material libary.

Once these has been set up we move on to solve by initializing the solution using standard initialization. After standard initialization, patch option is used to allocate the fluid domain with their respective fluids (i.e. water is allocated to the upper domain and air is allocated to the lower domain). In the patch volume fraction of waterin the upper domain is set to 1 and in the lower domain it is set to 0. Once the patching of the fluids is done we can check this using the contour of volume fraction of water as shown below:

Results:

Residual Plots:

Rayleigh-Taylor Contour:

Volume fraction of water contour animation:

Time step: 0.002

Case 3:

Number of elements : 80000

Element Size: 0.0001 m

The image of the mesh for the baseline case can be seen below:

Setup :

We will be using a pressure based solver and will be proceeding with transient simulation. We will be enabling the gravity and will be entering the value of -9.81 m/s as our water is in the upper domain. Then we will be selecting Multi-phase method as we have 2 different fluids. Once we select the multi-phase method we will be using a volume of fluid method. The VOF model is a surface tracking technique applied to a fixed Eulerian mesh. Then we select the material used for the fluid which would be the standard air and water from the ansys material libary.

Once these has been set up we move on to solve by initializing the solution using standard initialization. After standard initialization, patch option is used to allocate the fluid domain with their respective fluids (i.e. water is allocated to the upper domain and air is allocated to the lower domain). In the patch volume fraction of waterin the upper domain is set to 1 and in the lower domain it is set to 0. Once the patching of the fluids is done we can check this using the contour of volume fraction of water as shown below:

Results:

Residual Plots:

Rayleigh-Taylor Contour:

Volume fraction of water contour animation:

Time step: 0.001

Discussions:

For all the above cases, Rayleigh-Taylor Instability can be observed near the interface of the two fluids due to their difference in densities. Once the mesh gets finer the instability can be observed in detailed but takes a longer computing time as the numebr of elements is increased. Initially the lighter fluid tends to move up in a bulk manner (bubbles) as heavier fluid passes through them. The difference in the densities between the two fluids is large(i.e. Atwood Number is large), so a appropriate turbulence model should be used to observe the behaviour more accurately.

Why a steady state approach might not be suitable for these types of simulation.

• Transient solvers tends to capture the flows better compared to the steady state solvers.
• Transient solvers can be used to run at smaller timesteps and also can be run using varying time step.
• Since smaller timesteps are required to understand instabilities we will be using a transient solver.
• Also the transient solvers tends to achieve more accuracy for these type of simulations since it is a time marching solutions.
• The steady state solver ignores many cross terms and higher order terms with respect to time.

Atwood Number (A):

Atwood number (A) is a dimensionless number that describes density differences between two adjacent fluids with a common interface. This is an important parameter in the study of hydrodynamic instabilities indensity stratified flows. In Rayleigh-Taylor instability, the penetration distance of the heavy fluid bubbles into the lighter fluids is a function of acceleration time scale,

`A=(rho_1-rho_2)/(rho_1+rho_2)`

where,

`rho_1`= density of heavier fluid (kg/m^3)(i.e. density of water in our case= 997 kg/m^3)

`rho_2`= density of lighter fluid (kg/m^3)(i.e. density of air in our case= 1.225 kg/m^3)

Therefore in our case the atwood number = 0.9975.

How the variation in Atwood number affects the behaviour of the instability?

Since there are difference in densities between the mixing fluids, there are formation of spikes and bubbles. Bubble heads are shown to produce forward (reverse) transfer due to compression (extensional) straining at intermediate Atwood Number. At a higher atwood number the spikes are shown to produce forward transfer due to compression straining which is generated from larger spike pentration velocity. For our case we have a larger value of atwood number due to large difference in densities of the fluids. This results in formation of more bubbles and spikes.