Rayleigh Taylor Instability
In this tutorial, the setup of inviscid simulations of Rayleigh-Taylor instability in a two-dimensional domain are demonstrated. The following capabilities of HAMeRS will be shown:
- Conservative four-equation multi-species flow model
- Sixth-order shock capturing method WCNS6-LD
- Gravity as a source term
Rayleigh-Taylor instability (RTI) is a hydrodynamic instability that happens on an interface between fluids with different densities when an acceleration is directed from the heavier fluid to the lighter fluid. This tutorial is a single-mode Rayleigh-Taylor instability problem with domain size [0, 1/4] x [0, 1]. The initial interface is at y = 1/2. The heavier fluid with density ρ0 = 2 and specific heat capacity ratio γ0 = 5/3 is below the interface and the lighter fluid with density ρ1 = 1 and specific heat capacity ratio γ1 = 1.4 is above the interface. The gravitational acceleration (and body force) is pointing in the positive y direction. Pressure p is continuous across the interface and a small perturbation is given to the y-component of velocity, v. The initial conditions of the problem is:
The initial perturbed interface is located at
The file containing the initial conditions of this problem (RayleighTaylorInstability2D.cpp) can be found in problems/Euler/initial_conditions and the initial conditions should be set up by ln -sf <absolute path to RayleighTaylorInstability2D.cpp> src/apps/Euler/EulerInitialConditions.cpp. The code has to be re-compiled after the link to the actual initial condition file is set.
The configurations of the input file are discussed in this section. Only settings of several important blocks are discussed. The first thing to choose in the input file is whether it is an Euler (inviscid) or Navier-Stokes (diffusive and viscous) application. Since our problem is inviscid, we first set the application type to be Euler:
Application = "Euler"
The next thing to set is the block Euler. The key part is to set has_source_terms to be TRUE and specify the gravity vector in the block Source_terms:
Euler
{
// Name of project
project_name = "2D Rayleigh-Taylor instability"
// Number of species
num_species = 2
// Flow model to use
flow_model = "FOUR_EQN_CONSERVATIVE"
Flow_model
{
// Equation of state to use
equation_of_state = "IDEAL_GAS"
Equation_of_state_mixing_rules
{
species_gamma = 1.66666666667, 1.4
species_R = 1.0, 1.0
}
has_source_terms = TRUE
Source_terms
{
has_gravity = TRUE
gravity = 0.0, 1.0
}
}
// Convective flux reconstructor to use
convective_flux_reconstructor = "WCNS6_LD_HLLC_HLL"
Convective_flux_reconstructor{}
Boundary_data
{
// Set the boundary conditions for edges
boundary_edge_xlo
{
boundary_condition = "REFLECT"
}
boundary_edge_xhi
{
boundary_condition = "REFLECT"
}
boundary_edge_ylo
{
boundary_condition = "DIRICHLET"
partial_densities = 2.0, 0.0
velocity = 0.0, 0.0
pressure = 1.0
}
boundary_edge_yhi
{
boundary_condition = "DIRICHLET"
partial_densities = 0.0, 1.0
velocity = 0.0, 0.0
pressure = 2.5
}
// Set the boundary conditions for nodes
boundary_node_xlo_ylo
{
boundary_condition = "XREFLECT"
}
boundary_node_xhi_ylo
{
boundary_condition = "XREFLECT"
}
boundary_node_xlo_yhi
{
boundary_condition = "XREFLECT"
}
boundary_node_xhi_yhi
{
boundary_condition = "XREFLECT"
}
}
}
Dirichlet boundary conditions are specified for the top and bottom boundaries. The values of the primitive variables at the boundaries are set. The boundary conditions of the left and right boundaries are set to "REFLECT".
The settings for Main is very similar to those of tutorial 1 and are not discussed.
The following settings are used for CartesianGeometry:
CartesianGeometry
{
// Lower and upper indices of computational domain
domain_boxes = [(0, 0), (59, 239)]
x_lo = 0.0, 0.0 // Lower end of computational domain
x_up = 0.25, 1.0 // Upper end of computational domain
// Periodic dimension. A non-zero value indicates that the direction is periodic
periodic_dimension = 0, 0
}
No settings are used for ExtendedTagAndInitialize since this is not an AMR simulation:
ExtendedTagAndInitialize
{
}
The max_level in PatchHierarchy is set to 1 since this is not an AMR simulation:
PatchHierarchy
{
// Maximum number of levels in hierarchy
max_levels = 1
ratio_to_coarser
{
// Vector ratio to next coarser level
level_1 = 2 // all finer levels will use same values as level_1...
}
largest_patch_size
{
level_0 = 1000, 1000 // all finer levels will use same values as level_0...
}
smallest_patch_size
{
level_0 = 4, 4 // all finer levels will use same values as level_0...
}
}
The following settings are used for RungeKuttaLevelIntegrator and TimeRefinementIntegrator:
RungeKuttaLevelIntegrator
{
cfl = 0.5e0 // Max cfl factor used in problem
cfl_init = 0.5e0 // Initial cfl factor
lag_dt_computation = FALSE
use_ghosts_to_compute_dt = TRUE
}
TimeRefinementIntegrator
{
start_time = 0.0e0 // Initial simulation time
end_time = 1.95e0 // Final simulation time
grow_dt = 1.0e0 // Growth factor for timesteps
max_integrator_steps = 10000 // Max number of simulation timesteps
}
To run the simulation with one core, first put the input file inside a directory named tutorial_06 under HAMeRS. Then, execute the built main executable inside build/src/exec/main with the input file:
../build/exec/main <input filename>
To run the simulation with multiple cores, you can try mpirun/mpiexec/srun, depending on the MPI library used for the compilation of HAMeRS. For example, if mpirun is used:
mpirun -np <number of processors> ../build/src/exec/main <input filename>
The data can be visualized by opening dumps.visit inside the viz folder with VisIt.The figure below shows the density field:
