Week 5 - Rayleigh Taylor Instability
Week 5 Rayleigh Taylor instability
AIM:
- To mention some practical CFD models that have been based on the mathematical analysis of Rayleigh Taylor waves and explain how these mathematical models have been adapted for CFD calculations.
- To perform the Rayleigh Taylor instability simulation for 2 different mesh sizes with the base mesh being 0.5 mm and compare the results by showing the animations.
- To run one more simulation with water and user-defined material (density = 400 kg/m3, viscosity = 0.001 kg/m-s) for refined mesh.
- To define the Atwood Number and find out the Atwood number for both the cases and explain how the variation in Atwood number in the above two cases affects the behavior of instability.
- What are some practical CFD models that have been based on the mathematical analysis of Rayleigh Taylor waves? In your own words, explain how these mathematical models have been adapted for CFD calculations.
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 behavior 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.
- Water suspended atop oil is an everyday example of Rayleigh–Taylor instability, and it may be modeled by two completely plane-parallel layers of immiscible fluid, the denser fluid on top of the less dense one and both subject to the Earth's gravity.
- 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.

Some practical CFD models that have been based on the mathematical analysis of Rayleigh Taylor waves:
The free surface flows are successfully modelled with one of the existing free surface models such as,
- Level - Set Methods (LSM)
Level-set methods are a conceptual framework for using level sets as a tool for numerical analysis of surfaces and shapes. The advantage of the level-set model is that one can perform numerical computations involving curves and surfaces on a fixed Cartesian grid without having to parameterize these objects (this is called the Eulerian approach). Also, the level-set method makes it very easy to follow shapes that change topology, for example, when a shape splits in two, develops holes, or the reverse of these operations. All these make the level-set method a great tool for modeling time-varying objects, like inflation of an airbag, or a drop of oil floating in water.
- Volume Of Fluid (VOF) method
In computational fluid dynamics, the volume of fluid (VOF) method is a free-surface modelling technique, i.e. a numerical technique for tracking and locating the free surface (or fluid–fluid interface). It belongs to the class of Eulerian methods which are characterized by a mesh that is either stationary or is moving in a certain prescribed manner to accommodate the evolving shape of the interface. As such, VOF is an advection scheme. The Navier–Stokes equations describing the motion of the flow have to be solved separately. The same applies for all other advection algorithms.
The method is based on the idea of a so-called fraction function C. It is a scalar function, defined as the integral of a fluid's characteristic function in the control volume, namely the volume of a computational grid cell. The volume fraction of each fluid is tracked through every cell in the computational grid, while all fluids share a single set of momentum equations, i.e. one for each spatial direction.
- Front Tracking Method:
The method is based on writing one set of governing equations for the whole computational domain and treating the different phases as one fluid with variable material properties. It is the use of surfaces or lower dimensional manifolds as computational degrees of freedom in a numerical algorithm. Its purpose is to improve the resolution of discontinuities or steep gradients in the solution variables or in the laws of physics which describe them.
Thermal or concentration discontinuities, and thermodynamic phase discontinuities, often poorly handled by Eulerian advection schemes, may benefit from the use of front tracking.
- Two – Fluid Model:
The two-fluid model is formulated by considering each phase separately. Thus, the model is expressed in terms of two sets of conservation equations governing the balance of mass, momentum and energy in each phase. However, since the averaged fields of one phase are not independent of the other phase, we have interaction terms appearing in these balance equations.
In the two-fluid model formulation, the transfer processes of each phase are expressed by their own balance equations, thus it is anticipated that the model can predict more detailed changes and phase interactions. However, this means that the two-fluid model is far more complicated.
The real importance of the two-fluid model is that it can take into account the dynamic and non-equilibrium interactions between phases. Thus, it is expected that the two-fluid model can be useful to the analyses of transient phenomena, wave propagations and of the flow regime changes. Particularly if the two phases are weakly coupled such that the inertia of each phase changes rapidly, the two-fluid model should be used to study these phenomena.
- k – L Turbulence Model:
The k−L turbulence model, where k is the turbulent kinetic energy and L represents the turbulent eddy scale length, is a two-equation turbulence model that has been proposed to simulate turbulence induced by Rayleigh-Taylor (RT) and Richtmyer Meshkov (RM) instabilities. The equations are based on the successful buoyancy-drag models for RT and RM flows, but constructed only with local parameters so that it can be applied to multidimensional flows with multiple shells of materials. The model has several unknown coefficients that are determined by comparing analytical and numerical solutions with RT and RM experiments.
SIMULATION PROCEDURE
In this challenge, we are going to simulate the Rayleigh Taylor Instability phenomenon of a multi – phase fluid for a coarse mesh and a refined mesh case and observe its flow properties.
We are going to set up 3 cases.
- Simulation for a base (coarse) mesh (element size = 0.5 mm)
- Simulation for a refined mesh (element size = 0.2mm)
- Simulation for water and user defined material.
IMPORTANT PARAMETERS
Surface Dimensions:20 mm x 20 mm
Fluid Material: Phase – 1 = Air (0) primary phase
Phase – 2 = Water (1) secondary phase
Solver: Pressure – based
Velocity Formulation: Absolute
Time Management: Transient state
Gravity: -9.81 (negative y – axis)
Viscous Model: Laminar
Multiphase: Volume of Fluids (VOF) method
Formulation: Implicit
Timestep Size:0.01
Number of Timesteps:500
- GEOMETRY
In this step, we are going to create a new ANSYS Space Claim geometry representing the two phases of the fluid, being air and water. The bottom layer of the geometry represents air and the top region represents water.
To form these layers, we create two boxes of size 20 mm x 20 mm each as shown below.


here we give the name of the 1st component as Air and create new component as water



confirming the share topology by clicking the design structure
We then use the pull tool to create a surface geometry of the same. Also, we use the share tool to share the topology between the two surfaces and obtain a conformal mesh. We thereby share 1 co-incident edge.
- MESHING
After the geometry step, we are going to create mesh for the created geometry.
CASE – 1: BASE MESH (element size = 0.5 mm)
In this case study, we are going to analyze the effects of base meshing on the simulation. We carry out the procedure in the following steps.
Here we have got only a small number of nodes (3321) elements (3200) because the geometry is too small.
- SETUP & SOLUTION
In this step, we set up the necessary parameters required to run the simulation.
We use a pressure – based solver with absolute velocity formulation.
We enable gravity on the negative y – axis with a standard g value of -9.81 m/s2.
We use a transient state time condition to run the simulation.
We use the laminar viscous model.
Fluent fluid materials – Air and Water (H2O)
We use a multiphase homogeneous model representing the volume of fluids (VOF)
2 Eulerian phases with implicit volume fraction formulation.

We consider the first phase (primary) as AIR and the second phase as WATER.

We standardize the simulation and create contours of volume fraction and solution animations for the same.
We patch the two surfaces together with the volume fraction value being ‘1’ for water and ‘0’ for air.
here to confirm the patching we are checking the contours

We run the simulation for 500 timesteps, for a timestep size of 0.01 to obtain the following results.
Results

1
2
Final contour
Animation

CASE – 2: REFINED MESH (element size = 0.2 mm)
In this case setup, we run the simulation for the same parameters by refining the mesh nodes (20301) elements (20000) size so as to observe the effects the meshing size has on the simulation.
For the same geometry, we obtain the refined meshing as shown below.

For this meshed geometry, we run the simulation with the same setup parameters for 500 timesteps to obtain the following outputs.
Results



Animation

CASE – 3: USER – DEFINED MATERIAL (REFINED MESH)
In this case setup, we are going to run the simulation for water and a user-defined material, for a refined mesh.
Specifications of the user – defined material:
Density ρ = 400 kg/m3
Viscosity μ = 0.001 kg/m−s
For this case study, we are going to set up a user-defined material along with water, for the same meshed geometry (element size 0.1 mm). The user-defined material setup is as follows.
As shown above, we create a user-defined material of density
400kg/m3 and viscosity 0.001kg/ms. Furthermore, we create the 2 Eulerian phases as water and user-defined material. We run the calculation for 500 timesteps for a timestep size of 0.001 to obtain the following results.
Results




