Designing An Experiment To Simulate A 2 Dimensional Rayleigh Taylor Instability Problem Using SpaceClaim & Ansys Fluent.

ABSTRACT:- The purpose of this paper is to simulate a two-dimensional Rayleigh-Taylor instability problem using the classical method of Finite Volume analysis of a multiphase model using ANSYS FLUENT 20. The study focuses on setting up an experiment to visualize the Rayleigh Taylor instability problem which is an instability in the interface between two fluids of different densities that occurs when the lighter fluid is pushing the heavier fluid in a gravitational field.  We start by designing the domain in SpaceClaim, which consists of two square boxes with dimensions of (20mm *20mm) one on top of the other separated by an "interface". We will be using this setup to perform a grid dependency test on the Rayleigh Taylor Instability problem for 3 different mesh sizes and compare the results for the variation in Atwood Number. The results were then validated with the observations from previously published literature. The simulation will be run using a pressure-based transient solver. 

________________________________________________________________________________________________________________________________________________________________________________

                                                                   keywords:-  Rayleigh Taylor Instability, Grid Dependency Test, Transient state, Multiphase Flow, Atwood Number

________________________________________________________________________________________________________________________________________________________________________________

  1. Theory

 The Rayleigh–Taylor instability or RT instability is the instability of an interface between two fluids of different densities (ρ) which occurs when the lighter fluid is pushing the heavier fluid under the influence of the earth's gravity. In Multiphase flow, dynamic variables like velocity, viscosity, pressure, and density are generally used to explain the movement of fluids under the gravitational field. Multiphase flows occur in nuclear geothermal power plants, food processing industries, and many others. One of the basic examples is Rayleigh Taylor Instability. For our experiment, we used air & water with densities of 1.225 `(kg)/m^3` &  998 `(kg)/m^3` respectively. Since water has a higher density it would be suspended atop the air domain & both fluids will be separated by an interface (virtual diaphragm). The interface (the boundary between two spatial regions occupied by different matter) or diaphragm here is unstable to any disturbances of the interface if some portion of the heavier fluid is displaced downward with an equal volume of lighter fluid displaced upwards, the potential energy of the configuration is lower than the initial state. Thus the disturbance will grow and lead to a further release of potential energy, as the more dense material moves down under the (effective) gravitational field, and the less dense fluid is further displaced upwards. 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). It is a dimensionless density ratio & can be defined as: 

                                                `=> (water-air) / (water + air)` `=> (998-1.225)/ (998+1.225) => 0.9975`, where   = density of the heavier fluid,  = density of the lighter fluid

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. The Atwood number (A) is a dimensionless number in fluid dynamics used in the study of hydrodynamics instabilities in density stratified flows.

  2. Geometry 

The 2D geometry of the domain was modeled using SpaceClaim. We start by designing two squared boxes with dimensions of (20mm*20mm) placed one above the other on the XY plane. Next, we use the "Pull" tool feature to convert the geometry into surfaces. The bottom box named "air" the top box named "water" as shown in Fig 2.1. Since one edge of the box is shared by both the box we use the "share topology" feature & set it to "share". This concludes the geometry designing part.

                                                                                                                                Fig 2.1:- Computational Domain 

  3. Mesh Generation

We will be running 3 simulations & each case will vary based on the mesh size. Case 1 will have the largest element size with the least number of mesh cells generated Case 3 will have the smallest element size with the maximum number of mesh cells generated. This process of refining the mesh quality is called Grid Dependency    Test. In a grid dependency test we begin with a coarse mesh & then slowly improve the mesh quality by refining it till a point where the simulation results don't vary drastically. A high-density mesh usually produces accurate results, however, refining the mesh too much requires high computational power with long run times. 

                                                         Case 1                                                                    Case 2                                                                      Case 3

  4. Setting Up The Flow Physics

In Ansys Fluent after creating the mesh & assigning the boundaries, the next stage of the process is to set up the simulation parameters. This step involves specifying certain properties to capture the physics of the problem perfectly, such as materials used, simulation time parameters, solver parameters, initial & boundary conditions to solve the Navier Stokes equation which is PDEs, defining body forces, type of fluid flow & species involved.

We start by setting up the Domain, here we will be using a pressure-based transient solver with absolute velocity formation. We can justify using a transient solver since we are interested in visualizing the interface behavior with function of time. We will also enable "gravity", which will be acting along the Y-axis. The next parameter that we need to set up is the Flow Physics, here we will be navigating to the viscous feature & use a laminar flow since the velocity magnitude will be very small.

For materials, along with "air" at the bottom phase, we will be adding "water" for the top phase. Since the domain consists of 2 phases, we will be making use of a "multiphase" feature where we will be running the simulation with the Volume of Fluid (VoF) method. Next, we assign the phases to corresponding material i.e phase 1 = air & phase2= water. The next step involves "patching" the water phase with the water surface & assign a value of 1 (this value represents the phase1) & for the air surface, we assign a value of 0 ( represents the phase 2).  

The final step involved the Post-processing features where we created the phase contours & animations & ran the simulation for a time step size of 0.005 & the number of time steps will vary for each case. 

  5. Results

     (A) Initial Conditions

                                                Case 1                                                                                    Case 2                                                                                   Case 3

     (B) Phase Disturbance Contour 

                                              Case 1                                                                                    Case 2                                                                                   Case 3

     (C) Animation 

     (D) Final Observations

We can observe from the attached pictures, Case 1 as expected had the least number of cells since it had the largest mesh size while Case 3 had the maximum number of cells generated since it had the smallest mesh size. Although we observed bubbles & spikes for all 3 cases but the finer mesh was able to capture the physics more accurately & had a higher resolution for visualizing purposes as compared to the other two cases. While running the simulation we also observed if we entered an incorrect time step it leads to "Floating point error" which can occur due to many factors such as, improper mesh size, defining some property close to zero. Furthermore, as we calculated the Atwood Number earlier, our value for (A) was very close to 1, which indicated that the fluid with lighter density will take the form of large bubble-like plumes & that is exactly what we observed from our experiment for all 3 d cases irrespective of the mesh size parameters. 

   6. Results (Part 2)

 We keep all the parameters the same except add air as a "user-defined" material with modified values of ρ=400 kg/m³ & Viscosity(V) = 0.001 kg/m-s. 

Atwood Number (A)= 

    `=> (water-air) / (water + air)` `=> (998-400)/ (998+400) => 0.4277` [where   = density of the heavier fluid,  = density of the lighter fluid]

________________________________________________________________________________________________________________________________________________________________________________