Animation
google drive link for animations of 3 cases
https://drive.google.com/drive/folders/16qQxBkdgRoIkVt7V766mIIqmDAB2JcL5?usp=sharing
ATWOOD NUMBER
The Atwood number is a dimensionless number in fluid dynamics used in the study of hydrodynamic instabilities in density stratified flows. It is a dimensionless density ratio related to the density of heavier fluid and the density of lighter fluid.
It is a dimensionless density ratio defined as,
A=ρ1–ρ2/ρ1+ρ2
Were,
ρ1 = Density of heavier fluid
ρ2 = Density of lighter fluid.
Calculating the Atwood number for both the cases,
Case – 1 & 2
Density of heavier fluid (water)
ρ1 = 997 kg/m3
Density of lighter fluid (air)
ρ2 = 1.225 kg/m3
According to the definition, Atwood number is given by,
A=997–1.225 / 997+1.225
A = 0.99754
Case – 3
Density of heavier fluid (water)
ρ1= 997 kg/m3
Density of lighter fluid (user-defined material)
ρ2= 400 kg/m3
According to the definition, Atwood number is given by,
A=997–400 / 997+400
A = 0.4273
Atwood number is an important parameter in the study of Rayleigh–Taylor instability. For Atwood number close to 0, RT instability flows take the form of symmetric “fingers” of fluid; for Atwood number close to 1, the much lighter fluid “below” the heavier fluid takes the form of larger bubble-like plumes.
Atwood numbers have a direct effect on the growth of the perturbations. A higher Atwood number leads to a rapid or non-linear accelerating growth of the interface whereas for a lower Atwood number (closer to 0), leads to a steady growth with bubble like perturbation formation of the interface.
The behavior of RT Instability on the variation of Atwood number:
Rayleigh Taylor Instability 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 longer time. In contrast, for large Atwood numbers, bubble and mushroom head shapes are less symmetrical, also it occurs very rapidly at an accelerating rate which leads to a nonlinear growth rate.
Why is a steady – state approach not suitable for this simulation?
Steady state models are based on the assumptions that all flow conditions and properties of the system are constant with respect to time. Transient models, however, can handle conditions that change with time and assess the time-dependent impact on thermal and flow predictions. Furthermore, transient multiphase flow simulations may become relevant as the conditions change over time. In this case, a steady state solver can be used for the simulation if only the end result is required. But in this case, we need to study the characteristics of the fluid over a period of time and hence, a time marching solution would be more suitable for this simulation.
CONCLUSIONS
The transient – state simulation for all the three cases, representing the coarse mesh, refined mesh and for a user – defined material has been successfully simulated and the results have been obtained.
From the results, we can observe that a finer meshing results in a much clearer and more defined fluid interface, as it captures the details of the instability much better, as compared to a coarse mesh. It can also be noted that the computational time required for a refined mesh is higher.
For this simulation, a transient – state solver is preferred over a steady – state solver, as it involves variables and instabilities that are captured over a period of time. This gives us an accurate representation of the instabilities forming at the fluid interface in real – time. As a result of this simulation, the time required for the denser fluid to settle down is longer, as compared to that of the lighter fluid.
For the first two cases, we can observe that as the instability forms, finger like petrubations are formed, due to the value of the Atwood number being closer to 1 (0.99754 for cases 1 & 2). Whereas for the third case, we can observe that bubble like plumes are formed at the fluid interface due to the lower Atwood number value (0.4273). Therefore, we can conclude that as the Atwood number value increases (closer to 1), the interface growth is more rapid with finger shaped formations. Whereas, for a smaller Atwood number (closer to 0), we obtain a steady growth with bubble like formations. Therefore, as the density ratio increases, the interface instability due to the mixing is more rapid.